Li Qiao
email: qiaoli@bit.edu.cn

Jun Zhang
email: buaazhangjun@vip.sina.com

Zhen Gao
email: gaozhen16@bit.edu.cn

Derrick Wing

Fellow, IEEE Kwan Ng

Fellow, IEEE Marco Di Renzo
email: marco.di-renzo@universite-paris-saclay.fr

Fellow, IEEE Mohamed-Slim Alouini
email: slim.alouini@kaust.edu.sa

J Qiao

Z Zhang

Gao

D W K Ng

M.-S Alouini

Massive Access in Media Modulation Based Massive Machine-Type Communications

Keywords: Massive access, media modulation, massive multi-input multi-output, massive machine-type communications

The massive machine-type communications (mMTC) paradigm based on media modulation in conjunction with massive multi-input multi-output base stations (BSs) is emerging as a viable solution to support the massive connectivity for the future Internet-of-Things, in which the inherent massive access at the BSs poses significant challenges for device activity and data detection (DADD). This paper considers the DADD problem for both uncoded and coded media modulation based mMTC with a slotted access frame structure, where the device activity remains unchanged within one frame. Specifically, due to the slotted access frame structure and the adopted media modulated symbols, the access signals exhibit a doubly structured sparsity in both the time domain and the modulation domain. Inspired by this, a doubly structured approximate message passing (DS-AMP) algorithm is proposed for reliable DADD in the uncoded case. Also, we derive the state evolution of the DS-AMP algorithm to theoretically characterize its performance. As for the coded case, we develop a bit-interleaved coded media modulation scheme and propose an iterative DS-AMP (IDS-AMP) algorithm based on successive inference cancellation (SIC), where the signal components associated with the detected active devices are successively subtracted to improve the data decoding performance. In addition, the channel estimation problem for media modulation based mMTC is discussed and an efficient data-aided channel state information (CSI) update strategy is developed to reduce the training overhead in block fading channels. Finally, simulation results and computational complexity analysis verify the superiority of the proposed DS-AMP algorithm over state-of-the-art algorithms in the uncoded case. Also, our results confirm that the proposed SIC-based IDS-AMP algorithm can enhance the data decoding performance in the coded case and verify the validity of the proposed data-aided CSI update strategy.

I. INTRODUCTION

W ITH the advent of the Internet-of-Things (IoT), mas- sive machine-type communications (mMTC) are expected to support various potential applications, including smart metering, surveillance and healthcare. In fact, mMTC are considered to be an indispensable component for future beyond 5G/6G networks [START_REF] Shariatmadari | Machine-type communications: Current status and future perspectives toward 5G systems[END_REF], [START_REF] Wong | Key Technologies for 5G Wireless Systems[END_REF]. In stark contrast with conventional human-type communications (HTC) that are characterized by downlink transmissions with long packets and high data rates, mMTC are characterized by uplink transmissions with short packets from massively deployed machine type devices (MTDs) whose data traffic is sporadic [START_REF] Durisi | Toward massive, ultrareliable, and low-latency wireless communication with short packets[END_REF], [START_REF] Bockelmann | Massive machine-type communications in 5G: Physical and MAC-layer solutions[END_REF]. The unique traffic characteristics of mMTC indicate that the existing access solutions designed for HTC can not effectively support the massive access of MTDs, which necessitates lowlatency and high-reliable massive access techniques with lowcomplexity signal processing algorithms.

A. Related Work

Existing access solutions can be mainly divided into two categories: grant-based access approaches and grant-free access approaches [START_REF] Shariatmadari | Machine-type communications: Current status and future perspectives toward 5G systems[END_REF]- [START_REF] Bockelmann | Massive machine-type communications in 5G: Physical and MAC-layer solutions[END_REF]. As for grant-based access solutions, allocating orthogonal radio resources to different active MTDs via some sophisticated scheduling algorithms is necessary before the uplink transmission. However, this costs extra signaling overhead for handling the granted signals. In practice, due to the limited orthogonal resources, it would be difficult to support a large number of MTDs by applying grant-based schemes, and the use of complicated resource scheduling algorithms can lead to a prohibitive signaling overhead and latency, leading even to congestion [START_REF] Shariatmadari | Machine-type communications: Current status and future perspectives toward 5G systems[END_REF]- [START_REF] Chen | Fully nonorthogonal communication for massive access[END_REF]. As an alternative solution, grant-free access schemes have recently emerged and have drawn significant attention [START_REF] Chen | Massive access for 5G and beyond[END_REF]- [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF]. In grant-free access systems, the MTDs can transmit in the uplink without waiting for permission, which significantly simplifies the uplink access procedure and hence reduces the access latency compared with grant-based access protocols.

Due to the limited radio resources but large number of potential MTDs, grant-free massive access schemes result in non-orthogonal transmissions, which makes the device activity and data detection (DADD) problem more difficult to handle. In particular, the DADD can be treated as an underdetermined linear problem. This implies that the application of classical linear least-squares (LS) and linear minimum mean square error (LMMSE) detectors to tackle the DADD problem would result in poor detection performance [START_REF] Shim | Multiuser detection via compressive sensing[END_REF]. Fortunately, thanks to the sporadic traffic characteristics of MTDs in mMTC, the number of active MTDs is usually much smaller than the number of total MTDs in any given time interval, which mo-antenna elements for each extra spatial modulated data bit [START_REF] Di Renzo | Spatial modulation for multipleantenna wireless systems: A survey[END_REF].

Fortunately, during the last years many improved transmission schemes have been proposed in order to enhance the spectral efficiency of spatial modulation without compromising its low-complexity and energy-efficiency. A recent overview and comparison of the most popular solutions, which includes generalized spatial modulation, media modulation, and, more recently, metasurface-based modulation, is available in [START_REF] Li | Single-RF MIMO: From spatial modulation to metasurface-based modulation[END_REF]. In this paper, we consider media modulation introduced in [START_REF] Khandani | Media-based modulation: A new approach to wireless transmission[END_REF]. In particular, media modulation employs a single RF chain, a single radiating element, and several low-cost RF mirrors [START_REF] Naresh | On media-based modulation using RF mirrors[END_REF]- [START_REF] Mao | Secure single-input-multipleoutput media-based modulation[END_REF]. The information bits are encoded into the active/inactive (or ON/OFF) status of the RF mirrors, which determines the resulting radiation pattern of the entire structure [START_REF] Naresh | On media-based modulation using RF mirrors[END_REF], [START_REF] Basar | Media-based modulation for future wireless systems: A tutorial[END_REF]. In contrast to spatial modulation, the number of spatial bits encoded in media modulation is larger and depends on the number of distinguishable radiation patterns that can be realized [START_REF] Basar | Index modulation techniques for next-generation wireless networks[END_REF]. Due to the promising advantages of media modulation, it was recently adopted in [START_REF] Qiao | Compressive sensing based massive access for IoT relying on media modulation aided machine type communications[END_REF], [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF], [START_REF] Zhang | Low-complexity multi-user detection for MBM in uplink large-scale MIMO systems[END_REF], [START_REF] Shamasundar | Media-based modulation for the uplink in massive MIMO systems[END_REF] for application to the uplink transmission of MTDs with mMIMO BSs. Specifically, as for grant-based access schemes, the authors of [START_REF] Zhang | Low-complexity multi-user detection for MBM in uplink large-scale MIMO systems[END_REF] and [START_REF] Shamasundar | Media-based modulation for the uplink in massive MIMO systems[END_REF] proposed an iterative interference cancellation detector and a message passing algorithm based detector, respectively. As for grant-free access schemes with a slotted access frame structure, a structured orthogonal matching pursuit (StrOMP) algorithm for activity detection and a successive interference cancellation (SIC)-based structured subspace pursuit (SSP) algorithm for data detection were proposed in [START_REF] Qiao | Compressive sensing based massive access for IoT relying on media modulation aided machine type communications[END_REF]. In addition, a prior-information aided media modulation subspace matching pursuit (PIA-MSMP) algorithm was proposed in [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF] for DADD. Additionally, the authors of [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF] proposed a prior-information aided adaptive media modulation subspace matching pursuit algorithm to accommodate the dynamic device activity. Although the structured sparsity of media modulated symbols was exploited for better performance in [START_REF] Qiao | Compressive sensing based massive access for IoT relying on media modulation aided machine type communications[END_REF] and [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF], the proposed greedy algorithms failed to fully exploit the finite alphabet constraint of the transmit symbols. Hence, there are plenty of opportunities for further improving the design and optimization of media modulation for application to grant-free massive access. It is important to note, in addition, that the grant-free massive access schemes in the aforementioned research works, i.e., [START_REF] Wang | Joint user activity and data detection based on structured compressive sensing for NOMA[END_REF]- [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF], [START_REF] Gao | Compressive-sensing based multiuser detector for the large-scale SM-MIMO uplink[END_REF], [START_REF] Meng | Multi-user detection for spatial modulation via structured approximate message passing[END_REF], [START_REF] Zhang | Low-complexity multi-user detection for MBM in uplink large-scale MIMO systems[END_REF], [START_REF] Shamasundar | Media-based modulation for the uplink in massive MIMO systems[END_REF], are developed for an uncoded scenario, while their extension to coded systems needs further research.

B. Our Contributions

This paper considers both uncoded and coded media modulation based mMTC, by assuming a slotted access frame structure. It is assumed that the activity of the MTDs remains unchanged in each frame, which yields a sparse transmission frame. To design a reliable DADD scheme, we introduce a doubly structured AMP (DS-AMP) algorithm for uncoded transmission. The state evolution (SE) of the proposed DS-AMP algorithm is derived in order to characterize its performance. As for the coded transmission, we introduce a bit-interleaved coded media modulation (BICMM) scheme and propose an SIC-based iterative DS-AMP (IDS-AMP) scheme for improving the decoding performance. In addition, the channel estimation (CE) problem for media modulation based mMTC is discussed and an effective channel state information (CSI) update strategy is developed to reduce the training overhead in block fading channels. With the aid of numerical simulations and computational complexity analysis, we verify the superiority of the proposed DS-AMP algorithm and the SIC-based IDS-AMP algorithm for uncoded and coded transmission, respectively, with respect to state-of-theart benchmark schemes. Moreover, the effectiveness of the CSI update strategy is verified. The main contributions of this paper can be summarized as follows:

