Compressive Sensing Based Adaptive Active User Detection and Channel Estimation: Massive Access Meets Massive MIMO

TL;DR

This paper introduces an adaptive compressive sensing framework for active user detection and channel estimation in massive MIMO systems, significantly reducing latency and improving accuracy by exploiting channel sparsity.

Contribution

It develops a novel GMMV-AMP algorithm and an adaptive scheme that jointly detect users and estimate channels efficiently in massive MIMO systems.

Findings

The proposed algorithms outperform baseline schemes in accuracy.

The adaptive scheme guarantees reliable access with unknown sparsity levels.

State evolution accurately predicts the algorithm's performance.

Abstract

This paper considers massive access in massive multiple-input multiple-output (MIMO) systems and proposes an adaptive active user detection and channel estimation scheme based on compressive sensing. By exploiting the sporadic traffic of massive connected user equipments and the virtual angular domain sparsity of massive MIMO channels, the proposed scheme can support massive access with dramatically reduced access latency. Specifically, we design non-orthogonal pseudo-random pilots for uplink broadband massive access, and formulate the active user detection and channel estimation problems as a generalized multiple measurement vector compressive sensing problem. Furthermore, by leveraging the structured sparsity of the uplink channel matrix, we propose an efficient generalized multiple measurement vector approximate message passing (GMMV-AMP) algorithm to realize simultaneous active user detection and channel estimation based on a spatial domain or an angular domain channel model. To jointly exploit the channel sparsity presented in both the spatial and the angular domains for enhanced performance, a Turbo-GMMV-AMP algorithm is developed for detecting the active users and estimating their channels in an alternating manner. Finally, an adaptive access scheme is proposed, which adapts the access latency to guarantee reliable massive access for practical systems with unknown channel sparsity level. Additionally, the state evolution of the proposed GMMV-AMP algorithm is derived to predict its performance. Simulation results demonstrate the superiority of the proposed active user detection and channel estimation schemes compared to several baseline schemes.

Tables

Table 1. TABLE I: Computational Complexity for AUD and CE

Algorithm	Number of complex multiplications in each iteration
GSP	$2\left(G+1\right)KMP+G{M}^{2}P+2\left(M+1\right)GP{K}_a^2+2P{K}_a^3$
SOMP	$\left(2G+1\right)KMP+G{M}^{2}P+\left(M+1\right)GP{o}^2+P{o}^3$
DSAMP	$\left(2G+3\right)KMP+MP+2\left(M+1\right)GP{u}^2+2P{u}^3$
GMMV-AMP	$4GKMP+3GKP+16GMP+20KMP$
Turbo-GMMV-AMP	$2T_{\mathrm{amp}}\left(4GKMP+3GKP+16GMP+20KMP\right)$

Table 2. TABLE II: Considered Schemes

Scheme		Model	Algorithm	$G$
Simultaneous AUD and CE	Scheme 1	(11)	GMMV-AMP	fixed
	Scheme 2	(12)	GMMV-AMP	fixed
Alternating AUD and CE	Scheme 3	(11) and (12)	Turbo-GMMV-AMP	fixed
	Scheme 4	(11) and (12)	Turbo-GMMV-AMP	adaptive
Baseline Schemes		(11)	GSP [32], SOMP [31], and DSAMP [21]	fixed

Equations

$${\mathbf{y}}_{p,k}^t={\mathbf{h}}_{p,k}{s}_{p,k}^t+{\mathbf{n}}_p^t,$$

$${\mathbf{y}}_p^t=\sum _{k=1}^K{\alpha }_k{\mathbf{h}}_{p,k}{s}_{p,k}^t+{\mathbf{n}}_p^t={\mathbf{H}}_p{\mathbf{s}}_p^t+{\mathbf{n}}_p^t,$$

$${\tilde{\mathbf{h}}}_{p,k}=\sum _{l=1}^L{\beta }_{k,l}{\mathbf{a}}_R\left({\varphi }_{k,l}\right)e^{-j2\pi {\varpi }_{k,l}\left(-\frac{B_s}{2}+\frac{B_s\left(pN/P-1\right)}{N}\right)},$$

$${\left|\mathrm{supp}\left\{{\left[{\mathbf{X}}_p\right]}_{:,m}\right\}\right|}_c=K_a\ll K.$$

$$\mathrm{supp}\left\{{\left[{\mathbf{X}}_p\right]}_{:,1}\right\}=\mathrm{supp}\left\{{\left[{\mathbf{X}}_p\right]}_{:,2}\right\}=\cdots =\mathrm{supp}\left\{{\left[{\mathbf{X}}_p\right]}_{:,M}\right\}.$$

$$\mathrm{supp}\left\{{\mathbf{X}}_1\right\}=\mathrm{supp}\left\{{\mathbf{X}}_2\right\}=\cdots =\mathrm{supp}\left\{{\mathbf{X}}_P\right\}.$$

$${\tilde{\mathbf{w}}}_{p,k}={\mathbf{A}}_R^{\mathrm{H}}{\tilde{\mathbf{h}}}_{p,k},$$

$${\left|\mathrm{supp}\left\{{\tilde{\mathbf{w}}}_{p,k}\right\}\right|}_c\ll M,$$

$$\mathrm{supp}\left\{{\tilde{\mathbf{w}}}_{1,k}\right\}=\mathrm{supp}\left\{{\tilde{\mathbf{w}}}_{2,k}\right\}=\cdots =\mathrm{supp}\left\{{\tilde{\mathbf{w}}}_{P,k}\right\}.$$

$$\mathrm{supp}\left\{{\mathbf{W}}_1\right\}=\mathrm{supp}\left\{{\mathbf{W}}_2\right\}=\cdots =\mathrm{supp}\left\{{\mathbf{W}}_P\right\}.$$

$${\mathbf{Y}}_p^G={\mathbf{S}}_p^G{\mathbf{X}}_p+{\mathbf{N}}_p,\forall p\in \left[P\right],$$

$${\mathbf{R}}_p^G={\mathbf{Y}}_p^G{\mathbf{A}}_R^{\ast }={\mathbf{S}}_p^G{\mathbf{W}}_p+{\tilde{\mathbf{N}}}_p,\forall p\in \left[P\right],$$

$${\hat{x}}_{k,m}=\int x_{k,m}p\left(x_{k,m}|\mathbf{Y}\right)d{x}_{k,m},\forall k,m.$$

$$p\left(x_{k,m}|\mathbf{Y}\right)=\int p\left(\mathbf{X}|\mathbf{Y}\right)d{\mathbf{X}}_{\backslash k,m},$$

$$p\left(\mathbf{X}|\mathbf{Y}\right)$$

$$=\frac{1}{{\tilde{Z}}_1}\prod _{m=1}^M\left[\prod _{g=1}^{G}p\left(y_{g,m}|\mathbf{X}\right)\prod _{k=1}^K{p}_0\left(x_{k,m}\right)\right],$$

$$p\left(y_{g,m}|\mathbf{X}\right)=\frac{1}{\pi \sigma }\exp \left(-\frac{1}{\sigma }{\left|y_{g,m}-{\sum }_k{s}_{g,k}{x}_{k,m}\right|}^2\right),$$

$$p_0\left(\mathbf{X}\right)$$

$$=\prod _{m=1}^M\prod _{k=1}^K\left[\left(1-{\gamma }_{k,m}\right)\delta \left(x_{k,m}\right)+{\gamma }_{k,m}f\left(x_{k,m}\right)\right],$$

$$p\left(x_{k,m}|\mathbf{Y}\right)$$

$$\approx \frac{1}{{\tilde{Z}}_2}{p}_0\left(x_{k,m}\right)\mathcal{C}\mathcal{N}\left(x_{k,m};C_{k,m}^q,D_{k,m}^q\right),$$

$$D_{k,m}^q$$

$$C_{k,m}^q$$

$$V_{g,m}^q$$

$$Z_{g,m}^q$$

$$p\left(x_{k,m}|C_{k,m}^q,D_{k,m}^q\right)$$

$$+{\pi }_{k,m}^q\mathcal{C}\mathcal{N}\left(x_{k,m};A_{k,m}^q,B_{k,m}^q\right),$$

$$A_{k,m}^q$$

$${\pi }_{k,m}^q$$

$$\mathcal{L}$$

$$g_a\left(C_{k,m}^q,D_{k,m}^q\right)$$

$$g_c\left(C_{k,m}^q,D_{k,m}^q\right)$$

$$Q\left(\mathbf{\theta },{\mathbf{\theta }}^q\right)$$

$${\mathbf{\theta }}^{q+1}$$

$${\mu }_{k,m}^{q+1}$$

Full text

I Introduction

Video streaming, social networking, and the emerging Internet-of-Things (IoT) accelerate the development of base stations (BSs) that enable connectivity for billions of user equipments (UEs) with massive data volumes [2]. However, reliable massive access for massive connectivity is not supported by the current wireless networks [3].