• DS-AMP algorithm for uncoded media modulation based mMTC: By utilizing the structured sparsity of media modulated symbols in the modulation domain and the discrete distribution of quadrature amplitude modulation (QAM) alphabets, the proposed DS-AMP algorithm can reliably perform DADD. Furthermore, by leveraging the structured sparsity of the slotted access frame structure in the time domain, the active/inactive status of the MTDs and the noise variance can be adaptively learned via the EM algorithm with enhanced accuracy. Besides, we derive the theoretical SE, which closely matches the simulated results for the DS-AMP algorithm. • BICMM designed for coded media modulation based mMTC: We integrate the bit-interleaved coded modulation (BICM) [START_REF] Caire | Bit-interleaved coded modulation[END_REF], [START_REF] Koca | Bit-interleaved coded spatial modulation[END_REF] into media modulation and develop a BICMM scheme, which can effectively mitigate the error bursts in fading channels and improve the data decoding performance. Particularly, we employ a bit-wise interleaver between the channel encoder and the media modulation module, where the bits of the media modulated symbols and the QAM symbols are collectively processed by channel coding and interleaving. • SIC-based IDS-AMP scheme for coded media modulation based mMTC: We design a dedicated data packet structure for MTDs and develop an IDS-AMP detector at the BSs. The developed data packet includes a signature sequence part that is known by the transceiver and a payload part for data transmission. Particularly, we consider the decoding quality of the signature sequence as a metric to determine the SIC order. Different from the DS-AMP algorithm that yields hard-decision estimates, the developed IDS-AMP detector yields soft log-likelihood ratio (LLR) for both the media modulated and the QAM symbol bits, which are fed to the channel decoding module for improving the decoding accuracy. Furthermore, the SIC operation in the IDS-AMP detector successively subtracts the signal components associated with the well decoded active devices to improve the data decoding performance.

Notation: Matrices and vectors are denoted by symbols in boldface. Operators (•) II. SYSTEM MODEL In this section, we first introduce the media modulation based mMTC scheme and then focus on the massive access at the BS for DADD.

A. Media Modulation Based mMTC Scheme

As illustrated in Fig. 1, we consider K MTDs that employ media modulation for enhancing the throughput and the BS employs mMIMO with N r ≫ 1 antenna elements for realizing a reliable massive access. Only K a out of K (K ≫ K a ) MTDs are active simultaneously. Specifically, each MTD is equipped with one RF chain, one transmit antenna, and N RF low-cost RF mirrors, where each RF mirror has a controllable binary ON/OFF status. Each device has N t = 2 NRF different kinds of mirror activation patterns (MAPs), i.e., N t different channel realizations, which can be used to encode log 2 N t = N RF extra information bits. As for the uplink access, therefore, the throughput is η = log 2 M + N RF bit per channel use (bpcu), which consists of log 2 M bits conveyed by, e.g., an M -ary QAM (M -QAM) scheme and N RF bits conveyed by media modulation. Furthermore, the use of an mMIMO offers improved detection performance for massive access [START_REF] Gao | Compressive-sensing based multiuser detector for the large-scale SM-MIMO uplink[END_REF], [START_REF] Liu | Massive connectivity with massive MIMO-Part I: Device activity detection and channel estimation[END_REF].

B. Massive Access in Media Modulation Based mMTC

As shown in Fig. 2, we consider a media modulation based mMTC system with a slotted access frame structure, where the active/inactive status of the K MTDs is invariable in a frame, i.e., J successive time slots. We refer to this property as the structured sparsity in the time domain. For convenience, we focus our attention on the massive access problem in any given frame without loss of generality. Specifically, the received signal at the BS in the j-th, ∀j ∈ [J], time slot, which is denoted by y j ∈ C Nr×1 , can be written as

y j = K k=1 a k s k,j H k d k,j +w j = K k=1 H k x k,j +w j = Hx j + w j , (1)
where the binary activity indicator in the j-th time slot is selected from the M -QAM set S, d k,j ∈ C Nt×1 and x k,j = a k s k,j d k,j ∈ C Nt×1 are the media modulated symbol and the effective uplink transmitted symbols of the k-th MTD in the j-th time slot, respectively, H k ∈ C Nr×Nt is the multiple-input multiple-output (MIMO) channel matrix corresponding to the k-th MTD, w j ∈ C Nr×1 is the Gaussian noise whose elements follow an independent and identically distributed (i.i.d.) complex Gaussian distribution CN ([w j ] n ; 0, KNt) and xj = [(x 1,j ) T , (x 2,j ) T , ..., (x K,j ) T ] T ∈ C (KNt)×1 are the aggregated channel matrix and uplink access signal of the j-th time slot, respectively.

a k ∈ {0, 1} is one (zero) if the k-th MTD is active (inactive), s k,j of the k-th MTD
σ 2 w ), ∀n ∈ [N r ], H = [H 1 , H 2 , ..., H K ] ∈ C Nr×(
According to the media modulation transmission scheme, only one entry of the media modulated symbol d k,j , ∀j ∈ [J] and ∀k ∈ [K], is one and the others are zeros, i.e.,

supp {d k,j } ∈ [N t ], d k,j 0 = 1, d k,j 2 = 1, (2)
where supp{d k,j } denotes the support set of d k,j . We refer to this property as the structured sparsity in the modulation domain. Additionally, we consider the commonly adopted Gauss-Markov block fading channel model [START_REF] Xiang | Deep-Learning-Aided joint channel estimation and data detection for spatial modulation[END_REF], [START_REF] Ishikawa | Differential-detection aided large-scale generalized spatial modulation is capable of operating in high-mobility millimeter-wave channels[END_REF], and the details of the CE problem will be discussed in Section V.

III. PROPOSED SOLUTION FOR UNCODED MEDIA MODULATION BASED MMTC

In this section, we present the proposed DS-AMP algorithm for DADD in uncoded media modulation based mMTC. First, we exploit the doubly structured sparsity in media modulation based mMTC and formulate an optimization problem for massive access. Then, by exploiting the doubly structured sparsity as a priori information, we propose a DS-AMP algorithm for effective DADD. Subsequently, the SE of the DS-AMP algorithm is derived to theoretically predict the performance. Finally, the computational complexity of the DS-AMP algorithm is analyzed.

A. Preliminaries

Although the total number of MTDs in the IoT is generally massive, their sporadic traffic behavior results in a sparse activity, i.e., at any given time slot, only a minority of MTDs are active [START_REF] Shariatmadari | Machine-type communications: Current status and future perspectives toward 5G systems[END_REF]- [START_REF] Bockelmann | Massive machine-type communications in 5G: Physical and MAC-layer solutions[END_REF]. Hence, we introduce an activity indicator vector a = [a 1 , a 2 , ..., a K ] T ∈ C K×1 , which is sparse as the number of active MTDs K a = a 0 ≪ K. Furthermore, we collectively refer to the structured sparsity in the time domain due to the slotted access frame structure and the structured sparsity in the modulation domain shown in (2) as the doubly structured sparsity.

To exploit the structured sparsity in the time domain, we rewrite the received signals of J successive time slots in a compact matrix form as

Y = HX + W, (3)
where we have Y = [y 1 , y 2 , ..., = min The optimization problem in (4) minimizes the mean square error between Y and HX, which is equivalent to estimating the a posteriori mean of the uplink access signal X [45] 1 . In (4), the a posteriori mean of

y J ] ∈ C Nr×J , H ∈ C Nr×(KNt) , X = [x 1 ,
{a k ,d k,j ,s k,j } J,K j=1,k=1 J j=1 y j - K k=1 a k s k,j H k d k,j 2
[x k,j ] i , ∀k ∈ [K], ∀j ∈ [J], ∀i ∈ [N t ],
can be expressed as

[ x k,j ] i = [x k,j ] i ∈S [x k,j ] i p [x k,j ] i |y j , (5)
where S = {S, 0}, p [x k,j ] i |y j is the marginal distribution of p (x j |y j ) that can be expressed as

p [x k,j ] i |y j = ∼{[x k,j ] i } p (x j |y j ) . (6)
Based on Bayes' theorem, the joint posterior distribution p (x j |y j ) can be expressed as

p xj |y j ; σ 2 w , a = p y j |x j ; σ 2 w p (x j ; a) p (y j ) = 1 p (y j ) Nr n=1 p [y j ] n |x j ; σ 2 w K k=1 p (x k,j ; a k ) , (7)
1 Due to space limitations, interested readers are referred to Eqs. (10.2)-(10.5) of Chapter 10 in [START_REF] Kay | Fundamentals of Statistical Signal processing, Volume I: Estimation Theory[END_REF] for the detailed derivations. where the likelihood function can be expressed as

p [y j ] n |x j ; σ 2 w = 1 πσ 2 w exp - 1 σ 2 w [y j ] n - K k=1 [H k x k,j ] n 2 .
(8) According to the sparse device activity, the structured sparsity in the modulation domain in (2), and the discrete distribution of the QAM alphabet, the a priori distribution p (x k,j ; a k ) in ( 7) is formulated as

p (x k,j ; a k ) = (1 -a k ) Nt i=1 δ [x k,j ] i + a k    1 N t Nt i=1   1 M s∈S δ [x k,j ] i -s g∈[Nt],g =i δ [x k,j ] g      , (9)
where M = |S| c and δ (•) is the Dirac delta function.

A factor graph representation of the joint posterior distribution p xj |y j ; σ2 w , a in [START_REF] Chen | Massive access for 5G and beyond[END_REF] is illustrated in Fig. 3. However, calculating the marginal distribution p [x k,j ] i |y j , ∀k ∈ [K], ∀j ∈ [J], and ∀i ∈ [N t ], from the joint posterior distribution p xj |y j ; σ 2 w , a is extremely complicated in massive access, due to the exceedingly large value of K. Fortunately, the AMP algorithm can provide an effective approximation of the marginal distributions at a low complexity, yet achieving near MMSE performance. In particular, by applying the AMP algorithm on the factor graph illustrated in Fig. 3, we propose a DS-AMP algorithm for handling the massive access problem. Specifically, we apply the DS-AMP algorithm to calculate the a posteriori mean of the media modulated signals x k,j (k ∈ [K], j ∈ [J]), and resort to the EM algorithm for estimating the activity indicators a k and the variance σ 2 w of the complex Gaussian noise.

2) Update Rules of the DS-AMP Algorithm: As discussed in [START_REF] Donoho | Message-passing algorithms for compressed sensing[END_REF]- [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF], in the large system regime 2 (i.e., K → ∞, λ = Ka K and κ = Nr K are fixed), the AMP algorithm decouples, in the asymptotic regime, the matrix estimation problem of (3) into KJN t uncoupled scalar problems

Y = HX + W → r l,j = [x k,j ] i + w l,j , ∀i, j, k, (10)
where l = (k-1)N t +i, r l,j is the mean of [x k,j ] i estimated by the AMP algorithm, w l,j ∼ CN ( w l,j ; 0, φ l,j ) is the equivalent noise, and φ l,j is its variance.

In addition, the joint posterior distribution ( 7) can be approximated as

p xj |y j ; σ 2 w , a ≈ q xj |y j ; σ 2 w , a = K k=1 Nt i=1 q [x k,j ] i |r l,j , φ l,j ; σ 2 w , a k . ( 11
)
Based on Bayes' theorem, the approximated marginal a posteriori distribution, denoted as q [x k,j ] i |r l,j , φ l,j ; σ 2 w , a k , ∀i, j, k, is given as follows

q [x k,j ] i |r l,j , φ l,j ; σ 2 w , a k = 1 q (r l,j ; σ 2 w , a k ) q r l,j | [x k,j ] i ; σ 2 w p [x k,j ] i ; a k , (12)
where

q r l,j | [x k,j ] i = 1 πφ l,j exp - 1 φ l,j r l,j -[x k,j ] i 2 , (13)
p [x k,j ] i ; a k = ∼{[x k,j ] i } p (x k,j ; a k ) = 1- a k N t δ [x k,j ] i + a k N t M s∈S δ [x k,j ] i -s , (14)
q r l,j ; σ 2 w , a k = [x k,j ] i ∈S q r l,j | [x k,j ] i ; σ 2 w p [x k,j ] i ; a k . (15)
Then, the a posteriori mean and variance of [x k,j ] i , ∀i, j, k, are given, respectively, as

[ x k,j ] i = f m (r l,j , φ l,j ) = [x k,j ] i ∈S [x k,j ] i q [x k,j ] i |r l,j , φ l,j ; σ 2 w , a k , (16)
[ v k,j ] i = f v (r l,j , φ l,j ) = -[ x k,j ] i 2 + [x k,j ] i ∈S [x k,j ] i 2 q [x k,j ] i |r l,j , φ l,j ; σ 2 w , a k , (17)
where l = (k -1)N t + i.

It is worth noting that r l,j , φ l,j , x k,j , and v k,j are updated iteratively by the AMP algorithm. In the factor graph in Fig. 3, in particular, r l,j and φ l,j , ∀l, j, are updated iteratively at the variable nodes. The update rules at the t-th iteration are expressed as

φ t l,j = Nr n=1 H [n,l] 2 σ 2 w + V t n,j -1 , (18)
r t l,j = x t k,j i + φ t l,j Nr n=1 H * [n,l] [y j ] n -Z t n,j σ 2 w + V t n,j , (19)
where V t n,j and Z t n,j , ∀n, j, are updated at the factor nodes of the factor graph as

V t n,j = K k=1 H k[n,:] 2 v t k,j , (20)
Z t n,j = K k=1 H k[n,:] x t k,j -V t n,j [y j ] n -Z t-1 n,j σ 2 w + V t-1 n,j . (21)
Due to space limitations, interested readers are referred to [START_REF] Meng | An expectation propagation perspective on approximate message passing[END_REF], [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF] for the detailed derivation of the AMP update rules in ( 18)- [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF].

3) Parameters Estimation: It is noteworthy that the activity indicator a k , ∀k ∈ [K], and the noise variance σ 2 w in ( 16) and ( 17) are unknown parameters to be determined. Hence, we will integrate the EM algorithm to adaptively learn those parameters in each DS-AMP iteration.