To guarantee the availability of resources and the quality of service in massive access scenarios, uplink systems have to provide ultra-reliable low-latency detection and channel estimation (CE) for active UEs [2]. Conventional grant-based random access (RA) protocols require control signaling and the scheduling of uplink access requests for granting of resources [4, 5, 6, 7]. The physical random access channel (PRACH) protocol of Long-Term Evolution is one example of grant-based protocols, and can be classified into two categories: contention-free RA and contention-based RA [4]. For contention-free RA, the BS first allocates UEs the dedicated preambles, which are then transmitted by the active UEs, and finally the BS responds to the requesting UEs without further contention resolution [5]. For contention-based RA, multiple active UEs first transmit preambles selected from a predefined sequence set to access the BS. Further contention resolution is required if multiple UEs choose the same preamble [6]. Unfortunately, for massive access, collisions are likely to occur as the number of potential UEs can be much larger than the number of available preambles. The authors in [7] proposed a strongest-user collision resolution protocol to resolve collisions in overloaded networks. However, such grant-based solutions generally suffer from high access latency and require complicated collision resolution schemes in massive access scenarios [8].

As a promising alternative, grant-free RA protocols have recently attracted significant attention, where each active UE directly transmits its pilots and data to the BS without waiting for permission [9, 10]. Allocating orthogonal channel resources (e.g., using orthogonal pilots as in [11]) can facilitate the detection of the active UEs at the BS and the estimation of their channels. However, for massive numbers of potential UEs, this approach fails due to the limited number of available orthogonal channels for a given channel coherence time. Fortunately, a key characteristic of massive access in future wireless networks is the sporadic traffic of the UEs, i.e., out of the many potential UEs, only a small number are activated and want to access the network in any given time interval [7].

Exploiting this sporadic traffic property, several compressive sensing (CS)-based grant-free RA schemes have been proposed, where the active user detection (AUD) is formulated as a sparse signal recovery problem [13, 12]. In [14, 15], two advanced CS-based multi-user detection schemes, which leverage the structured sparsity over multiple time slots for accurate support detection, were proposed to jointly detect the active UEs and to decode their data. To further exploit the a priori information about the transmitted discrete symbols for RA, the authors in [16] developed a joint approximate message passing (AMP) and expectation maximization (EM) algorithm to improve the detection performance of sparsely active UEs. Furthermore, the authors in [17] proposed a threshold aided block sparsity adaptive subspace pursuit algorithm, which enables improved sparse signal recovery. This work was also extended to cloud radio access networks (C-RANs) for efficient UE activity and data detection [18, 19]. However, the solutions in [14, 15, 16, 17, 18, 19] rely on the availability of perfect channel state information (CSI), which is difficult to obtain, especially for massive wireless-connected UEs.

To jointly perform AUD and CE for single-antenna BSs, blind detection of sparse code multiple access was proposed to support grant-free RA in massive access scenarios in [20]. For multi-antenna systems, the authors in [22] proposed a novel joint AUD and CE scheme, where the sparsity of delay-domain channel impulse response (CIR) is leveraged for facilitating CE. Particularly, in [22], by iteratively exchanging the active user information and CIR estimates, an identified user cancelation technique is further employed for enhanced performance. However, the delay-domain CIR sparsity highly depends on the channel environment and may be violated in some practical scenarios. Considering frequency-domain CE, a modified Bayesian compressive sensing (BCS)-based access scheme for uplink C-RANs was proposed in [21], where the structured sparsity over multiple receive antennas are exploited. To reduce the computational complexity, the authors in [23] and [24] developed an AMP-based scheme for AUD and CE for massive access in massive multiple-input multiple-output (MIMO) systems. However, the solutions in [23] and [24] require the full knowledge of the a priori distribution of the channels and the noise variance, which might not be available in practice. Besides, the work in [20, 21, 22, 23, 24] considers a narrow-band massive access scenario assuming single-carrier transmission.

Besides the aforementioned work, there are also related solutions that address the massive access problem with sparse user activity from information-theoretical perspectives [25, 26]. Especially in [26], the authors provided an elegant solution to support massive access with fully non-coherent detection.

In this paper, we consider massive access for the more challenging enhanced mobile broadband (eMBB) scenario, and investigate AUD and CE for uplink massive MIMO orthogonal frequency division multiplexing (OFDM) systems. By exploiting the sporadic traffic of the UEs and the virtual angular domain sparsity of massive MIMO channels, we develop a CS-based adaptive AUD and CE scheme. Specifically, a pilot design based on distributed CS (DCS) theory is proposed for broadband massive access. Moreover, the AUD and CE problems at the BS are formulated as a generalized multiple measurement vector (GMMV) CS problem [27]. By leveraging the structured sparsity of the uplink channel matrix, we propose a GMMV-AMP algorithm for efficient simultaneous AUD and CE based on a spatial domain or an angular domain channel model. To further improve performance, a Turbo-GMMV-AMP algorithm is proposed for detecting the active UEs and estimating their channels in an alternating manner. This forms the basis for an adaptive access scheme which allows the adaptation of the required access latency111For grant-free massive access, the time slots consumed for transmitting the access pilot sequence contribute to the major part of the access latency. In this paper, we focus on the time slot overhead for pilot transmission. to the sparsity level of the uplink channel matrix, i.e., the maximum number of non-zero entries of its columns. Additionally, the state evolution (SE) of the proposed GMMV-AMP algorithm is derived to characterize its performance. Our main contributions can be summarized as follows.

• •

DCS theory-based pilot design tailored for multi-carrier systems: Previous work [14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24] mainly focuses on frequency-flat narrow-band massive access with single-carrier transmission. In contrast, we consider the more challenging massive access problem for eMBB, where OFDM is employed. Based on DCS theory [27], we design pseudo-random pilot sequences tailored for multi-carrier systems. Thereby, the structured channel sparsity for the different subcarriers is further leveraged to improve AUD and CE performance.

• •

GMMV-AMP algorithm: Exploiting the structured sparsity of the massive access channel matrix observed at multiple receive antennas and multiple subcarriers, the proposed GMMV-AMP algorithm facilitates efficient simultaneous AUD and CE. In particular, exploiting the EM algorithm, the GMMV-AMP algorithm can learn the unknown hyper-parameters of the a priori distribution of the channels and the noise variance.

• •

Turbo-GMMV-AMP algorithm: To jointly leverage the channel sparsity presented in both the spatial and the angular domains, this algorithm performs AUD and CE in an alternating manner for further enhanced performance. Compared with GMMV-AMP-based simultaneous processing methods and the state-of-the-art solutions, the proposed alternating approach will reap a significant reduction of access latency for massive access.

• •

CS-based adaptive AUD and CE: Most prior work [14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24] leverages the UEs’ sporadic traffic only to provide a fixed access latency. In contrast, the access latency for the proposed adaptive access scheme can be adapted to the actual sparsity level of the uplink massive access channel matrix for reliable AUD and CE.

Notations: Throughout this paper, scalar variables are denoted by normal-face letters, while boldface lower and upper-case letters denote column vectors and matrices, respectively. $\left[{\bf{X}}_{p}\right]_{k,m}$ is the $(k,m)$-th element of matrix ${{\bf{X}}_{p}}{\in\mathbb{C}^{K\times M}}$; $\left[{\bf{X}}_{p}\right]_{k,:}$ and $\left[{\bf{X}}_{p}\right]_{:,m}$ are the $k$-th row vector and the $m$-th column vector of matrix ${{\bf{X}}_{p}}{\in\mathbb{C}^{K\times M}}$, respectively. The transpose, complex conjugate, and conjugate transpose operators are denoted by ${(\cdot)^{\rm{T}}}$, ${(\cdot)^{*}}$, and ${(\cdot)^{\rm{H}}}$, respectively. $\left|{\cal K}\right|_{c}$ is the number of elements in set ${\cal K}$, $\left[K\right]$ denotes the set $\left\{{1,\cdots,K}\right\}$, and $\rm{supp}\left\{\cdot\right\}$ is the support set of a vector or a matrix. ${\mathbb{E}}\left[\cdot\right]$ denotes the statistical expectation. $\emptyset$ denotes the empty set and ${\bf{0}}_{K\times M}$ is the ${K\times M}$ zero matrix. Finally, ${\cal C}{\cal N}\left({x;\mu,v}\right)$ denotes the complex Gaussian distribution of a random variable $x$ with mean $\mu$ and variance $v$.

II System Model

In this section, we introduce the system model of the uplink massive access in massive MIMO-OFDM systems. Moreover, the sparsity properties of the massive access channel matrix presented in the spatial and the angular domains are further explained, respectively.

II-A Uplink Massive Access in Massive MIMO Systems

We consider the typical uplink massive access scenario in massive MIMO systems, as illustrated in Fig. 1. There are one BS equipped with an $M$-antenna uniform linear array (ULA) and $K$ potential UEs, where $K$ is usually large (e.g., $K=10^{3}$ in [24]). OFDM with $N$ subcarriers is adopted to combat time dispersive channels, and $P$ pilots are uniformly allocated across the $N$ subcarriers. For the subchannel of the $p$-th pilot subcarrier $\left(1\leq p\leq P\right)$, the signal ${\bf y}_{p,k}^{t}\in\mathbb{C}^{M\times 1}$ received at the BS from the $k$-th UE in the $t$-th time slot (or equivalently the $t$-th OFDM symbol) can be expressed as

[TABLE]

where ${\bf h}_{p,k}\in\mathbb{C}^{M\times 1}$ is the subchannel associated with the $k$-th UE, $s_{p,k}^{t}$ is the uplink access pilot of the $k$-th UE, and ${\bf n}_{p}^{t}$ denotes the additive white Gaussian noise (AWGN) at the BS for the $p$-th pilot subcarrier and the $t$-th time slot. Here, without loss of generality, we consider single-antenna UEs. For a typical massive access scenario, within a given time interval, only a small number of UEs are activated to access the BS. The UE activity indicator is denoted as $\alpha_{k}$, and is equal to 1 when the $k$-th UE is active and 0 otherwise. Meanwhile, we define the set of active UEs as ${\cal A}=\left\{k|{\alpha_{k}}=1,\;1\leq k\leq K\right\}$, and the number of active UEs is denoted by $K_{a}=\left|\cal A\right|_{c}$. Hence, the signal received at the BS from all active UEs for the $p$-th pilot subcarrier and the $t$-th time slot is given as follows

[TABLE]

where ${\bf H}_{p}=\left[{\alpha_{1}}{{\bf h}_{p,1}},\cdots,{\alpha_{K}}{{\bf h}_{p,K}}\right]\in\mathbb{C}^{M\times K}$ and ${\bf s}_{p}^{t}=\left[s_{p,1}^{t},\cdots,s_{p,K}^{t}\right]^{\rm T}\in\mathbb{C}^{K\times 1}$. By considering both the large-scale and the small-scale fading, we can model ${\bf h}_{p,k}$ as ${\bf h}_{p,k}={\rho_{k}}{\widetilde{\bf h}_{p,k}}$, where $\rho_{k}$ is the large-scale fading caused by path loss and shadowing, and $\widetilde{\bf h}_{p,k}$ is the small-scale fading. For the $p$-th pilot subcarrier, the subchannel of the $k$-th UE is modeled as follows [32]

[TABLE]

where $N/P$ is an integer, $L$ denotes the number of multi-path components (MPCs), $\beta_{k,l}$ and $\varpi_{k,l}$ are the complex path gain and the path delay of the $l$-th MPC, respectively, and $B_{s}$ is the two-sided bandwidth. The array response vector ${{\bf a}_{R}}\left(\phi_{k,l}\right)$ is given by ${{\bf a}_{R}}\left(\phi_{k,l}\right)=\left[1,e^{-j2\pi{\phi_{k,l}}},\cdots,e^{-j2\pi\left(M-1\right){\phi_{k,l}}}\right]^{\rm T}$, where $\phi_{k,l}=\frac{d}{\lambda}{\sin\left(\varphi_{k,l}\right)}$. Here, $\varphi_{k,l}$ is the angle of arrival (AOA) of the $k$-th UE’s $l$-th MPC, $\lambda$ is the wavelength, and $d=\lambda/2$ is the antenna spacing.

II-B Space-Frequency Structured Sparsity in Massive Access

Due to the sporadic traffic of the UEs, only a small number of UEs are active, i.e., $K_{a}\ll K$. Thus, by defining ${\bf X}_{p}={\bf H}_{p}^{\rm T}$, the channel vector $\left[{\bf X}_{p}\right]_{:,m}\in\mathbb{C}^{K\times 1}$ observed at the $m$-th receive antenna for the $p$-th pilot subcarrier is sparse as

[TABLE]

Moreover, all BS antennas exhibit the same sparsity,

[TABLE]

We refer to this property as the spatial domain structured sparsity of massive access. Since the $\alpha_{k},\forall k$, are identical for all subchannels, the $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$ also exhibit a common sparsity pattern in the frequency domain as follows

[TABLE]

The joint structured sparsity in (5) and (6) is referred to as the space-frequency structured sparsity of $\left\{{{{\mathbf{X}}_{p}}}\right\}_{p=1}^{P}$. To illustrate this structured sparsity in Fig.2(a), as an example, we assume that $K_{a}=3$ active UEs out of $K=10$ total UEs access the BS, which is equipped with $M=16$ antennas.

II-C Angular-Frequency Structured Sparsity in Massive MIMO

On the other hand, the BS is usually at high elevation with few scatterers around, whereas the UEs are typically located at low elevation in a rich local scattering environment far from the BS [27]. This scenario can be modeled using the classical one-ring channel model [28]. For a UE which is located at a distance of $R$ from the BS and surrounded by rich scatterers located within a radius of $r$ around the UE, the angular spread $\Delta\approx{\rm arctan}\left(r/R\right)$ seen from the BS is expected to be very small as usually $R\gg r$. This leads to sparsity of massive MIMO channels in the virtual angular domain [29, 30, 31]. Specifically, the virtual angular domain massive MIMO channel associated with the $k$-th UE for the $p$-th pilot subcarrier can be represented as

[TABLE]

where the transformation matrix ${\bf A}_{R}\in\mathbb{C}^{M\times M}$ at the BS side is a unitary matrix. Here, ${\bf A}_{R}$ depends on the geometry of the adopted array, and becomes the discrete Fourier transform matrix for a ULA with $d=\lambda/2$ [32]. Due to the small $\Delta$ and large $M$, the channel vector $\widetilde{\bf w}_{p,k}$ is sparse, i.e.,

[TABLE]

and this sparsity is clustered, as illustrated in Fig. 2(b). Moreover, since the spatial propagation characteristics of all wireless channels within the total bandwidth are similar, all subchannels associated with different subcarriers are affected by the same scatterers [27]. Consequently, the $\left\{\widetilde{\bf w}_{p,k}\right\}_{p=1}^{P},\forall k$, have a common sparsity pattern in the frequency domain, i.e.,

[TABLE]

We refer to the jointly structured sparsity in (8) and (9) as the angular-frequency structured sparsity of massive MIMO channels. Additionally, we further define the virtual angular domain channel matrix as ${\bf W}_{p}={{\bf X}_{p}}{{\bf A}_{R}^{*}}=\left[{\alpha_{1}}{{\bf w}_{p,1}},\cdots,{\alpha_{K}}{{\bf w}_{p,K}}\right]^{\rm T}$, where ${\bf w}_{p,k}={\rho_{k}}{\widetilde{\bf w}_{p,k}}$. Considering the space-frequency structured sparsity described in (4)-(6), we further have $\left|{\rm supp}\left\{\left[{\bf W}_{p}\right]_{:,m}\right\}\right|_{c}\ll K_{a}$, and

[TABLE]

The illustration in Fig. 2(b) takes both the space-frequency structured sparsity of massive access and the angular-frequency structured sparsity of massive MIMO channels into account. These sparsity properties will be exploited in the remainder of this paper to achieve low-latency and highly-reliable AUD and CE performance.

III CS-Based Active User Detection and Channel Estimation Schemes

In this section, we detail the proposed CS-based AUD and CE schemes for massive access. First, a DCS-based pilot design is proposed for broadband massive access. Then, two simultaneous AUD and CE schemes are developed based on a spatial domain and an angular domain channel model, respectively. On this basis, the alternating AUD and CE schemes are further proposed for enhanced performance. Finally, the computational complexity of the proposed schemes will be analyzed.

The frame structure of the uplink signals is illustrated in Fig. 3. A frame consists of $T$ time slots, where the first $G$ time slots include both pilots and data, and the remaining $(T-G)$ time slots are reserved for data transmission only. Here, we assume $T$ is smaller than the channel coherence time, and the activity of the UEs during the $T$ time slots remains unchanged. At the BS, the received signals in $G$ successive time slots for the $p$-th pilot subcarrier are collected as

[TABLE]

where ${\bf Y}_{p}^{G}=\left[{\bf y}_{p}^{1},\cdots,{\bf y}_{p}^{G}\right]^{\rm T}\in{\mathbb{C}}^{G\times M}$, ${\bf S}_{p}^{G}=\left[{\bf s}_{p}^{1},\cdots,{\bf s}_{p}^{G}\right]^{\rm T}\in{\mathbb{C}}^{G\times K}$, ${\bf X}_{p}={\bf H}_{p}^{\rm T}\in{\mathbb{C}}^{K\times M}$, and ${\bf N}_{p}=\left[{\bf n}_{p}^{1},\cdots,{\bf n}_{p}^{G}\right]^{\rm T}$. To avoid complicated scheduling protocols and the associated latencies, in RA, the active UE set (AUS) $\cal A$ and the corresponding channel vectors $\left\{{\bf h}_{p,k}\right\}_{p=1}^{P}$, $k\in{\cal A}$, have to be reliably estimated based on the noisy measurements $\left\{{\bf Y}_{p}^{G}\right\}_{p=1}^{P}$ and the known pilot matrices $\left\{{\bf S}_{p}^{G}\right\}_{p=1}^{P}$, which is equivalent to estimating $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$ based on (11). Due to the space-frequency structured sparsity of $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$, the AUD and CE based on (11) can be formulated as a CS problem with $G\ll K$. Moreover, considering the virtual angular domain sparsity shown in (8), we can further transform (11) as

[TABLE]

where ${\widetilde{\bf N}_{p}}={{\bf N}_{p}}{{\bf A}_{R}^{*}}$. Based on this, we develop two categories of AUD and CE schemes:

• •

Simultaneous AUD and CE: The estimate of $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$, denoted by $\{\widehat{\bf X}_{p}\}_{p=1}^{P}$, can be directly obtained based on (11). Alternatively, we can first estimate $\left\{{\bf W}_{p}\right\}_{p=1}^{P}$ based on (12), and then obtain $\{{\widehat{\bf X}_{p}}\}_{p=1}^{P}$ according to (7).

• •

Alternating AUD and CE: Compared with ${\bf X}_{p}$ in (11), the sparser ${\bf W}_{p}$ in (12) will yield a better CE performance, but the common sparsity pattern across multiple columns of ${\bf W}_{p}$ is destroyed. Therefore, the proposed alternating scheme leverages (11) for AUD and (12) for CE, i.e., (11) and (12) are alternately exploited to reap both the structured sparsity of ${\bf X}_{p}$ and the enhanced sparsity of ${\bf W}_{p}$ for further improved performance.

In the following, we will first discuss the pilot design and then explain the proposed AUD and CE schemes.

III-A DCS-Based Pilot Design for Broadband Massive Access

For the CS problems in (11) and (12), the RA pilot matrices ${\bf S}_{p}^{G}$, $\forall p$, serve as measurement matrix. The properties of the measurement matrix are crucial for guaranteeing reliable recovery of sparse channel matrices $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$. Hence, the pilot signals should be carefully designed to guarantee reliable AUD and CE. The sparse signal recovery algorithms proposed in this paper are based on the family of AMP algorithms, which usually require independent and identically distributed (i.i.d.) Gaussian measurement matrices [33]. Hence, for the $p$-th pilot subcarrier, the pilot associated with the $k$-th UE in the $t$-th time slot is generated from a standard complex Gaussian distribution, i.e., $s_{p,k}^{t}\sim{\cal CN}(s_{p,k}^{t};0,1)$. Furthermore, the ${\bf S}_{p}^{G}$ should be different for different pilot subcarriers to achieve diversity, which means that (11) and (12) are GMMV-CS models [27]. Compared to the conventional multiple measurement vector (MMV) problem [34], where identical pilots would be allocated to all pilot subcarriers, employing different pilot matrices across different pilot subcarriers can improve AUD and CE performance according to DCS theory [27].

III-B Simultaneous AUD and CE Schemes

For massive access, reliable inference of AUS $\cal A$ and the corresponding $\left\{{\bf h}_{p,k}\right\}_{p=1}^{P}$, $k\in{\cal A}$, from $\left\{{\bf Y}_{p}^{G}\right\}_{p=1}^{P}$ is challenging. In this subsection, we propose two simultaneous AUD and CE schemes, where a GMMV-AMP algorithm is developed to solve the CS problems in (11) and (12), respectively. Without loss of generality, we consider (11) first and focus on the $p$-th pilot subcarrier. The obtained results can be easily extended to the model (12) and multiple pilot subcarriers.

1) Spatial Domain Simultaneous AUD and CE (Scheme 1): Define $x_{p,k,m}=\left[{\bf X}_{p}\right]_{k,m}$. The minimum mean square error (MMSE) estimate of ${\bf X}_{p}$ is the posterior mean, which can be expressed as

[TABLE]

In (13), the superscript $G$ and index $p$ in $x_{p,k,m}$, ${\bf X}_{p}$, and ${\bf Y}_{p}^{G}$ are dropped to simplify the notation, and the marginal posterior distribution is given by

[TABLE]

where ${\bf X}_{\backslash k,m}$ denotes the collection of the $\left\{x_{i,j}\right\}_{1\leq j\leq M,j\neq m}^{1\leq i\leq K,i\neq k}$. The joint posterior distribution in (14) can be computed according to the Bayesian rule as

[TABLE]

where ${\widetilde{Z}_{1}}=\iint p\left({\bf Y}|{\bf X}\right)p_{0}\left({\bf X}\right)d{\bf X}d{\bf Y}$ is a normalization factor and $p_{0}\left({\bf X}\right)$ is the a priori distribution of ${\bf X}$. Under the assumption of AWGN, the likelihood function in (15) is

[TABLE]

where $\sigma$ is the variance of the complex AWGN. In this paper, to characterize the sparsity of ${\bf X}$, we consider a flexible spike and slab a priori distribution for ${\bf X}$, i.e.,

[TABLE]

which can effectively capture the actual prior knowledge of the channel matrix $\bf X$ [29, 30, 31]. In (17), $0<\gamma_{k,m}<1$ is the sparsity ratio, i.e., the probability of $x_{k,m}$ being non-zero, $\delta\left(\cdot\right)$ is the Dirac delta function, and $f\left(\cdot\right)$ is the distribution of the non-zero entries. This distribution arises from the literature of AMP algorithm [35, 36], and has been successfully employed in various AMP-based channel estimation schemes [29, 30, 31].

The factorization in (15) can be represented by a bipartite graph, which consists of variable nodes, factor nodes, and the corresponding edges [37]. This suggests the use of message passing algorithms [37] to realize the MMSE estimator. As the messages for marginal posterior probabilities are difficult to compute for massive access, we resort to the AMP algorithm [33], which employs low-complexity heuristics for approximating $p\left(x_{k,m}|{\bf Y}\right)$.

Proposition 1

In the large system limit, i.e., $K\rightarrow\infty$, while $\gamma=K_{a}/K$ and $\kappa=G/K$ are fixed, the AMP algorithm decouples the matrix estimation problem based on (11) into $KM$ scalar estimation problems. Considering this, the posterior distributions of $x_{k,m}$, $\forall k,m$, are approximated as

[TABLE]

where $q$ denotes the $q$-th iteration, and ${\widetilde{Z}_{2}}$ is a normalization factor. In (18), $D_{k,m}^{q}$ and $C_{k,m}^{q}$ are updated at the variable nodes of the bipartite graph as

[TABLE]

where $V_{g,m}^{q}$ and $Z_{g,m}^{q}$ are updated at the factor nodes of the bipartite graph as

[TABLE]

Proof:

Please refer to the Appendix A. ∎

It is worth noticing that, although the proposed algorithms are developed from the large system limit ($K\rightarrow\infty$), in practice, they perform well even in the medium size problems, such as massive access with hundreds even thousands of UEs, which has been discussed in the literature of AMP algorithms [33, Sec. II], [41]. Moreover, we consider the widely used Gaussian a priori distribution for the channel gains, i.e., $f\left(x_{k,m}\right)={\cal CN}\left(x_{k,m};\mu,\tau\right)$ [29]. By exploiting this a priori model in (18), the posterior distribution of $x_{k,m}$ is obtained as

[TABLE]

where

[TABLE]

and $\pi_{k,m}^{q}$ is referred to as the belief indicator. The posterior mean (62) and variance (63) can now be explicitly calculated as

[TABLE]

respectively.

Therefore, for the $p$-th pilot subcarrier, the MMSE estimate of ${\bf X}_{p}$ can be acquired by iteratively calculating (19)-(28) instead of solving the high-dimensional integrals in (14). The resulting procedure is referred to as the basic MMV-AMP algorithm. However, the basic MMV-AMP algorithm requires full knowledge of the a priori distribution of the channels $\left\{\mu,\tau,\gamma_{k,m},\forall k,m\right\}$ and the noise variance $\sigma$, which may be difficult to obtain in practice. Hence, the EM algorithm is exploited to learn the unknown hyper-parameters, i.e., $\bm{\theta}=\left\{\mu,\tau,\sigma,\gamma_{k,m},\forall k,m\right\}$. The EM algorithm involves two steps

[TABLE]

where ${\mathbb{E}}\left[\cdot|{\bf Y};{\bm{\theta}}^{q}\right]$ denotes the expectation conditioned on measurements ${\bf Y}$ with parameters ${\bm{\theta}}^{q}$, i.e., the expectation is with respect to the posterior distribution $p\left({\bf X}|{\bf Y};{\bm{\theta}}^{q}\right)$. There are two challenges in harnessing the EM algorithm: (a) the computation of $p\left({\bf X}|{\bf Y};{\bm{\theta}}^{q}\right)$ is of high complexity and (b) the joint optimization of all elements of ${\bm{\theta}}$ is difficult. Fortunately, in the large system limit with $K\to\infty$, the high complexity of calculating $p\left({\bf X}|{\bf Y};{\bm{\theta}}^{q}\right)$ can be considerably reduced by using the approximation $p\left({\bf X}|{\bf Y};{\bm{\theta}}^{q}\right)\approx\prod\limits_{m}{\prod\limits_{k}{p\left(x_{k,m}|C_{k,m}^{q},D_{k,m}^{q}\right)}}$ according to (18). Moreover, the incremental EM algorithm [38] can be used to simplify the joint optimization of all elements of ${\bm{\theta}}$, where ${\bm{\theta}}$ is updated one element at a time and the other parameters are held constant. By setting the derivative of (29) with respect to one element of ${\bm{\theta}}$ to zero, the update rules of the hyper-parameters are obtained as, $\forall k,m$:

[TABLE]

As the EM algorithm may converge to a local extremum of the likelihood function, the proper initialization of the hyper-parameters is crucial. Here, we use the following initialization [39], $\forall k,m$:

[TABLE]

Here, $\Phi(-c)$ and $\phi(c)$ are the cumulative distribution function and the probability density function of the standard normal distribution, respectively, and ${\rm SNR}^{0}\!=\!100$ is suggested in [39].

Equations (19)-(28) and (31)-(37) are the main EM steps incorporated into the MMV-AMP algorithm to learn the unknown hyper-parameters. However, the resulting overall algorithm is limited to the AUD and CE of a single pilot subcarrier. Hence, we extend the MMV-AMP algorithm to the GMMV-AMP algorithm (summarized in Algorithm 1), where the subchannel matrices ${\bf X}_{p}$, $\forall p$, for all pilot subcarriers are jointly estimated with different measurement matrices ${\bf S}_{p}^{G},\forall p$. Specifically, in lines 3-6, the messages are updated independently for all pilot subcarriers; in lines 4 and 5, a damping parameter $\rho=0.3$ is used to prevent the algorithm from diverging according to [40]; line 7 uses the incremental EM algorithm to learn the unknown hyper-parameters $\bm{\theta}$; line 8 refines the update rule for the sparsity ratio ${\gamma}_{p,k,m}$ to leverage the structured sparsity of the channel matrix for improved CS recovery. By contrast, the state-of-the-art AMP-based estimators in [23] and [24] require $\bm{\theta}$ as a priori information.

In Algorithm 1, the sparsity ratio $\gamma_{p,k,m}$ is the probability that the $\left(k,m\right)$-th element of ${\bf X}_{p}$ is non-zero. In line 7, $\gamma_{p,k,m}$ is updated independently for all $p$, $k$, and $m$ according to (32), which indicates that the common sparsity pattern described in (5) and (6) is not exploited. To fully exploit the structured sparsity of the channel matrix, as discussed in Section II-B and illustrated in Fig. 2(a), we assume that the channel elements associated with the same UE have a common sparsity ratio, and further propose to refine $\gamma_{p,k,m}$ as in line 8 of Algorithm 1, where we use

[TABLE]

With the estimate of $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$, the AUS and the corresponding channel vectors can be simultaneously acquired. Specifically, for AUD, we develop two UE activity detectors based on the $\{\widehat{\bf X}_{p}\}_{p=1}^{P}$ and the belief indicators $\pi_{p,k,m}$, $\forall p,k,m$, respectively, as follows. We first define a threshold function $r\left(x;\varepsilon\right)$, where $r\left(x;\varepsilon\right)=1$ if $\left|x\right|>{\varepsilon}$, otherwise $r\left(x;\varepsilon\right)=0$.

Definition 1

Based on $\{\widehat{\bf X}_{p}\}_{p=1}^{P}$, a channel gain-based activity detector (CG-AD) is proposed for AUD as follows

[TABLE]

where $\varepsilon_{\rm cg}={\widetilde{\alpha}}\max\left\{\left|\hat{x}_{p,k,m}\right|,\;\forall p,k,m\right\}$, $\widetilde{\alpha}=0.01$ is the ratio of the minimum and maximum amplitudes of the channel coefficients, see [23, Sec. IV], and $p_{\rm cg}=0.9$222If more than $90\%$ of the elements of $\{[\widehat{\bf X}_{p}]_{k,:}\}_{p=1}^{P}$ are decided to be non-zero, the $k$-th UE is declared active..

Proposition 2

In the large system limit, if a reliable estimate of ${\bf X}_{p}$ is acquired after the convergence of the GMMV-AMP algorithm,

[TABLE]

Proof:

Please refer to the Appendix B. ∎

Definition 2

Since the belief indicator $\pi_{p,k,m}$ tends to be 1 for ${\hat{x}}_{p,k,m}\neq 0$ and 0 for ${\hat{x}}_{p,k,m}=0$, we further design a belief indicator-based activity detector (BI-AD) as

[TABLE]

For spatial domain simultaneous AUD and CE, we set $\varepsilon_{\rm bi}=\varepsilon_{\rm bi}^{\rm spa}=0.5$ to make the missed detection and false alarm probabilities identical, and the same as $p_{\rm cg}$, we set $p_{\rm bi}=p_{\rm bi}^{\rm spa}=0.9$ based on empirical experience. Nevertheless, we note that the decisions of the CG-AD and BI-AD are not sensitive to the values of $p_{\rm cg}$ and $p_{\rm bi}$333We found empirically via simulations for a wide range of system parameters that $p_{\rm cg}=p_{\rm bi}^{\rm spa}=0.9$ yields a high AUD performance.. Finally, if the $k$-th UE is declared active, its channel is estimated as ${\widehat{\bf h}_{p,k}}=\left[{{{\widehat{\mathbf{X}}}_{p}}}\right]_{k,:}^{\text{T}}$.

2) Angular Domain Simultaneous AUD and CE (Scheme 2): The GMMV-AMP algorithm designed for CS model (11) can be directly applied to model (12) for angular domain simultaneous AUD and CE by replacing ${\bf Y}_{p}^{G}$ and $\widehat{\bf X}_{p}^{G}$ with ${\bf R}_{p}^{G}$ and $\widehat{\bf W}_{p}^{G}$, respectively. Actually, (11) and (12) are equivalent signal models. Meanwhile, it is worth noticing that the channel matrix ${\bf X}_{p}$ and ${\bf W}_{p}$ exhibit two different forms of structured sparsity. Hence, the difference between the spatial domain and the angular domain schemes mainly lies in the different update rules (i.e., (38) and (42)) for refining the sparsity ratio $\gamma_{p,k,m}$ in line 8 of Algorithm 1. Different from (11), both the space-frequency and the angular-frequency structured sparsity of the channel matrix are considered in (12), as has been discussed in Section II-C and illustrated in Fig. 2(b). Define the neighbors of $w_{p,k,m}$ as

[TABLE]

Due to the clustered sparsity and the structured sparsity described in (8)-(10), $w_{p,k,m}$ and the elements of ${\widetilde{\cal N}}_{p,k,m}$ tend to be simultaneously either zero or non-zero444Here, we assume that no power leakage caused by the discrete Fourier transformation in (7) [29]. Otherwise, the channel would be approximately sparse.. Hence, when the GMMV-AMP algorithm is applied to (12), ${\cal N}_{p,k,m}$ in line 8 is replaced by ${\widetilde{\cal N}}_{p,k,m}$. Based on the estimate of $\left\{{\bf W}_{p}\right\}_{p=1}^{P}$, the $\{\widehat{\bf X}_{p}\}_{p=1}^{P}$ can be acquired according to (7). Moveover, the active UEs can be detected via CG-AD in (39) based on $\{\widehat{\bf X}_{p}\}_{p=1}^{P}$ or via BI-AD in (41) based on the belief indicators of $w_{p,k,m}$, $\forall p,k,m$. However, as the common sparsity across multiple BS antennas is destroyed, it is challenging to find a suitable $p_{\rm bi}$ for BI-AD. Hence, we set $p_{\rm bi}$ to $p_{\rm bi}^{\rm ang}=p_{\rm bi}^{\rm spa}{S_{\rm min}^{a}}/M$, where $S_{\rm min}^{a}$ is the minimum number of non-zero elements of $\left[{\bf W}_{p}\right]_{k,:}$, $\forall p,k$. For AUD, our simulations in Section V reveal that, for Scheme 1 based on model (11), BI-AD is more reliable than CG-AD, but for Scheme 2 based on model (12), CG-AD is more reliable than BI-AD .

III-C Alternating AUD and CE Schemes

The simultaneous AUD and CE based on (11) or (12) can not fully exploit the enhanced sparsity of ${\bf W}_{p}$ and the common sparsity pattern across the multiple columns of ${\bf X}_{p}$. Hence, we propose a Turbo-GMMV-AMP algorithm that performs AUD and CE in an alternating manner, where model (11) and model (12) are alternately exploited for enhanced performance. This facilitates the development of a CS-based adaptive AUD and CE scheme for practical massive access scenarios with unknown channel sparsity level.

1) Turbo-GMMV-AMP Algorithm (Scheme 3): This algorithm is summarized in Algorithm 2, and consists of module A and module B, as illustrated in Fig. 4. In the first turbo iteration ($j=1$), module A determines a rough estimate of the AUS, where the GMMV-AMP algorithm is applied based on model (11) to obtain the belief indicators $\pi_{p,k,m}^{j}$, $\forall p,k,m$. Even for small $G$ with poor estimation quality of ${\bf X}_{p}$, i.e., $v_{k,m}^{\infty}$ is smaller than a specific value rather than $v_{k,m}^{\infty}\!\to\!0$, the $0\leq\pi_{p,k,m}^{j}\leq 1$ tend to be 0 for $x_{p,k,m}=0$ but 1 for $x_{p,k,m}\neq 0$ after convergence of the GMMV-AMP algorithm, the related proof is similar to Appendix B. Hence, we use the BI-AD rather than CG-AD to acquire two AUS estimates with different reliability: a rough AUS $\Omega$ with a lower threshold ${\varepsilon_{\rm{det}}}$ and a reliable AUS $\Xi^{j}$ with a higher threshold ${\varepsilon_{\rm{rel}}}$, as shown in lines 6-13 in Algorithm 2, so that $\Xi^{j}\subseteq\Omega$. Here, we set $\varepsilon_{\rm det}$ to 0.4 to obtain a low missed detection probability, and $\varepsilon_{\rm rel}$ to 0.9 to avoid false alarm and to achieve reliable detection of the active UEs, i.e., the UEs in $\Xi^{j}$ are active with high probability. These two AUSs, $\Omega$ and $\Xi^{j}$, are passed on to module B.