Note that the EM algorithm is an iterative approach that finds the maximum likelihood solutions for probabilistic models with unknown parameters [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF]. By defining θ = σ 2 w , a k , k ∈ [K] , the EM algorithm updates the parameter set θ as follows

Q θ, θ t = E ln p (X, Y; θ) |Y; θ t , (22)
θ t+1 = arg max θ Q θ, θ t , (23)
where θ t is the parameter set estimated at the t-th iteration, E {•|Y; θ t } denotes the expectation conditioned on the received signal Y under θ t . According to the AMP algorithm, the a posteriori distribution p xj |y j ; σ 2 w , a can be approximated as q xj |y j ; σ 2 w , a , which reduces the complexity of the EM estimation. Hence, we can obtain the update rules of the noise variance σ 2 w and the activity indicator a k , ∀k, as follows

σ 2 w t+1 = 1 JN r J j=1 Nr n=1    [y j ] n -Z t n,j 2 1 + V t n,j (σ 2 w ) t 2 + σ 2 w t V t n,j V t n,j + (σ 2 w ) t    , (24)
a t+1 k = f a r t l,j , φ t l,j ; a t k = 1 J J j=1 x k,j ∈Γ0 Nt i=1 q [x k,j ] i |r t l,j , φ t l,j ; a t k , (25)
where l = (k -1)N t + i and Γ 0 is the set of all possible x k,j when the k-th MTD is active 3 . The detailed derivation of the EM update rules as indicated in ( 24) and ( 25) is available in the Appendix.

4) Arithmetic Flow of the Proposed DS-AMP Algorithm:

Based on ( 12)-( 25), we summarize our proposed DS-AMP algorithm in Algorithm 1.

The details of Algorithm 1 are explained as follows. Since the active/inactive status of the MTDs is unknown, we initialize the activity indicator a 1 k = 0.5, ∀k ∈ [K], in line 1. We note that lines 4-5 are the key steps of the AMP 3 For example, if Nt = 2 and S = {+1, -1}, then

Γ 0 = +1 0 , - 1
0 , 0 +1 , 0 -1 . ( 26
)
Algorithm 1: Proposed DS-AMP Algorithm KNt ) , and the maximum iteration number T0. Output: The set of active MTDs Ω and the reconstructed media modulation signal X ∈ C KNt×J . 1: ∀i, j, k, n: We initialize the iterative index t=1, the activity indicator

Input: The received signals Y=[y1, ..., yJ ] ∈ C Nr×J , the channel matrix H = [H1, ..., HK ] ∈ C Nr×(
a 1 k = 0.5, Z 0 n,j = [yj ] n , V 0 n,j = 1, the noise variance σ 2 w 1 = 100, the reconstructed signal X = 0KN t ×J , x 1 k,j i = a 1 k s∈S s/M Nt, and v 1 k,j i = a 1 k s∈S |s| 2 /M Nt -x 1 k,j i 2 ;
2: for t = 1 to T0 do 3:

%AMP operation:

4:
∀i, j, k, n: Compute V t n,j , Z t n,j , φ t l,j , and r t l,j by using ( 20), ( 21), [START_REF] Choi | Compressed sensing for wireless communications: Useful tips and tricks[END_REF], and [START_REF] Donoho | Message-passing algorithms for compressed sensing[END_REF], respectively, where l = (k -1)Nt + i; {Decoupling step} 5:

∀i, j, k, n: Compute x t+1 k,j i and v t+1 k,j i by using ( 16) and [START_REF] Shim | Multiuser detection via compressive sensing[END_REF], respectively; {Denoising step} 6:

%EM operation:

7: ∀k: Compute (σ 2
w ) t+1 and a t+1 k by using ( 24) and ( 25); 8: end for 9: %Min-max normalization:

Ω = {k| [ã] k > 0.5}; 13: ∀k, j: η * = arg max η∈[Nt] x T 0 k,j η ; 14: ∀k ∈ Ω, ∀j: The reconstructed signal is X [(k-1)Nt+η * ,j] = x T 0 k,j η * .
operation, which consists of two parts. For the first part, a decoupling operation is performed in the same way as the original AMP algorithm in order to decouple the superimposed received signal into uncoupled scalar elements based on [START_REF] Jeong | MAP-based active user and data detection for massive machine-type communications[END_REF] and ( 18)-( 21) [START_REF] Donoho | Message-passing algorithms for compressed sensing[END_REF]- [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF]. For the second part, the denoising step computes the a posteriori mean and variance of each scalar element by using ( 16) and [START_REF] Shim | Multiuser detection via compressive sensing[END_REF], where the structured sparsity in the modulation domain and the discrete distribution of the QAM alphabet are exploited in the a priori probability. In line 7, based on ( 24) and ( 25), the EM algorithm updates the noise variance σ 2 w and the activity indicators a k , ∀k ∈ [K], by exploiting the structured sparsity in the time domain. Furthermore, the iteration stops when a predefined maximum number T 0 is reached, where T 0 is chosen to guarantee the convergence of the algorithm. Line 10 linearly transforms the estimated activity indicator vector a to ã by using the min-max normalization. In line 12, if the k-th (k ∈ [K]) element of ã is larger than 0.5, the k-th MTD is considered to be active. Line 13 selects the possible MAPs based on the structured sparsity of the media modulated symbol as indicated in (2). Finally, line 14 reconstructs the media modulated signals.

Note that all the estimated activity indicators can be small, i.e., [ a] k < 0.5, ∀k ∈ [K], in poor conditions, e.g., at low signal-to-noise ratio (SNR), where the set of active MTDs given in [START_REF] Wei | Approximate message passing-based joint user activity and data detection for NOMA[END_REF] (i.e., {k|[ a] k > 0.5}, ∀k ∈ [K]) results in the degraded active MTDs detection error. Fortunately, we find that the elements of a are linearly separable in most cases, even when all the elements are smaller than 0.5. Hence, in line 10, we preprocess a by using the min-max normalization to enlarge the separability of its elements, and the output is ã. After the preprocessing, in line 12, we can obtain the set of active MTDs based on ã. The simulation results illustrated in Section V verify the advantage of our improved active MTDs detection method over that in [START_REF] Wei | Approximate message passing-based joint user activity and data detection for NOMA[END_REF].

C. State Evolution of the DS-AMP Algorithm

SE is a tool for analyzing the performance of AMP algorithms in the large system limit, i.e., KN t → ∞, by tracking the mean-square errors (MSE) of each iteration [START_REF] Donoho | Message-passing algorithms for compressed sensing[END_REF], [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF]. In particular, capitalizing on the SE, we can characterize the performance of the proposed DS-AMP algorithm theoretically.

To start with, the MSE and average variance of the estimated signals at the t-th iteration are respectively defined as

e t = 1 KJN t K k=1 J j=1 Nt i=1 [ x t k,j ] i -[x k,j ] i 2 , ( 27
)
v t = 1 KJN t K k=1 J j=1 Nt i=1 [ v t k,j ] i . (28)
We define a scalar random variable X 0 that obeys the a priori distribution in [START_REF] Ma | Two-level sparse structure based compressive sensing detector for uplink spatial modulation with massive connectivity[END_REF]. Based on [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF], it can be shown, in large system limit and if the elements of the measurement matrix follow an i.i.d. distribution with zero mean and variance γ, that the estimated mean of x 0 at the t-th iteration, denoted as r t 0 , can be expressed as

r t 0 = x 0 + σ 2 w + γKN t e t N r γ z, (29)
where x 0 is a realization of X 0 , z ∼ CN (z; 0, 1). In addition, the estimated variance of x 0 at the t-th iteration, denoted as φ t 0 , can be expressed as

φ t 0 ≈ σ 2 w + γKN t v t N r γ . ( 30
)
By substituting ( 29) and ( 30) into ( 16) and ( 17), we obtain the approximated MMSE estimation of x 0 for the DS-AMP algorithm at the (t + 1)-th iteration. Hence, the MSE and the average variance of the estimated signals at the (t + 1)-th iteration can be expressed as

e t+1 = dx 0 p 0 (x 0 ) Dz f m (r t 0 , φ t 0 ) -x 0 2 , ( 31
)
v t+1 = dx 0 p 0 (x 0 ) Dzf v (r t 0 , φ t 0 ), (32)
where p 0 (x 0 ) is the a priori distribution in [START_REF] Ma | Two-level sparse structure based compressive sensing detector for uplink spatial modulation with massive connectivity[END_REF], Dz = e -|z| 2 /πdz, f m (r t 0 , φ t 0 ) and f v (r t 0 , φ t 0 ) are defined in ( 16) and [START_REF] Shim | Multiuser detection via compressive sensing[END_REF], respectively.

As listed in Algorithm 2, we present the detailed procedures of the SE of the proposed DS-AMP algorithm. In particular, we adopt Monte Carlo simulations to generate a large number of realizations of the transmit signals, where the sporadic traffic and the doubly structured sparsity are fully embodied. Simulation results in Section V-B will demonstrate that the SE can accurately predict the simulated results of the proposed DS-AMP algorithm. w , the sparsity level λ = Ka K , the number of MAPs Nt, the frame length J, the order of the QAM modulation, the variance γ of the elements in the measurement matrix, the number of Monte Carlo simulations NMC, the maximum SE iterations TSE, and the terminal threshold ε. Output: The theoretically predicted MSE e.

for m = 1 to NMC do 5: ∀i, j, k: r m,t l,j = x m k,j i + σ 2 w +γKNt e t Nr γ z, φ m,t l,j = σ 2 w +γKNt v t Nr γ ; 6:
∀i, j, k: x m k,j i =fm(r m,t l,j , φ m,t l,j ), v m k,j i =fv(r m,t l,j , φ m,t l,j );

7:

∀k: a t+1 k,m = fa(r m,t l,j , φ m,t l,j ; a t k,m );

8:
Calculating e t+1 m and v t+1 m referring to ( 27) and ( 28), respectively;

e t+1 = 1 N MC N MC m=1 e t+1 m , v t+1 = 1 N MC N MC m=1 v t+1 m ;
11:

e = e t+1 ;

12:

if e t+1 -e t < ε then

D. Computational Complexity of the DS-AMP Algorithm

The computational complexity of the proposed DS-AMP algorithm mainly arises from the complex-valued matrix multiplications 4 of the following operations at each iteration [START_REF] Meng | Multi-user detection for spatial modulation via structured approximate message passing[END_REF].

AMP decoupling step: The complexity of performing the AMP decoupling step, i.e., ( 18)-( 21), is O( 52 JKN t N r ). AMP denoising step: The complexity of performing the AMP denoising step, i.e., ( 15)- [START_REF] Shim | Multiuser detection via compressive sensing[END_REF], is O JKN t (|S| c + 1 4 ) . Furthermore, simulation results demonstrate that the predefined maximum number of iterations T 0 can be small to guarantee the convergence of the proposed DS-AMP algorithm. Hence the overall complexity is of the order of

O T 0 JKN t ( 5 2 N r + |S| c + 1 4
) , which scales linearly with the number of MTDs, the number of MAPs in media modulation, the order of the QAM modulation, and the number of receive antennas at the BSs. This linear complexity is appealing for efficiently processing the massive access of future IoT.

IV. PROPOSED SOLUTION FOR CODED MEDIA MODULATION BASED MMTC

In this section, we propose a solution for coded media modulation based mMTC scheme, whose block diagram is illustrated in Fig. 4. Specifically, we propose a dedicated data packet structure as well as a BICMM scheme at the MTDs, and develop an IDS-AMP detector at the receiver, so that the DADD performance can be further improved. ) can be used to determine the order for the subsequent SIC processing. By contrast, the conventional SIC order metric, e.g., the SNR or the signal-tointerference-plus-noise ratio (SINR), is difficult to be acquired in practical massive access processing [START_REF] Cho | MIMO-OFDM Wireless Communications With MAT-LAB[END_REF]. Particularly, for 0 ≤ D(b s , b s ) ≤ L s , a smaller D(b s , b s ) indicates a better decoding quality of the signature sequence, which usually implies a higher decoding quality of the associated payload data and a higher priority to be subtracted during the SIC stage. As a result, the error propagation in the SIC processing can be mitigated with the aid of the proposed signature sequence.