In module B, with the rough AUS estimate $\Omega$, the angular domain channel vectors of the UEs in $\Omega$, i.e., $\left[{\bf W}_{p}\right]_{\Omega,:}$, are estimated based on the model in (12) as

[TABLE]

where $\left[{\bf S}_{p}^{G}\right]_{:,\Omega}\in\mathbb{C}^{G\times\left|\Omega\right|_{c}}$ and $\left[{\bf W}_{p}\right]_{\Omega,:}\in\mathbb{C}^{\left|\Omega\right|_{c}\times M}$ are sub-matrices of ${\bf S}_{p}^{G}$ and ${\bf W}_{p}$, respectively, ${\overline{\bf N}_{p}}=\left[{\bf S}_{p}^{G}\right]_{:,{{\cal K}-\Omega}}\left[{\bf W}_{p}\right]_{{{\cal K}-\Omega},:}+{\widetilde{\bf N}_{p}}$, ${\cal K}$ is the set of all potential UEs, and ${\cal K}-\Omega$ denotes the difference set of sets $\cal K$ and $\Omega$. To reduce the power of ${\overline{\bf N}_{p}}$, ${\cal A}\subseteq\Omega$ is desirable, i.e., a low missed detection probability. The dimension of the uplink channel matrix for CE is reduced by considering only the UEs in $\Omega$. Furthermore, the low-dimensional channel matrix $\left[{\bf W}_{p}\right]_{\Omega,:}$ is still sparse due to the angular-frequency structured sparsity of massive MIMO channels. Hence, we can estimate $\left[{\bf W}_{p}\right]_{\Omega,:},\forall p$, by applying the GMMV-AMP algorithm to (43), see line 16 of Algorithm 2. Moreover, the signals received from the UEs in $\Gamma$, a subset of $\Xi^{j}$, are removed to enhance the sparsity of the uplink massive access channel matrix for AUD. The residual received signals $\left(\widetilde{\bf Y}_{p}^{G}\right)^{j}$ are computed in lines 17 and 18, and are passed on to module A.

In the following turbo iterations ($j>1$), the AUD problem in module A is to recover $\left({\bf X}_{p}^{\rm re}\right)^{j}$ based on the following model

[TABLE]

where $(\widetilde{\bf Y}_{p}^{G})^{j}$ contains the residual received signals in the $j$-th turbo iteration, $\left({\bf X}_{p}^{\rm re}\right)^{j}={\bf X}_{p}-{\widetilde{\bf X}_{p}^{j}}$, and ${\widetilde{\bf X}_{p}^{j}}\in\mathbb{C}^{K\times M}$ is defined as, $[{\widetilde{\bf X}_{p}}^{j}]_{\Gamma,:}=[\widehat{\bf X}_{p}^{j-1}]_{\Gamma,:}$, while $[{\widetilde{\bf X}_{p}^{j}}]_{{\cal K}-\Gamma,:}={\bf 0}_{\left|{\cal K}-\Gamma\right|_{c}\times M}$. To prevent the GMMV-AMP algorithm from diverging, we only remove the signals received from a part of the UEs in $\Xi^{j}$, i.e., $\lambda_{\rm aus}<1$ (e.g. $\lambda_{\rm aus}=0.8$). Modules A and B will be executed iteratively. Since the $\left({\bf X}_{p}^{\rm re}\right)^{j}$ become sparser and the channels of the UEs in $\Omega$ are iteratively re-estimated as the turbo iterations proceed, the $\Omega$ and the corresponding channels are constantly refined. Therefore, compared to the simultaneous processing approaches in Section III-B, the proposed alternating approach facilitates more reliable AUD and CE with significantly smaller $G$, which means a dramatic reduction of access latency.

2) CS-Based Adaptive AUD and CE (Scheme 4): For practical systems, the UE activity and the channel environment are time-varying. As a result, the sparsity level of the uplink massive access channel matrix may change over time. If the channel matrix is relatively sparse, a small time slot overhead $G$ is sufficient to acquire accurate AUS and CSI estimates, while if the channel matrix is relatively dense, a large $G$ is required to guarantee reliable sparse signal recovery. This motivates us to propose a CS-based adaptive AUD and CE scheme, as shown in Algorithm 3, to adaptively adjust $G$ to facilitate low-latency and high-reliability AUD and CE. Algorithm 3 can be summarized as follows.

• •

Step 1: In each time slot, all active UEs transmit non-orthogonal RA pilots to the BS. The pilots are pre-designed based on DCS theory and known to the system.

• •

Step 2: The BS collects the received signals over multiple successive time slots, and the Turbo-GMMV-AMP algorithm is utilized to alternately acquire the AUS and CSI estimates, see lines 3-5 of Algorithm 3. Besides, the estimation reliability is evaluated based on a pre-specified criterion, see line 8 of Algorithm 3.

• •

Step 3: If the criterion is met, the outcome of the evaluation is informed to all UEs, and the active UEs stop transmitting pilots and start to transmit their data to the BS without scheduling permission. Otherwise, Step 1 is repeated until the received signals collected at the BS are sufficient to meet the evaluation criterion.

Here, based on CS theory, $G_{0}\geq{\mathbb{E}}[|{\rm supp}\{\left[{\bf W}_{p}\right]_{:,m}\}|_{c}]$ is desirable, and $\epsilon=0.8$ is suggested in [27, Sec. V].

III-D Computational Complexity Analysis

For the practical implementation of the algorithms, the related computational complexity determines the hardware cost and the power consumption for processing. Hence, the complexity analysis for the proposed algorithms is an important issue, especially for massive connectivity with large scale systems.

Table I compares the complexity of the proposed GMMV-AMP algorithm and Turbo-GMMV-AMP algorithm, as well as the conventional greedy CS recovery algorithms, i.e., generalized subspace pursuit (GSP) [43], simultaneous orthogonal matching pursuit (SOMP) [42], and distributed sparsity adaptive matching pursuit DSAMP [27, 44], in terms of the number of the required complex multiplications in each iteration for AUD and CE. Obviously the matrix inversion implemented in three greedy CS recovery algorithms for least square operation contributes to most of the computational complexity. Hence, three greedy methods have the same order of computational complexity, i.e., in order of cubic magnitude of the number of active UEs $K_{a}$. By contrast, the complexity of the proposed GMMV-AMP and Turbo-GMMV-AMP algorithms increases linearly with $K,G,M$, and $P$. Hence, for massive access scenarios with large $K_{a}$, the proposed algorithms can be more computationally efficient.

IV State Evolution

SE is a framework for analyzing the performance of AMP algorithms in the large system limit where $K\rightarrow\infty$ [41]. In this section, we harness SE to characterize the mean square error (MSE) performance of the proposed GMMV-AMP algorithm. The MSE of the estimation and the variance of the estimated channels are defined as

[TABLE]

respectively. Based on the derivations in Appendix A, the GMMV-AMP algorithm can be explained intuitively. For the $p$-th pilot subcarrier, the proposed GMMV-AMP algorithm decouples the matrix estimation problem in (11) into $KM$ independent scalar estimation problems, as

[TABLE]

where index $p$ and superscript $G$ are omitted for notational simplicity, $C_{k,m}^{q}\sim{\cal CN}\left(C_{k,m}^{q};x_{k,m},D_{k,m}^{q}\right)$ is the equivalent measurement of $x_{k,m}$ in the $q$-th iteration, and ${\tilde{n}_{k,m}^{q}}\sim{\cal CN}\left({\tilde{n}_{k,m}^{q}};0,D_{k,m}^{q}\right)$ denotes the effective noise.

Proposition 3

Define a scalar random variable ${X_{0}}\!\sim\!p_{0}\left(X\right)$. Then, we have $C^{q}=x_{0}+{\tilde{n}^{q}}$, ${\tilde{n}^{q}}\sim{\cal CN}\left({\tilde{n}^{q}};0,D^{q}\right)$, and the posterior distribution of $x_{0}$ can be expressed as

[TABLE]

where

[TABLE]

and $z\!\sim\!{\cal CN}\left(z;0,1\right)$. Hence, $e^{q+1}$ and $\vartheta^{q+1}$ are updated as

[TABLE]

where ${\cal D}x_{0}=p_{0}\left(x_{0}\right)d{x_{0}}$ and ${\cal D}z=e^{-\left|z\right|^{2}}/\pi dz$. Therefore, defining a scalar random variable following the same distribution as the channel coefficient, i.e., $X_{0}\sim p\left(X\right)$, $C_{k,m}^{q}$ and $D_{k,m}^{q}$ in the GMMV-AMP algorithm can be calculated as in (49), and $e^{q+1}$ and $\vartheta^{q+1}$ can be obtained from (50) and (51).

Proof:

Please refer to the Appendix C. ∎

Since the a priori distribution $p_{0}\left(X\right)$ does not take the structured sparsity of channel matrix into consideration, the scalar SE in (49)-(51) can not accurately analyze the MSE performance of the proposed GMMV-AMP algorithm. Hence, we use Monte Carlo simulation to carry out the SE, so that (50) and (51) are simplified and the exploitation of the structured sparsity is also taken into into account. In contrast to the conventional AMP algorithms in [23] and [24], which assume full knowledge of the a priori distribution of the channels and the noise variance, the SE for the proposed GMMV-AMP algorithm also needs to track the update rules of the hyper-parameters ${\bm{\theta}}$ in (31)-(33), and [29]

[TABLE]

where $\sigma_{0}$ is the actual noise variance. The SE of the proposed GMMV-AMP algorithm is summarized in Algorithm 4.

V Simulation Results

For the presented simulation results, we assume the BS employs a ULA with $M$ antennas, $K=500$ potential UEs are randomly distributed in the cell with radius 1 km, and $K_{a}=50$ ($K_{a}\ll K$) UEs are active unless otherwise specified. Furthermore, the carrier frequency is 2 ${\rm GHz}$, the bandwidth is $B_{s}=10$ ${\rm MHz}$, and the received SNR is 30 ${\rm dB}$. The system adopts OFDM for massive access in an eMBB scenario, where $N=2048$ subcarriers and a cyclic prefix of length $N_{\rm CP}=64$ are employed. $P=N_{\rm CP}$ pilots are uniformly allocated to the $N$ subcarriers [27]. We consider the one-ring channel model with limited angular spread [28], so that each UE’s massive MIMO channel exhibits clustered sparsity in the virtual angular domain. The large scale fading $\rho_{k}$ is given by the standard Log-distance path loss model as $\rho_{k}=128.1+37.6{\rm lg}(d_{k})$ with distance $d_{k}$ measured in km. The small scale fading channel is generated by (3), where $L$ varies from 8 to 40 [29], the related AOAs are generated within an angular spread $\Delta$ varying from [math] to [math] so that the effective sparsity level in virtual angle domain $S_{a}$ varies from 8 to 14, ${\beta_{k,l}}\sim{\cal CN}({\beta_{k,l}};0,1)$, and $\varpi_{k,l}$ is randomly and uniformly selected from $[0,N_{CP}/B_{s}]$. Furthermore, $T_{\rm amp}=200$, $T_{\rm tur}=10$, $\eta=10^{-5}$, and the simulation results are obtained by averaging over $N_{\rm sim}=3\times 10^{3}$ simulation runs unless otherwise specified.

An overview of the proposed Schemes 1-4 and the considered baseline schemes, GSP [43], SOMP [42], and DSAMP [27], is presented in Table II. For performance evaluation, we consider the detection error probability $P_{e}$ for AUD and the MSE for CE, which are respectively defined as

[TABLE]

Here, to reduce computational complexity, for Schemes 1-4, only $\widetilde{P}$ out of $P$ pilot subcarriers are used to estimate AUS $\cal A$ and the corresponding channels. The remaining $(P-\widetilde{P})$ pilot subchannels of the active UEs can be easily estimated by applying the GMMV-AMP algorithm to (43) given $\Omega={\widehat{\cal A}}$.

Fig. 5 verifies the superiority of the proposed DCS-based pilot design for broadband massive access based on Scheme 1. The quality of the pilots is evaluated in terms of the success rate, which is defined as the ratio of the number of simulation runs with $P_{e}=0$ to the total number of simulation runs. Fig. 5 shows that, as expected, employing different ${\bf S}_{p}^{G}$ for different $p$ improves the success rate. Moreover, the performance is further improved when massive MIMO and larger $\widetilde{P}$ are employed, since the structured sparsity of $\left\{{\bf X}_{p}\right\}_{p=1}^{P}$ is leveraged.

Fig. 6 examines the AUD performance of Scheme 1, where the performances of three state-of-the-art GMMV-CS algorithms are shown as benchmarks. As can be observed, the GMMV-AMP algorithm outperforms the other three algorithms. Hence, the access latency can be considerably reduced for a given target $P_{e}$. For example, for $M=16$ and ${\widetilde{P}}=1$, the GSP-based scheme requires $G=72$ to achieve $P_{e}\leq 10^{-5}$, whereas Scheme 1 with BI-AD needs only $G=58$, which indicates a reduction of approximately $19\%$ in the access latency. Moreover, Scheme 1 can achieve a better detection performance by equipping more antennas at the BS and/or utilizing larger $\widetilde{P}$, since a larger $M$ and/or $\widetilde{P}$ can enhance the space-frequency structured sparsity shown in Fig. 2(a). However, this improvement becomes negligible when $M$ and $\widetilde{P}$ are sufficiently large.

In Fig. 6, we further compare the performance of CG-AD and BI-AD. For $G<55$ or ${\widetilde{P}}\geq 8$, CG-AD and BI-AD have similar performance. However, when $G\geq 56$ and $\widetilde{P}\leq 8$, BI-AD outperforms CG-AD, as CG-AD suffers from a detection error floor. The reason for this behavior is that when there are sufficiently many measurements for reliable CS recovery, the belief indicator $\pi_{p,k,m}$ takes values of 0 and 1, but the estimated channel gain ${\hat{x}_{p,k,m}}$ takes the true value of $x_{p,k,m}$. Hence, based on the threshold function $r\left(x;\varepsilon\right)$, BI-AD can reliably determine whether $x_{p,k,m}$ is zero or not, while for CG-AD, missed detections and false alarms can not be avoided, which may lead to a detection error floor. Clearly, for Scheme 1, BI-AD is more reliable than CG-AD for AUD.

Fig. 7 depicts the CE MSE performance of Scheme 1 and the three state-of-the-art GMMV-CS algorithms also considered in Fig. 6. The oracle least square (LS) estimator with known AUS $\cal A$ is used as performance upper bound [27]. When $G$ is sufficiently large, both Scheme 1 and the three baseline schemes approach the oracle LS performance bound, since ${\cal A}$ is accurately estimated in this case, and the CE problem is reduced to an oracle LS problem. However, for $G<70$, Scheme 1 outperforms the three baseline algorithms, and its performance improves when $M$ and/or $\widetilde{P}$ increase. Besides, the MSE performance of Scheme 1 is accurately predicted by SE. Here, an important observation is that when $G<K_{a}$, both Scheme 1 and the oracle LS estimator can not perform reliable CE. This suggests that $G\geq K_{a}$ is required for reliable CE in (11). Hence, the reduction of $G$ is limited to $K_{a}$, which is identical to the sparsity level of the column vectors in ${\bf X}_{p}$. This motivates Scheme 2 for simultaneous AUD and CE, where the sparsity level of the channel matrix, defined as ${\widetilde{K}_{a}}={\rm max}\left\{|{\rm supp}\{\left[{\bf W}_{p}\right]_{:,m}\}|_{c},\forall p,m\right\}$, is less than $K_{a}$.

Fig. 8 compares the AUD performance of Scheme 1 and Scheme 2. For $M=16$, Scheme 1 outperforms Scheme 2 for both BI-AD and CG-AD, respectively. In contrast, for $M=32$, Scheme 2 achieves a much better AUD performance than Scheme 1 when $G>28$. This is because ${\bf W}_{p}$ is sparser than ${\bf X}_{p}$, i.e., ${\widetilde{K}_{a}}<K_{a}$, and the required $G$ for reliable AUD and CE in Scheme 2 mainly depends on ${\widetilde{K}_{a}}$ rather than $K_{a}$. However, the virtual angular domain sparsity weakens the common sparsity pattern across multiple columns of ${\bf W}_{p}$. Therefore, if the BS has a small number of antennas, e.g., $M=16$, Scheme 1 outperforms Scheme 2 by leveraging the common sparsity observed at different BS antennas. However, when $M$ becomes large, Scheme 2 can considerably reduce the required $G$ for reliable AUD and CE compared to Scheme 1. In addition, BI-AD in Scheme 2 suffers from an obvious detection error floor, which indicates that for Scheme 2, CG-AD is more reliable than BI-AD for detection of the active UEs. Fig. 9 compares the CE MSE performance of Scheme 1 and Scheme 2, which again verifies the superiority of Scheme 2 for massive MIMO systems. The theoretical SE accurately predicts the MSE.

Fig. 10 and Fig. 11 compare the MSE and $P_{e}$ performance of Schemes 1-4, respectively. The details of Scheme 4 are shown in Algorithm 3. Given $\gamma=K_{a}/K=0.1$ and $M=64$, an initial time slot overhead of

[TABLE]