In addition to the data packet structure, we develop a BICMM scheme. As shown in Fig. 4, the BICMM scheme consists of an encoder, a bit-wise interleaver, and a media modulation module. Specifically, after channel coding, the length of one data packet is expanded from L bits to L ′ bits. Then, the L ′ -bit data packet is delivered to a bit-wise interleaver module. We consider a block interleaver with η columns and J = L ′ /η rows, where we assume L ′ = ηJ without loss of generality. The L ′ -bit data packet is fed into the interleaver by rows and read out by columns 5 . At the transmitter, every η bits of the interleaved L ′ -bit data packet are sequentially modulated into J media modulation symbols that are transmitted in J successive time slots (i.e., a frame). Hence, after interleaving, the bits originally associated with the same media modulated symbol are separated by (J -1) other bits, which makes the coded bits to the media modulation module more dispersed to combat the spatial-selective channel fading 6 . In particular, each media modulation-based MTD has N t unique radiation patterns associated with N t different channels with different spatial channel fading. Typically, active MTD selects one of the channels based on its input spatial data bits to transmit a media modulated symbol with η bits. Burst errors may occur if the selected channel suffers from deep fading, which results in the spatial-selective channel fading. Thanks to the bit-wise interleaver combined with the encoder, η successive bits originally associated with the same media modulated symbol are separated into different media modulated symbols with different spatial channel fading, which mitigates the spatial-selective channel fading and improves the data decoding performance [START_REF] Caire | Bit-interleaved coded modulation[END_REF], [START_REF] Koca | Bit-interleaved coded spatial modulation[END_REF], [START_REF] Goldsmith | Wireless Communications[END_REF].

B. Proposed IDS-AMP Detector at the BS

As shown in Fig. 4, the proposed IDS-AMP detector has 8 modules, including a DS-AMP algorithm module, an LLR calculation module, a deinterleaver, a soft decoder, a decoding quality judgement module, an encoder, an interleaver, and an interference cancellation module. The main procedure is summarized in Algorithm 3, and the details are illustrated as follows.

DS-AMP algorithm module (lines 4-10 of Algorithm 3): Firstly (line 4), we obtain the approximated a posteriori distribution q x k,j |y j ; σ 2 w , a k , ∀k, j, and the activity indicator vector a = [ a 1 , ..., a K ] T , which can be acquired by calling the DS-AMP algorithm (i.e., performing lines 1-8 of Algorithm 1). Secondly (lines 5-8), if the iteration index is i = 0, we acquire the indices of the active MTDs that are detected, denoted as Ω 0 , by executing lines 10-12 of Algorithm 1. The index set of the MTDs that remain to be iteratively decoded, denoted by Ω 1 , is assigned to be equivalent to Ω 0 in the first iteration (i.e., the iteration with index i = 0) and is updated in the subsequent SIC. Thirdly (lines 9, 10), we update the iteration index i = i + 1, and then select the N MTDs most likely to be active based on the activity indicators for the subsequent SIC, where the index set of the N MTDs is denoted as Ω 2 =Θ([ a] Ω1 ,N ) and we consider the MTD with a larger activity indicator is more likely to be active. If |Ω 1 | c is smaller than the predefined constant N , we set Ω 2 = Ω 1 .

LLR calculation module (line 12 of Algorithm 3):

According to the approximated a posteriori distribution q x k,j |y j ; σ 2 w , a k , ∀k, j, obtained from the DS-AMP algorithm module, we calculate the LLR information of the MTDs whose indices are in Ω 2 . For any media modulation symbol x k,j , ∀k, j, the LLR of the associated media modulated bit B MED k,j,b , ∀b ∈ [log 2 N t ], and the LLR of the associated M -QAM bit B QAM k,j,d , ∀d ∈ [log 2 M ], can be expressed, respectively, as

LLR B MED k,j,b = log x k,j ∈Φ b 0 q x k,j |y j ; σ 2 w , a k x k,j ∈Φ b 1 q (x k,j |y j ; σ 2 w , a k ) , (33)
LLR B QAM k,j,d = log x k,j ∈Ψ d 0 q x k,j |y j ; σ 2 w , a k x k,j ∈Ψ d 1 q (x k,j |y j ; σ 2 w , a k ) , (34)
where

Φ b 0 (Φ b 1
) is the set of x k,j for which the media modulated bit B MED k,j,b , ∀b, is equal to zero (one), and

Ψ d 0 (Ψ d 1 )
is the set of x k,j for which the M -QAM bit B QAM k,j,d , ∀d, is equal to zero (one) 7 .

Deinterleaver and soft decoder modules (line 14 of Algorithm 3): Firstly, the LLR information of the MTDs with indices in Ω 2 is deinterleaved. Note that the deinterleaver is an array identical to the interleaver described in Section IV-A, with the exception that the L ′ -bit data packet is read in by columns and read out by rows. Secondly, for each MTD whose index is in Ω 2 , the soft decoder decodes the deinterleaved LLR information to get the data packet of L s + L d bits. If there exist elements equal to zero in m, on the other hand, this indicates that there are one or several MTDs with almost perfect decoding quality. In this case, we apply the signals reconstruction module and the interference cancellation module (i.e., lines [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF][START_REF] Gao | Compressive-sensing based multiuser detector for the large-scale SM-MIMO uplink[END_REF][START_REF] Meng | Multi-user detection for spatial modulation via structured approximate message passing[END_REF][START_REF] Di Renzo | Spatial modulation for multipleantenna wireless systems: A survey[END_REF][START_REF] Di Renzo | Spatial modulation for generalized MIMO: Challenges, opportunities, and implementation[END_REF][START_REF] Xiao | Generalized space time block coded spatial modulation for open-loop massive MIMO downlink communication systems[END_REF][START_REF] Xiao | Differentially-encoded rectangular spatial modulation approaches the performance of its coherent counterpart[END_REF][START_REF] Renzo | Spatial modulation based on reconfigurable antennas -A new air interface for the IoT[END_REF][START_REF] Phan-Huy | Single-carrier spatial modulation for the Internet of Things: Design and performance evaluation by using real compact and reconfigurable antennas[END_REF][START_REF] Liu | Massive connectivity with massive MIMO-Part I: Device activity detection and channel estimation[END_REF]. %DS-AMP algorithm module: 4:

Φ 1 0 = +1 0 , - 1
0 , Φ 1 1 = 0 +1 , 0 -1 , (35)
Ψ 1 0 = +1 0 , 0 +1 , Ψ 1 1 = -1 0 , 0 -1 . (36)
Obtain a = [ a 1 , ..., a K ] T and q x k,j |y j ; σ 2 w , a k , ∀k, j, by performing lines 1-8 of Algorithm 1; 5:

if i = 0 then 6:

Acquire the index set of active MTDs detected, denoted as Ω 0 , by performing lines 10-12 of Algorithm 1; 7:

The index set of MTDs remaining to be iteratively decoded is denoted as Ω 1 , and

Ω 1 = Ω 0 in
Ω 1 = {Ω 1 \ Ω 3 }, Λ = {Ω 0 \ Ω 1 }; 28: H = H 0 [:,{[KNt]\ Λ}]
; {Update the measurement matrix} 29:

if Ω 1 = ∅ then 30:

break; {End the algorithm} 31:

end if 32:

end if 33: end while 34: Results: The estimated index set of active MTDs is Ω 0 and the associated demodulated payload data is

B [Ω 0 ,Ls+1:L] .
Signals reconstruction module (lines 22-24 of Algorithm 3): Firstly (line 22), we acquire the indices of the zero elements in m, denoted by

Ω 3 = {Ω 2 | n | [m] n = 0,
∀n}, which is a fine-grained index set of the MTDs with higher decoding quality. Secondly (lines 23), for each MTD with indices in Ω 3 , to reconstruct the signal components, we sequentially perform encoding, interleaving, and media modulation (i.e., we repeat the BICMM scheme at MTD) according to the de- subtract the index set Ω 3 of the MTDs cancellated in the current iteration from the index set Ω 1 of the MTDs to be decoded in the following iterations, which is updated as Ω 1 = {Ω 1 \ Ω 3 }. Also, we obtain the index set of the MTDs that are already subtracted in the previous iterations, which is denoted as Λ = {Ω 0 \ Ω 1 }. Thirdly (line 28), we update the measurement matrix as

H = H 0 [:,{[KNt]\ Λ}]
, where H 0 ∈ C Nr×(KNt) is the input channel matrix, Λ denotes the index set of the MAPs of the MTDs whose indices are in Λ, and the relationship between Λ and Λ is defined in line 1 of Algorithm 3. Finally (lines 29-31), if Ω 1 is empty, the algorithm is terminated; otherwise, the updated Y, H, and Ω 1 are fed back to line 4 to continue the next IDS-AMP iteration. Since the sparsity level in the next iteration is reduced with the unaltered dimension N r of the observations (i.e., the number of antennas at the BS), the proposed IDS-AMP detector is capable of achieving the improved decoding performance by using the aforementioned interference cancellation module.

V. CHANNEL ESTIMATION FOR MEDIA MODULATION BASED MMTC

In this section, we first present the Gauss-Markov block fading channel model. Then, the CE is discussed in two stages. In the first stage, we discuss the initial CE stage to acquire the CSI of all the MTDs based on preambles. In the following stage, we propose a data-aided CSI update strategy to reduce the training overhead.

Firstly, the Gauss-Markov block fading channel of the kth MTD in the (t + 1)-th frame (block), ∀k ∈ [K], can be expressed as

H t+1 k = √ αH t k + √ 1 -αV t k , (37)
where H t k ∈ C Nr×Nt is the MIMO channel matrix corresponding to the k-th MTD in the t-th time slot, elements in H t k obey the i.i.d. complex Gaussian distribution, and V t k ∈ C Nr×Nt is the uncorrelated channel aging error. The elements in V t k also follow the i.i.d. complex Gaussian distribution with zero mean and unit variance [START_REF] Xiang | Deep-Learning-Aided joint channel estimation and data detection for spatial modulation[END_REF]. Furthermore, α is the autoregressive (AR) coefficient and the time-domain

correlation is α τ /2 = E[[H t k ] * n,i [H t+τ k ]n,i] E[[H t k ] * n,i [H t k ]n,i] , where E [•] is the expectation operation, n ∈ [N r ], i ∈ [N t ] [44]
. The discrete time-lag τ ∈ Z is equivalent to the number of symbols within a frame (block) [START_REF] Ishikawa | Differential-detection aided large-scale generalized spatial modulation is capable of operating in high-mobility millimeter-wave channels[END_REF] and a case study will be presented in the simulations.

As for the initial CE at the first stage, the frame structure of the access signals consists a unique preamble allocated to each MTD and the following payload data, as shown in Fig. 5. The active MTDs transmit their preambles and data to the BS without any grant. At the BS, joint activity detection and channel estimation (JADCE) as well as data detection (DD) is conducted in each frame. Due to the sporadic traffic, the JADCE-DD process can be realized by using CS techniques. The details of CS-based algorithms for JADCE-DD applied to media modulation aided grant-free massive access can be found in e.g., Section VI of [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF]. In particular, the authors of [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF] employed greedy CS algorithms, while other kinds of CS algorithms can also be utilized, e.g., the AMP-based algorithms according to [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF], which is beyond the scope of this paper. Owing to the slowly time-varying CSI and the sporadic device activity [START_REF] Chen | Massive access for 5G and beyond[END_REF], the BS can obtain the CSI of all the MTDs after several frames.

At the following stage, as indicated in Fig. 5, only data symbols are transmitted in successive frames and DADD is performed at the BS. Hence, we refer to this stage as the DADD stage. Note that the CSI utilized in the first frame of the DADD stage is obtained at the JADCE-DD stage. Furthermore, to combat the channel aging and to reduce the training overhead, we propose a CSI update strategy with the aid of the detected data symbols, as shown in Fig. 5. Specifically, by providing

H t = [H t 1 , H t 2 , ..., H t K ] ∈ C Nr×(KNt)
Y t ≈ H t X t + W t , (38)
where X t = X t [ Ω t ,:] ∈ C NaNt×J and H t = H t [:, Ω t ] ∈ C Nr×NaNt are the estimated media modulated signals and the associated CSI for active MTDs, respectively, and the approximately equal sign "≈" can be written as an equal sign "=" if perfect DADD is achieved. Due to the slowly time-varying IoT channels, the frame length J can be longer than N a N t , where N a ≤ K a and K a is only a fraction of the total MTDs. Hence, we can refine the CSI H t based on X t . For example, by using the MMSE estimator 8 , we can obtain the refined CSI of active MTDs, denoted as

H t = Y t ( X t ) H R t H X t + N a N t σ 2 w I -1 ( X t ) H R t H ∈ C Nr×NaNt , where R t H = E ( H t
) H H t is the channel's covariance matrix [START_REF] Sugiura | Effects of channel estimation on spatial modulation[END_REF]. Then, in the (t + 1)-th frame, we set H t+1 = H t ∈ C Nr×(KNt) and update part of it by performing H t+1 [:, Ω t ] = H t . Hence, we can adopt the updated CSI H t+1 as an input of the proposed Algorithm 1 or Algorithm 3 for DADD in the (t + 1)-th frame.

VI. PERFORMANCE EVALUATION A. Simulation Setup

In this section, an extensive simulation investigation is carried out to evaluate the activity detection error rate (ADER) of the MTDs, the symbol error rate (SER), and the bit error rate (BER) of the proposed massive access solution. The ADER, SER, and BER are defined as

ADER = E m + E f K , (39)
SER = JE m + E symbol JK a , (40)
BER = ηJE m + E MED + E QAM ηJK a , (41)
respectively, where E m is the number of active MTDs missed to be detected, E f is the number of inactive MTDs falsely detected to be active, E symbol is the number of error symbols 9 of the detected active MTDs, JK a is the total number of symbols transmitted by the K a active MTDs within one frame, E MED and E QAM are the numbers of the errors of the media modulated bits and the quadrature amplitude modulated bits of the detected active MTDs within one frame, respectively, and ηJK a is the number of total bits transmitted by the K a active MTDs within one frame. As for the CSI update strategy, given the t-th frame, ∀t, NMSE H denotes the normalized mean-square errors (NMSE) between the updated CSI H t ∈ C Nr×(KNt) utilized for DADD and the true CSI H t ∈ C Nr×(KNt) . In addition, the NMSE between the reconstructed signal matrix X t and the true signal matrix X t is denoted as NMSE X . Specifically, NMSE H and NMSE X can be respectively expressed as

NMSE H = H t -H t F H t F , (42)
NMSE X = X t -X t F X t F . (43)
8 Denote Pu as the probability when the estimated media modulated signal matrix X t

[ Ω t u ,:] ∈ C Nt×J of given MTD, ∀u ∈ [Na]
, is a singular matrix.

Then, considering the typical case, i.e., one row of X t

[ Ω t u ,:]
is all zeros, we have Pu ≈ JNt(1 -1 Nt ) J . If J = 200 and Nt = 4, we have Pu ≈ 8.2 × 10 -23 . Hence, it is suitable to assume that X t is a non-singular matrix, which enables us to adopt the MMSE estimator. 9 We consider that one symbol is in error if either its media modulated bits or quadrature amplitude modulated bits are in error.

In the simulations, the number of MTDs is K = 500 with K a = 50 active MTDs, where each MTD adopts N RF = 2 RF mirrors for media modulation and 4-QAM (M = 4). Hence, the throughput is η = N RF + log 2 M = 4 bpcu. Moreover, the number of receive antennas is N r = 256, the maximum iteration number is set to T 0 = 15. Additionally, the frame length J is set to 12 for uncoded media modulation based mMTC. For the SE of the DS-AMP algorithm, the number of Monte Carlo simulations is N MC = 500, the maximum number of iterations is T SE = 50, and the terminal threshold is ε = 10 -5 . Since we can obtain the a posteriori estimation of the media modulation signals x k,j , ∀k, j, in each Monte Carlo simulation (i.e., line 6 in Algorithm 2), the ADER, BER, and SER of the theoretical SE can be calculated in the same way as those for the DS-AMP algorithm, and they can be then averaged over all Monte Carlo simulations. For the proposed coded media modulation based mMTC scheme, we consider a Turbo code with 1/3 rate and 12 tail bits. The length of the data packet is L = 120 with the length of the signature sequence being L s = 20. Hence, after channel encoding, the length of the data packet is L ′ = 3L+12 = 372 and the frame length is J = L ′ /η = 93 for coded media modulation based mMTC. In addition, N is set to 5 in line 10 of Algorithm 3. Without loss of generality, we investigate the ADER, SER, and BER performance of different algorithms in any given frame with perfect CSI at the BS, and a Rayleigh MIMO channel is considered.

B. Performance of the Proposed DS-AMP Algorithm

For comparison, we consider the following benchmarks. Benchmark 1: LMMSE multi-user detector for a traditional uplink mMIMO system [START_REF] Gao | Compressive-sensing based multiuser detector for the large-scale SM-MIMO uplink[END_REF], where K a single-antenna users (after the grant-based scheduling) adopting 16-QAM (for achieving the same throughput of 4 bpcu) are supported by an mMIMO BS with N r = 256 receive antennas. Benchmark 2: The StrOMP algorithm (i.e., algorithm 1 in [START_REF] Qiao | Compressive sensing based massive access for IoT relying on media modulation aided machine type communications[END_REF]) is used for the activity detection and the SIC-SSP algorithm (i.e., algorithm 2 in [START_REF] Qiao | Compressive sensing based massive access for IoT relying on media modulation aided machine type communications[END_REF]) is used for data detection, where the terminal threshold P th for the StrOMP algorithm is set to 1.5. Benchmark 3: A modified DS-AMP algorithm without executing the min-max normalization (i.e, the min-max normalization in line 10 of Algorithm 1 is replaced by ã = a), where the activity detection method is the same as that in [START_REF] Wei | Approximate message passing-based joint user activity and data detection for NOMA[END_REF]

(i.e., {k|[ a] k > 0.5}, ∀k ∈ [K]). AMP:
The conventional AMP algorithm [START_REF] Meng | An expectation propagation perspective on approximate message passing[END_REF] (i.e., only perform lines 1-5 and 8 in Algorithm 1), where the sparsity level is λ = Ka K and the noise variance σ 2 w are perfectly known in advance, and the a priori probability in [START_REF] Ma | Two-level sparse structure based compressive sensing detector for uplink spatial modulation with massive connectivity[END_REF]

is replaced by p [x k,j ] i = (1 -λ)δ [x k,j ] i + λ M s∈S δ [x k,j ] i -s . TLSSCS:
The cutting-edge TLSSCS algorithm from [START_REF] Ma | Two-level sparse structure based compressive sensing detector for uplink spatial modulation with massive connectivity[END_REF], where the scaling factor α = 4 (i.e., α in (6) of literature [START_REF] Ma | Two-level sparse structure based compressive sensing detector for uplink spatial modulation with massive connectivity[END_REF]). PIA-MSMP: The state-of-the-art PIA-MSMP algorithm (i.e., algorithm 1 in [START_REF] Ma | Improved compressed sensing-based joint user and symbol detection for media-based modulation-enabled massive machine-type communications[END_REF]) with the perfectly known sparsity level. scheme with the proposed DS-AMP algorithm outperforms the traditional uplink mMIMO scheme (benchmark 1) when the same throughput is considered. This verifies the advantages of the proposed media modulation based mMTC scheme over the traditional scheme. Since benchmark 1 perfectly knows the indices of the active MTDs, comparing the proposed scheme with benchmark 1 is actually unfair. Secondly, it can be observed that the proposed DS-AMP algorithm outperforms the TLSSCS algorithm, the PIA-MSMP algorithm, and benchmark 2 in terms of ADER, SER, and BER performance, which verifies the superiority of the proposed algorithm. Moreover, thanks to the exploitation of the doubly structured sparsity, our proposed DS-AMP algorithm outperforms the conventional AMP algorithm in ADER, SER, and BER. Furthermore, it is apparent that the proposed DS-AMP algorithm outperforms benchmark 3 in the low SNR regime (i.e., 0 dB∼2 dB), which verifies the effectiveness of the proposed min-max normalization (i.e., line 10 in Algorithm 1) in the DS-AMP algorithm. Finally, we observe that the SE offers a good tightness compared with the proposed DS-AMP algorithm in terms of ADER, SER, and BER performance, which can be observed from Fig. 7, Fig. 8, and Fig. 9 as well.

From Fig. 7, we observe that the proposed DS-AMP algorithm outperforms the conventional AMP algorithm, the TLSSCS algorithm, the PIA-MSMP algorithm, and benchmark 2 versus different frame lengths at SNR = 5 dB in terms of ADER, SER, and BER performance. From Fig. 7(a), owing to the exploitation of the structured sparsity in the time domain, it can be observed that the advantage of the proposed DS-AMP algorithm over other algorithms in terms of ADER performance becomes more evident upon increasing J. From Fig. 7(b) and Fig. 7(c), if J is small (i.e., J < 5), we can obtain improved BER and SER performance as well as improved ADER performance. If J is large (i.e., J > 9), on the other hand, the BER and SER performance almost stays unaltered for different values of J. Since the conventional AMP algorithm reconstructs the signals of each time slot separately, its performance remains unchanged against different values of J. under different system parameters, i.e., the sparsity level or the number of receive antennas, in typical IoT scenarios. Fig. 10 depicts the ADER and BER performance of the proposed DS-AMP algorithm versus the maximum iteration number T 0 . We observe that the ADER and BER performance of the proposed DS-AMP algorithm converges fast at various SNRs (usually fewer than 15 iterations are needed), which indicates that we can adopt the maximum iteration number T 0 = 15 for Algorithm 1. In particular, the SER performance of the DS-AMP algorithm versus T 0 is similar to that of the BER and ADER performance.

The computational complexity comparison of different solutions is provided in Table I. It is apparent that the computational complexity of the proposed DS-AMP algorithm is an order of magnitude lower than that of the TLSSCS algorithm, the PIA-MSMP algorithm, and benchmark 2 under the considered simulation parameters. Moreover, the complexity of the DS-AMP algorithm scales linearly with the number of receive antennas N r , whereas the computational complexity of the TLSSCS and PIA-MSMP algorithms can be approximately proportional to the square of N r . Hence, the proposed DS-AMP is more attractive than other state-of-the-art algorithms for solving the massive access problem in mMIMO systems. Furthermore, compared with the conventional AMP algorithm and benchmark 3, the proposed DS-AMP algorithm achieves better performance without substantially increasing the computational complexity.

C. Performance of the Proposed IDS-AMP Scheme

In this subsection, we compare the SER and BER performance of the proposed IDS-AMP scheme with the following benchmarks. Benchmark 4: The proposed DS-AMP algorithm adopting uncoded media modulation and a hard decision (i.e., it performs a hard decision according to the output signal X ∈ C KNt×J from Algorithm 1 to get the demodulated binary bits). Benchmark 5: The proposed DS-AMP algorithm adopting coded media modulation and soft decision, while the processing of interleaving/deinterleaving and SIC in Fig. 4 is not adopted. Benchmark 6: The proposed DS-AMP algorithm adopting coded media modulation, interleaving/deinterleaving, and soft decision, while the SIC processing is not adopted. Benchmark 7: The proposed IDS-AMP scheme with the exception of the proposed decoding quality judgement (i.e., lines 15-17 of Algorithm 3) is removed and Ω 3 in line 22 is equal to Ω 2 . Fig. 11 compares the SER and BER performance of benchmarks 4-7 and the proposed IDS-AMP scheme. The worst performance is achieved by benchmark 4, which indicates the necessity of adopting channel coding and soft decoding for improving the data decoding performance. The superiority of benchmark 6 over benchmark 5 verifies the efficiency of 1 The order of complex-valued multiplications is obtained under parameters J = 12, Nt = 4, K = 500, Ka = 50, T0 = 15, |S|c = 4. the proposed BICMM in overcoming the effect of burst error in spatial-selective channel fading. Besides, the superiority of benchmark 7 over benchmark 6 in the high SNR regime (i.e., greater than -2 dB) verifies the effectiveness of SIC processing in improving the data decoding performance. However, benchmark 7 is observed to be inferior to benchmark 5 in the low SNR regime (i.e., -3.5 dB∼-2.5 dB), since the SIC at a low SNR can degrade the performance. By contrast, the proposed IDS-AMP scheme shows a consistent superiority over the four considered benchmarks. Particularly, the superiority of the IDS-AMP scheme over benchmark 7 verifies the data decoding gain achieved by the proposed decoding quality judgement module.

D. Performance of the CSI Update Strategy

To investigate the data-aided CSI update strategy, we consider successive N f (N f ≫ 1) frames (blocks), where the Gauss-Markov block fading channel model ( 37) is considered and a Rayleigh MIMO channel model is utilized for the first frame. The coherence time T c is defined by the duration that the time-domain correlation function is above 0.5 [START_REF] Ishikawa | Differential-detection aided large-scale generalized spatial modulation is capable of operating in high-mobility millimeter-wave channels[END_REF] and we have T c = 0.423/f m , where f m = vf c /v c is the maximum Doppler shift, v is the maximum velocity of the MTDs, and v c is the speed of light [START_REF] Ishikawa | Differential-detection aided large-scale generalized spatial modulation is capable of operating in high-mobility millimeter-wave channels[END_REF]. According to 3GPP Narrowband Internet-of-Things (NB-IoT) specifications, W = 15 kHz bandwidth for single-tone uplink transmission with carrier frequency f c = 900 MHz can be supported [START_REF] Ratasuk | NB-IoT system for M2M communication[END_REF]. In this case, the time-lag τ = ⌊T c /T s ⌋ = 213, if v = 35.6 km/h, where T s = 1/W denotes the symbol duration. Since the timedomain correlation α τ /2 has to be above 0.5 between the time indices i = 1 and i = τ + 1, we obtain the AR coefficient α = 0.9935 when α τ /2 is set to 0.5 [START_REF] Ishikawa | Differential-detection aided large-scale generalized spatial modulation is capable of operating in high-mobility millimeter-wave channels[END_REF]. Furthermore, the frame length after channel coding is J = 213, and the associated number of bits transmitted by a frame is 280 bits. 12. The CSI update process of the t-th frame, 0 < t < N f . The AR coefficient is 0.99. The encoder and other system parameters are the same as that in Section VI-A, and the proposed DS-AMP algorithm with channel coding is used for DADD in each frame.

For comparison, we consider two strategies, both of which assume that the CSI of all the MTDs has been acquired at the BS in the first frame. In particular, for the CSI nonupdate strategy: the CSI used for DADD in each frame is the same as that in the first frame. As for the proposed CSI update strategy: the CSI of correctly detected active MTDs are refined by the MMSE estimator at the BS in the t-th frame, 0 < t < N f , as shown in [START_REF] Basar | Index modulation techniques for next-generation wireless networks[END_REF], and then the updated CSI matrix of all the MTDs is used for DADD in the (t+1) frame. In particular, by decoding, re-encoding, and media modulating the bits associated with the detected active MTDs, as shown in Fig. 12, we can obtain the signal matrix X t for refining the CSI matrix H t . Then, H t+1 is updated as shown in [START_REF] Basar | Index modulation techniques for next-generation wireless networks[END_REF].

It can be seen from Fig. 13, as the time increases, the NMSE H and NMSE X performance of the CSI non-update strategy decreases monotonically due to the channel aging. On the contrary, by exploiting the proposed CSI update strategy, the NMSE H and NMSE X performance dramatically outperforms that of the non-update strategy, especially, after a relatively large number of frames. Hence, by using the proposed CSI update strategy in slowly time-varying IoT channel scenarios, it is unnecessary for the BS to estimate the CSI of all the MTDs in every frame, which can reduce the training overhead significantly.

VII. CONCLUSIONS

This paper investigated the DADD problem for media modulation based mMTC, where both uncoded and coded transmissions were considered. By exploiting the doubly structured sparsity of media modulation signals, we first proposed a DS-AMP algorithm to solve the DADD problem for the uncoded case. Also, we derived the SE of the DS-AMP algorithm to theoretically predict its performance. Furthermore, for the coded case, we developed a BICMM scheme and proposed an IDS-AMP scheme based on SIC for improving the data decoding performance. In particular, a signature sequence part for facilitating the SIC processing was proposed to be embedded into the dedicated data packet and then the IDS-AMP detector was developed, whereby the estimated signal components were iteratively subtracted from the received signals for enhancing the data decoding performance. In addition, we discussed the CE problem and developed a data-aided CSI update strategy to reduce the training overhead in block fading channels. Finally, extensive simulation results demonstrated the superiority of the proposed DS-AMP algorithm over cutting-edge algorithms in terms of ADER, SER, and BER performance while ensuring a lower computational complexity in the uncoded transmission. Improved data decoding performance of the proposed SICbased IDS-AMP scheme was verified in the coded case as well. Also, simulation results verified the efficiency of the proposed CSI update strategy.

APPENDIX

Here we provide the detailed derivation of the EM update rules in the proposed algorithms. According to ( 22) and ( 23

In order to find the parameter set θ that maximizes Q (θ, θ t ) in the (t + 1)-th iteration, we differentiate Q (θ, θ t ) with respect to θ and let it equal zero.

First, according to ( 9) and ( 44

where Γ = {Γ 0 , 0 Nt×1 } and Γ 0 is defined in [START_REF] Di Renzo | Spatial modulation for generalized MIMO: Challenges, opportunities, and implementation[END_REF]. The function q [x k,j ] i |r t k,j , φ t k,j ; σ 2 w t , a t k can be calculated according to [START_REF] Wang | Dynamic compressive sensing-based multi-user detection for uplink grant-free NOMA[END_REF]

and C is the normalization factor, V t n,j and Z t n,j are calculated in [START_REF] Meng | An expectation propagation perspective on approximate message passing[END_REF] and [START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF], respectively.

Hence, according to (49), (50), (52), and (53), we obtain the EM update rule of the noise variance as (54)

Fig. 1 .Fig. 2 .

 12 Fig. 1. Media modulation based mMTC scheme: The MTDs adopt media modulation to access the mMIMO BSs.

2 s

 2 .t. (2), a 0 ≪ K, and s k,j ∈ S, k ∈ [K], j ∈ [J]. (4) B. Proposed DS-AMP Algorithm for DADD 1) Problem Formulation Based on Factor Graph:

Fig. 3 .

 3 Fig. 3. Factor graph of the joint posterior distribution p xj |y j ; σ 2 w , a , where the activity of the k-th MTD, ∀k ∈ [K], remains unchanged in J successive time slots and is denoted by a k . The circles represent the variable nodes and the squares represent the factor nodes.

10 :

 10 Let ã = a-min( a) max( a)-min( a) , where a=[ a1, ..., aK] T = a T 0 1 , ..., a T 0 K T , min(•) and max(•) are the minimum value and maximum value of the arguments, respectively; 11: %Extract the active MTDs and their MAPs: 12: ∀k: The set of active MTDs

Algorithm 2 :

 2 State Evolution of DS-AMP AlgorithmInput: The noise variance σ 2

1 :

 1 ∀m ∈ [NMC]: Generate NMC realizations of the transmit signals X m ∈ C KNt×J , according to the a priori distribution in (9). 2: ∀m, k: Define e 1 = 0N MC ×1 and v 1 = 0N MC ×1 to record the predicted MSE and average variance of the m-th Monte Carlo realization. We initialize the iteration number t = 1, the predicted MSE e 1 = 1, the average variance v 1 = 1, and the activity indicators for the m-th signal realization a 1 k,m = 0.5; 3: for t = 1 to TSE do 4:

Fig. 4 .

 4 Fig. 4. Communication process of the proposed massive access solution for coded media modulation based mMTC.A. Dedicated Data Packet Structure and BICMM at the MTDsThe proposed data packet structure is tailored for the implementation of an SIC-based IDS-AMP detector at the receiver. Specifically, the data packet is composed of two parts, i.e., the L s -bit signature sequence part and the L dbit payload data. We consider that all the MTDs share the same binary signature sequence b s ∈ N Ls×1 , which is a predefined pseudo-random 0/1 sequence known at the transceiver. At the receiver, we consider the Hamming distance D(b s , b s ) between b s and b s as a metric to evaluate the decoding quality of the associated payload data part of the MTDs, where b s ∈ N Ls×1 is the estimated binary signature sequence of any detected MTD. Hence, D(b s , b s ) can be used to determine the order for the subsequent SIC processing. By contrast, the conventional SIC order metric, e.g., the SNR or the signal-tointerference-plus-noise ratio (SINR), is difficult to be acquired in practical massive access processing[START_REF] Cho | MIMO-OFDM Wireless Communications With MAT-LAB[END_REF]. Particularly, for 0 ≤ D(b s , b s ) ≤ L s , a smaller D(b s , b s ) indicates a better decoding quality of the signature sequence, which usually implies a higher decoding quality of the associated payload data and a higher priority to be subtracted during the SIC stage. As a result, the error propagation in the SIC processing can be mitigated with the aid of the proposed signature sequence.In addition to the data packet structure, we develop a BICMM scheme. As shown in Fig.4, the BICMM scheme consists of an encoder, a bit-wise interleaver, and a media modulation module. Specifically, after channel coding, the length of one data packet is expanded from L bits to L ′ bits. Then, the L ′ -bit data packet is delivered to a bit-wise interleaver module. We consider a block interleaver with η columns and J = L ′ /η rows, where we assume L ′ = ηJ without loss of generality. The L ′ -bit data packet is fed into the interleaver by rows and read out by columns5 . At the transmitter, every η bits of the interleaved L ′ -bit data packet are sequentially modulated into J media modulation symbols that are transmitted in J successive time slots (i.e., a frame). Hence, after interleaving, the bits originally associated with the same media modulated symbol are separated by (J -1) other bits, which makes the coded bits to the media modulation module more dispersed to combat the spatial-selective channel

  Finally, for the |Ω 2 | c MTDs detected as active, the decoded L|Ω 2 | c bits are recorded in a binary data packet matrix B ∈ N K×L , i.e., B [Ω2,1:L] . In particular, the first L s columns of B store the decoded signature sequence bits and the remaining L d columns record the decoded payload data bits. Decoding quality judgement module (lines 15-17 of Algorithm 3): Firstly (lines 16), for each MTD whose index is in Ω 2 | n , ∀n ∈ [|Ω 2 | c ], we calculate the Hamming distances between the decoded signature sequence B [Ω2|n,1:Ls] T and the true signature sequence b s , denoted as [m] n = D b s , B [Ω2|n,1:Ls] T . Secondly (line 17), based on the Hamming distances m ∈ N |Ω2|c×1 , whether and how to perform the SIC processing is judged. In particular, if no element of m is equal to zero, we skip the interference cancellation module, decode the remaining MTDs indexed in {Ω 1 \ Ω 2 } in line 18, and stop the algorithm in line 19.

Algorithm 3 :

 3 Proposed IDS-AMP Detector Input: The received signals Y = [y 1 , ..., y J ] ∈ C Nr×J and the channel matrix H = H 0 =[H 1 , ..., H K ] ∈ C Nr×(KNt) . Output: The estimated index set of active MTDs Ω 0 , and the associated demodulated payload data bits B [Ω 0 ,Ls+1:L] with B ∈ N K×L defined as the binary data packet matrix of all MTDs. 1: Initialization: The iteration index i=0, the binary data packet matrix B = 0 K×L , and for possible active MTDs given their temporary index set Λ, their MAPs' index set is defined as Λ = { Λn} |Λ|c n=1 , where Λn = {Nt(Λ|n-1)+u} Nt u=1 is the MAPs' index set of the n-th MTD in Λ for n ∈ [|Λ|c]; 2: while 1 do 3:

Fig. 5 .

 5 Fig. 5. Schematic diagram of the two-stage massive access scheme: the JADCE-DD stage and the following DADD stage, where the frame structure at the transmitter and the processing procedure at the receiver for each stage are presented in detail. coded B [Ω3,1:L] . Thirdly (line 24), we record the reconstructed signals in X[ Ω3,:] , where Ω 3 is the index set of the MAPs of the MTDs whose indices are in Ω 3 . Interference cancellation module (lines 26-31 of Algorithm 3): Firstly (line 26), we subtract the signal components X[ Ω3,:] from the received signals Y. Secondly (line 27), we

  and the received signal Y t = [y t 1 , ..., y t J ] ∈ C Nr×J as the input of the proposed Algorithm 1 or Algorithm 3, we can obtain the set of correctly detected active MTDs Ω t and the reconstructed signal X t ∈ C KNt×J of the t-th frame, ∀t. The number and the set of MAPs for the correctly estimated active MTDs are denoted as N a = |Ω t | c (N a ≤ K a ) and Ω t = { Ω t u } Na u=1 , respectively, where Ω t u = {N t (Ω t | u -1) + i} Nt i=1 is the set of MAPs for the u-th MTD in Ω t for u ∈ [N a ]. Hence, we can transform the system model in (3) as

Fig. 6 (Fig. 6 .

 66 a), Fig. 6(b), and Fig. 6(c) compare the ADER, SER, and BER performance of different solutions versus the SNR, respectively. Firstly, we observe that the SER and BER performance of the media modulation based mMTC Performance comparison of different solutions versus SNR: (a) ADER performance comparison; (b) SER performance comparison; (c) BER performance comparison.

Fig. 7 .

 7 Performance comparison of different solutions versus the numbers of time slots J within a frame at SNR = 5 dB: (a) ADER performance comparison; (b) SER performance comparison; (c) BER performance comparison.

Fig. 8 Fig. 8 .Fig. 9 .

 889 provides the ADER, SER, and BER performance comparisons of different algorithms versus the sparsity level (λ = Ka K given K = 500 and SNR = 3 dB). From Fig. 8(a), Fig. 8(b), and Fig. 8(c), we observe that the proposed DS-AMP algorithm outperforms the conventional AMP algorithm, the TLSSCS algorithm, the PIA-MSMP algorithm, and benchmark 2 in terms of ADER, SER, and BER performance, respectively. Fig. 9 depicts the ADER, SER, and BER performance of different algorithms versus the numbers of receive antennas N r at SNR = 5 dB. From Fig. 9, similar conclusion as observed in Fig. 8 can be obtained. In particular, both figures verify the superiority and robustness of the proposed DS-AMP algorithm Performance comparison of different solutions versus the sparsity level λ = Ka K , given K = 500 and SNR = 3 dB: (a) ADER performance comparison; (b) SER performance comparison; (c) BER performance comparison. Performance comparison of different solutions versus the numbers of receive antennas Nr at SNR = 5 dB: (a) ADER performance comparison; (b) SER performance comparison; (c) BER performance comparison.

Fig. 10 .

 10 Fig. 10. Performance of the proposed DS-AMP algorithm versus the maximum iteration number T 0 .

Fig. 11 .

 11 Fig. 11. Performance of the proposed SIC-based massive access solution in comparison with the benchmarks: (a) SER performance comparison; (b) BER performance comparison.

Fig. 13 .

 13 Fig. 13. Performance comparison of the proposed CSI update strategy versus time (i.e., the number of frames) at SNR = 30 dB and the AR coefficient α = 0.99: (a) NMSE H performance comparison; (b) NMSE X performance comparison.

  ), ln p (X, Y; θ) can be expressed as ln p (X, Y; θ) = ln p Y|X; σ 2 w + ln p (X; a) = J j=1 ln p y j |x j ; σ 2 w + ln p (x j ; a) y j ] n |x j ; σ 2 w + K k=1 ln p (x k,j ; a k ) .

  Ω denote the l p norm, m-th element of x, and the elements of x indexed by Ω, respectively. Θ(x, n) is the operator that selects the indices of the top n largest elements of x. For an ordered set Γ and its subset Ω, |Γ| c , Γ| i , and Γ \ Ω denote the cardinality of Γ, the i-th element of Γ, and the complement of subset Ω in Γ, respectively. For the binary codewords c 1 and c 2 , we define their Hamming distance D(c 1 , c 2 ) as the number of elements in which they differ. [K] is the set {1, 2, ..., K} and 0 m×n denotes a zero matrix with size m × n. If h is the joint distribution of variable x = [x 1 , x 2 , x 3 ] T , then the marginal distribution of x 2 is denoted by ∼{x2} h(x) = x1∈A x3∈A h(x), where x i , ∀i ∈ {1, 2, 3}, belongs to a finite domain A. Finally, CN (x; µ, ν) denotes the complex Gaussian distribution of a random variable x with mean µ and variance ν.