is adopted. For different simulation runs, the varying $S_{a}$ yields a different sparsity level for $\left\{{\bf W}_{p}\right\}_{p=1}^{P}$, i.e., different ${\widetilde{K}_{a}}$. The results and the consumed $G$, i.e., the overhead, are recorded after the pre-defined criterion in line 8 of Algorithm 3 is met. For very low time slot overheads, i.e., $G\in\left[10,18\right]$, both Scheme 1 and Scheme 2 have a poor performance. This is because for these schemes, $G\ll K_{a}$ leads to extremely insufficient measurements. By contrast, Scheme 3 using the proposed Turbo-GMMV-AMP algorithm can achieve a much better AUD and CE performance than Scheme 1 and Scheme 2, which confirms its superiority in reducing access latency. Finally, Scheme 4, i.e., the proposed CS-based adaptive AUD and CE scheme, adaptively adjusts the overhead to achieve satisfactory AUD and CE performance. In Fig. 11, for Scheme 4, the percentage of simulation runs requiring the given time slot overhead $G$ is provided. This reveals that $88.6\%$ of the simulation runs require overheads of $G\in\left[12,18\right]$, and the corresponding MSE performance is much better than those for Schemes 1-3.

Fig. 12 further compares the AUD and CE performance of Scheme 3 and Scheme 4 for different number of active UEs $K_{a}$. The proposed Scheme 4 adaptively adjusts the time slot overhead to guarantee reliable AUD and CE for different $K_{a}$. However, Scheme 3 employs a fixed overhead and suffers from poor performance when $K_{a}$ becomes large. This means that some of the active UEs will not be able to access the network. Hence, for practical systems with time-varying UE activity, the superiority of the proposed CS-based adaptive AUD and CE scheme is evident.

VI Conclusion

This paper investigates new methods for facilitating massive access in massive MIMO systems, which leverage the sporadic traffic of the UEs and the virtual angular domain sparsity of massive MIMO channels to dramatically reduce the access latency. The space-frequency structured sparsity of the channel matrix in the spatial domain improves the AUD performance, while the angular-frequency structured sparsity of the channel matrix in the angular domain improves the CE performance. Therefore, joint AUD and CE schemes exploiting only spatial domain or only angular domain channel model can not take full advantage of the sparsity properties of massive access in massive MIMO systems. This motivates the derivation of the proposed Turbo-GMMV-AMP algorithm, which achieves a significant performance improvement by performing AUD based on a spatial domain channel model and CE based on an angular domain channel model in an alternating manner. Furthermore, for practical systems, where the number of active UEs is not known, the proposed CS-based adaptive AUD and CE scheme can adjust the time slot overhead to realize ultra-reliable low-latency massive access.

Appendix A Proof of the Proposition 1

The factorization in the joint posterior probability (15) can be represented by a bipartite graph, which motivates the application of the sum-product (SP) algorithm to realize the MMSE estimator [37]. As the bipartite graph consists of $M$ independent subgraphs, we only discuss the $m$-th subgraph in the following derivations, and the antenna index $m$ is dropped for notational simplicity. For the $m$-th subgraph, we define variable nodes $V=\left[K\right]$, factor nodes of the likelihood function $F=\left[G\right]$, and edges $E=\left[K\right]\times\left[G\right]=\left\{\left(k,g\right):k\in\left[K\right],g\in\left[G\right]\right\}$. The update rules for the messages associated to the edges are [33]

[TABLE]

where $b\in\left[G\right]$, $j\in\left[K\right]$, $q$ denotes the $q$-th iteration, and $\propto$ denotes equality up to a constant scale factor.

One practical hurdle for the large-scale implementation of the SP algorithm lies in the required evaluation of high-dimensional integrals for calculation of messages ${\hat{\xi}_{g\to k}^{q}}\left(x_{k}\right)$. This leads to an unacceptably high complexity. However, a key observation is that, in the large system limit with $K\to\infty$, the messages ${\hat{\xi}_{g\to k}^{q}}\left(x_{k}\right)$ can be approximated by Gaussian distributions [33]. Since the random variables $\left\{x_{j}\right\}_{j=1}^{K}$ are independently complex Gaussian distributed, random variable $Z_{g,k}=\sum\nolimits_{j\neq k}s_{g,j}x_{j}$ follows a complex Gaussian distribution $Z_{g,k}\sim{\cal CN}\left(z_{g,k};Z_{g\to k}^{q},V_{g\to k}^{q}\right)$, with mean $Z_{g\to k}^{q}$ and variance $V_{g\to k}^{q}$ given as

[TABLE]

where ${\hat{x}_{j\to g}^{q}}$ and $v_{j\to g}^{q}$ are the mean and variance of message ${\xi_{j\to g}^{q}}\left(x_{j}\right)$, respectively. Hence, the messages ${\hat{\xi}_{g\to k}^{q}}\left(x_{k}\right)$ can be approximated as

[TABLE]

Further, the posterior distribution of $x_{k}$ is calculated as

[TABLE]

where

[TABLE]

It is convenient to introduce a family of densities

[TABLE]

where ${\widetilde{Z}_{4}}=\int{p_{0}\left(x\right){\cal CN}\left(x;C,D\right)}dx$ is a normalization constant. The corresponding mean and variance are

[TABLE]

respectively. Define ${\xi_{k}^{q+1}}(x_{k})\propto{\cal CN}\left(x_{k};{\hat{x}_{k}^{q+1}},v_{k}^{q+1}\right)$, with ${\hat{x}_{k}^{q+1}}=g_{a}\left(C_{k}^{q},D_{k}^{q}\right)$ and $v_{k}^{q+1}=g_{c}\left(C_{k}^{q},D_{k}^{q}\right)$. The messages ${\xi_{k\to g}^{q+1}}\left(x_{k}\right)$ can be approximated as

[TABLE]

where

[TABLE]

At this point, the messages ${\xi_{k\to g}^{q+1}}\left(x_{k}\right)$ and ${\hat{\xi}_{g\to k}^{q}}\left(x_{k}\right)$ have been approximated as Gaussian densities. However, the computational complexity is still high when the system is large, since the number of messages scales with the number of potential UEs $K$. In order to reduce the number of messages in the $q$-th iteration, we can further simplify the update rules by making some approximations. Defining $Z_{g}^{q}=\sum\nolimits_{j=1}^{K}{s_{g,j}{\hat{x}_{j\to g}^{q}}}$, $V_{g}^{q}=\sum\nolimits_{j=1}^{K}{\left|s_{g,j}\right|^{2}v_{j\to g}^{q}}$, we can rewrite (57) as

[TABLE]

Substituting (67) into (64) and ignoring terms approximated as 0 in the large system limit $K\to\infty$ [33], the update rules at the variable nodes and the factor nodes can be approximated as in (19)-(22). Finally, the posterior distribution of $x_{k}$ is given by (18).

Appendix B Proof of the Proposition 2

If a reliable estimate of ${\bf X}_{p}$ is acquired after the convergence of the GMMV-AMP algorithm, the variance of the posterior distribution of $x_{k,m}$ tends to be zero, i.e., $v_{k,m}^{\infty}\to 0$, thus,

[TABLE]

Therefore, $D_{k,m}^{\infty}$ in (19) is calculated as

[TABLE]

Here, approximation $(a)$ is because the pilots are generated from an i.i.d standard complex Gaussian distribution, i.e., $s_{g,k}\sim{\cal CN}\left(s_{g,k};0,1\right)$, thus $\sum\nolimits_{g}{\left|s_{g,k}\right|^{2}}\approx\sum\nolimits_{g}{\mathbb{E}}\left[\left|s_{g,k}\right|^{2}\right]=G$. In the large system limit, as $G\to\infty$, $D_{k,m}^{\infty}\to 0^{+}$. For a given realization of the massive access channel matrix ${\bf X}_{p}$, defining $r^{\infty}=\sum\nolimits_{g}{\frac{s_{g,k}^{*}\left(y_{g,m}-Z_{g,m}^{\infty}\right)}{\sigma}+V_{g,m}^{\infty}}$, $C_{k,m}^{\infty}$ in (20) is given as

[TABLE]

Hence, substituting (69) and (70) into (26), for $x_{k,m}=0$, we have

[TABLE]

while for $x_{k,m}\neq 0$,

[TABLE]

which yields

[TABLE]

Appendix C Proof of the Proposition 3

Here, for simplicity of derivation, we focus on the $m$-th subgraph only and drop the antenna index $m$, as in Appendix A. The derivation is based on (11), and can be easily extended to model (12), thus we have

[TABLE]

Substituting (57) and (74) into (60), $C_{k}^{q}$ is computed as

[TABLE]

In (75), approximation $(b)$ is because the term $\sigma+\sum\nolimits_{j\neq k}{\left|s_{g,j}\right|^{2}}v_{j\to g}^{q}$ is approximately independent of $g$ in the large system limit with $K\to\infty$, which has been proven in [33]. Define

[TABLE]

Since $s_{g,k}\sim{\cal CN}\left(s_{g,k};0,1\right)$, $n_{g}\sim{\cal CN}\left(n_{g};0,\sigma\right)$, and $s_{g,k}$ is independent of $n_{g}$, $r_{k}$ follows a complex Gaussian distribution according to the central limit theorem as $K\to\infty$. Moreover, by substituting (45) into (76), we find the mean and the variance of $r_{k}^{q}$ are 0 and $G(\sigma+Ke^{q})$, respectively. Hence,

[TABLE]

where $z\sim{\cal CN}\left(z;0,1\right)$. Meanwhile, by substituting (46) into $D_{k}^{q}$ in (60), we can obtain

[TABLE]

Hence, the Proposition 2 is proven.