T , (•) * , and (•) -1 represent transpose, conjugate, and inverse, respectively. For a matrix A and an ordered set Ω, A [m,n] denote the m-th row and n-th column element of A, A [Ω,:] (A [:,Ω] ) is the sub-matrix containing the rows (columns) of A indexed by Ω, and A [Ω,m:n] is the sub-matrix containing from the m-th to the n-th columns of A [Ω,:] . For a vector x and an ordered set Ω, x p , [x] m , and [x]

  x2 , ..., xJ ] ∈ C (KNt)×J , and W = [w 1 , w 2 , ..., w J ] ∈ C Nr×J . Hence the massive access problem can be formulated as the following optimization problem

			J
	min X	Y -HX 2 F = min {xj} J j=1	j=1	y j -Hx j	2 2

Table I :

 I Computational complexity comparison of different algorithms for uncoded media modulation based mMTC TLSSCS O{(JNrKa + 2NrKa 2 + Ka 3 ) + (Ka + 1)[Nr 2 (KNt + J) + NrJKNt] +

	Algorithms	Computational complexity	Complex-valued multiplications 1 Nr = 128 Nr = 256
	Benchmark 1	O(JNrKa + 2NrKa 2 + Ka 3 )	0.84 × 10 6	1.56 × 10 6
	DS-AMP AMP Benchmark 3 PIA-MSMP	O[T0JKNt( 5 2 Nr + |S|c + 1 4 )] O[T0JKNt( 5 2 Nr + |S|c + 1 4 )] O[T0JKNt( 5 2 Nr + |S|c + 1 4 )] Ka+1 s=1 [Nr 2 + 2Nr(sNt) 2 + (sNt) 3 ]} O{3JKaNr(Nt + 1) + (Ka + 1)[Nr 2 (KNt + J) + NrJKNt] + Ka 2Nr(sNt) 2 + (sNt) 3 ]} s=1 [Nr 2 +	1.17 × 10 8 1.17 × 10 8 1.17 × 10 8 2.14 × 10 9 2.12 × 10 9	2.32 × 10 8 2.32 × 10 8 2.32 × 10 8 7.53 × 10 9 7.50 × 10 9
	Benchmark 2	O{KaJKNtNr + Ka s=1 [JNr(s + 2s 2 + 2(sNt) 2 ) + J(s 3 + (sNt) 3 )]} Ka s=1 [JNr(s + 2s 2 + 2(sNt) 2 ) + J(s 3 + (sNt) 3 )] +	4.82 × 10 9	8.16 × 10 9

  ), we differentiate Q (θ, θ t ) with respect to a k , ∀k, as follows∂Q (θ, θ t ) ∂a k = ∂E {ln p (X, Y; θ) |Y; θ t } ∂a k

	=	J j=1	E	d da k	lnp (x k,j ; a k ) |Y; θ t
	=	J j=1		x k,j ∈Γ	p x k,j |Y; θ t d da k	lnp (x k,j ; a k )
	≈	J j=1 x k,j ∈Γ	q x k,j |r t k,j , φ t k,j ; σ 2 w	t , a t k	d da k	lnp (x k,j ; a k )
	= j=1 x k,j ∈Γ J i=1 Nt q [x k,j ] i |r t k,j ,φ t k,j ; σ 2 w	t ,a t k	d da

k lnp (x k,j ; a k ) ,

  Hence, the EM update rule for a k , ∀k, can be expressed as Next, definez j = Hx j (∀j ∈ [J]), then we have [z j ] n = K k=1 [H k x k,j ] n (n ∈ N r ). According to (8) and (44), we differentiate Q (θ, θ t ) with respect to σ 2 [y j ] n -[z j ] n 2 |Y; θ t , (49) where for n ∈ [N r ], E [y j ] n -[z j ] n 2 |Y; θ t = [y j ] n -E [z j ] n |y j ; θ t 2 + Var [z j ] n |y j ; θ t ,(50)and Var {•|y j ; θ t } is the variance conditioned on y j with parameter θ t .Note that the distribution of [z j ] n can be expressed asf [z j ] n |y j ; θ t = 1 C f [y j ] n | [z j ] n ; θ t f [z j ] n ; θ t = 1 C CN [y j ] n ; [z j ] n , σ 2 w t CN [z j ] n ; Z t n,j , V t n,j = CN [z j ] n ; E [z j ] n |y j ; θ t , Var [z j ] n |y j ; θ t ,(51)whereE [z j ] n |y j ; θ t = σ 2Var [z j ] n |y j ; θ t = σ 2

	Hence, we have
	σ 2 w	t+1 =	1 JN r	J j=1	Nr n=1	E w	t Z t n,j + [y j ] n V t n,j (σ 2 w ) n,j t + V t	,	(52)
					w w ) (σ 2 t + V t t V t n,j n,j	,
										and
										d da k	lnp (x k,j ; a k ) =	-	Nt i=1 δ [x k,j ] i p (x k,j ; a k )	+
										1 N t	Nt i=1	1 M s∈S	δ [x k,j ] (46)
										a t+1 k	=	1 J	J j=1 x k,j ∈Γ0	Nt i=1	q [x k,j ] i |r t l,j , φ t l,j ; a t k .	(47)
										w as follows
										∂Q (θ, θ t ) ∂ (σ 2 w )	=	∂E {ln p (X, Y; θ) |Y; θ t } ∂ (σ 2 w )
										=	J j=1	Nr n=1	E	d d (σ 2 w )	lnp [y j ] n |x j ; σ 2 w |Y; θ t
										=	J j=1	Nr n=1	E	d d (σ 2 w )	-ln σ 2 w -	[y (48)

is n∈[Nt],n =i δ [x k,j ] n p (x k,j ; a k ) . j ] n -[z j ] n 2 (σ 2 w ) |Y; θ t .

Although developed in the large system regime, the AMP algorithm performs well even for medium size problems, such as massive access with hundreds or thousands of MTDs[START_REF] Donoho | Message-passing algorithms for compressed sensing[END_REF],[START_REF] Ke | Compressive sensing based adaptive active user detection and channel estimation: Massive access meets massive MIMO[END_REF].

For the sake of simplicity, a real-valued multiplication is assumed to have a complexity that equals a quarter of a complex-valued multiplication.

Note that the deinterleaver is an array identical to the interleaver, while the L ′ bits data packet is read in by columns and read out by rows.

Practical slowly time-varying fading channels may exhibit burst error which is beyond the error correction capability of the codes of reasonable complexity[START_REF] Goldsmith | Wireless Communications[END_REF]. Therefore, in fading channels, coding is typically combined with interleaving to mitigate the effect of burst error[START_REF] Goldsmith | Wireless Communications[END_REF].

For example, supposing that Nt = 2, M = 2, S = {+1, -1}, b ∈ [1], and d ∈ [1], then we can get

1]. The work was supported by NSFC under Grants 62071044 and 61827901, the BJNSF under Grant L182024. D. W. K. Ng is supported by funding from the UNSW Digital Grid Futures Institute, UNSW, Sydney, under a crossdisciplinary fund scheme and by the Australian Research Council's Discovery Project (DP210102169). The codes and some other materials about this work may be available at https://gaozhen16.github.io.