Message Passing-Based Joint User Activity Detection and Channel Estimation for Temporally-Correlated Massive Access

To accomplish joint user activity detection and channel estimation for massive access, a variety of advanced sparse signal recovery algorithms have been proposed for different wireless communication systems. For the single-cell scenario, the works [9, 10] utilize the standard CS algorithms of orthogonal matching pursuit (OMP) and basis pursuit denoising (BPDN). However, these two kinds of algorithms usually suffer high computational complexity due to the extremely large amount of devices. Then the works [11, 12] propose the computationally efficient AMP algorithms with Bayesian denoiser for user activity detection. In the work [13], the message-scheduling generalized AMP algorithm is developed to further reduce the computational complexity without degrading the detection and estimation performance. The AMP-based algorithms are also extended to perform data detection along with the activity detection in [14, 15]. Besides the AMP-based algorithms, covariance matching pursuit [16, 17], dimension reduction-based optimization [18] and deep learning methods [19, 20, 21] are also investigated for performance enhancement. In particular, the covariance matching pursuit algorithm [16, 17] and dimension reduced-based optimization algorithm [18] are especially designed for the massive multiple-input multiple-out (MIMO) systems, which can significantly outperform the conventional CS-based methods. The deep learning methods [19, 20, 21] employ the artificial neural networks and learn the network parameters based on the pre-collected data, which are able to outperform the traditional hand-designed algorithms and to enjoy low computational complexity. Moreover, a sparsity-constrained method is proposed for the non-ideal scenario where different users have unknown frequency offsets in [22]. In the multi-cell systems, the cooperative user activity detection and channel estimation is studied in [23, 24] based on the AMP-based algorithms. Recently, unsourced random access as a new random access protocol is also studied in [25, 16], where all users share a common codebook and the BS only needs to detect the transmit codewords.

The aforementioned works all detect the user activity in each frame individually by assuming the user activity is independent for each frame. Due to the low-rate transmission, the information transmission of one active user may occupy multiple consecutive frames. Such temporal correlation of user activity can be exploited for performance enhancement. Assuming the user sparsity level and the channel state information (CSI) are available, the work [26] proposes a DCS-based algorithm where the set of the detected active user in the last slot is used as the initial set to be detected in the current slot. By adaptively exploiting the prior support based on the corresponding support quality information, the work [27] proposes a prior-information-aided adaptive subspace pursuit (PIA-ASP) algorithm, which can always outperform the DCS-based algorithm in [26]. These algorithms include the matrix inversion calculation in each iteration, and thus suffer high computational complexity in the large-scale problems. Moreover, both the works [26, 27] assume that the CSI is perfectly known in advance, which is unrealistic in massive access. By leveraging the system statistics, the authors in [28] propose a sequential AMP (S-AMP) algorithm to sequentially perform user activity detection and channel estimation in each frame. Then the work [29] proposes to extract the side information (SI) from the estimation results in the previous frame to enhance the user activity detection performance in the current frame based on the AMP framework. Different from [29] which only utilizes the historical estimation results as SI, our previous work [30] proposes to extract the double-sided information by considering the estimation from the next frame as well for further exploiting the temporal correlation.

Compared with these prior works [26, 27, 29, 28, 30] that exploit the temporal correlation for use activity detection and channel estimation, this paper differs mainly in the following two aspects. First, while all the exiting works perform activity detection in a frame-by-frame manner, (i.e., the user activity is still detected sequentially frame by frame, though the temporal correlation among adjacent frames has been exploited), this paper performs multi-frame activity detection in a block-by-block manner. Second, this paper proposes to adaptively learn the parameters of the system statistics by using the EM technique, while the previous algorithms either do not consider the system statistics [26, 27] or assume the system statistics are perfectly known [29, 28, 30].

In this work, we consider the problem of joint user activity detection and channel estimation in multiple consecutive frames for the temporally-correlated massive access systems. The main contributions and distinctions of this work are summarized as follows:

1. 1.

   We propose to perform the user activity detection and channel estimation jointly in multiple consecutive frames and formulate it as a DCS problem. Based on the probabilistic model, not only the sparse activity pattern is considered in the problem, but also the statistical relationships of the activities in these consecutive frames are exploited.

2. 2.

   In contrast to only making use of the historical estimations in [28, 29], this paper proposes a HyGAMP-DCS algorithm to fully exploit the temporally-correlated user activities in the multiple consecutive frames. The HyGAMP-DCS algorithm combines the computationally efficient GAMP algorithm for channel estimation and the standard MP algorithm for the activity likelihood update. In particular, the activity likelihood values in each frame are updated by aggregating both the forward messages from the previous frame and backward messages from the next frame at each iteration. Numerical results show that HyGAMP-DCS can achieve superior performance to the existing DCS-based algorithms [26, 27, 28, 29] thanks to the bidirectional message propagation.

3. 3.

   To make the proposed algorithm applicable for the practical system with imperfect system statistics, this paper incorporates the EM algorithm in HyGAMP-DCS to adaptively learn the hyperparameters during the estimation procedure. Compared with the traditional EM-AMP algorithm [31] for the CS problem, the proposed EM-HyGAMP-DCS algorithm additionally learns the statistical dependencies between the activities in the neighboring frames. Simulation results demonstrate that EM-HyGAMP-DCS realizes very similar performance to HyGAMP-DCS which requires the perfect system statistics.

4. 4.

   We provide the performance and complexity analysis of the proposed algorithms. We first introduce the state evolution (SE) to predict the asymptotic performance of HyGAMP-DCS and EM-HyGAMP-DCS. In particular, we also point out that the SE is quite essential to find the appropriate hyperparameter initialization for EM-HyGAMP-DCS, which is validated by the numerical results. Then the computational complexity comparison with the state-of-the-art methods is given to illustrate the computational efficiency of our proposed algorithms.

The remaining part of this paper is organized as follows. Section II introduces the system model of temporally-correlated massive access. In Section III, we present the probabilistic model of the considered system and then introduce the Bayesian inference methods. In Section IV, the HyGAMP-DCS algorithm is proposed, then the EM algorithm is employed for hyperparameter update in Section V. Next, we analyze the performance and the computational complexity of the proposed algorithms in Section VI. The simulation results of the proposed algorithms are shown in Section VII. Finally, we conclude this paper in Section VIII.

In this paper, upper-case and lower-case letters denote random variables and their realizations, respectively. Letters $\mathbf{x}$, $\mathbf{X}$ denote vector and matrix, respectively. Superscripts $(\cdot)^{T},(\cdot)^{*}$ denote transpose and conjugate, respectively. Further, $\mathbb{E}[\cdot]$ and $\mathbb{V}[\cdot]$ denotes the expectation operation and variance operation, respectively; $|\cdot|$ denotes the magnitude of a variable. In addition, a random vector $\mathbf{x}\in\mathbb{C}^{M\times 1}$ drawn from the complex Gaussian distribution with mean $\mathbf{x}_{0}\in\mathbb{C}^{M\times 1}$ and covariance matrix $\boldsymbol{\Sigma}\in\mathbb{C}^{M\times M}$ is characterized by the probability density function $\mathcal{CN}(\mathbf{x};\mathbf{x}_{0},\boldsymbol{\Sigma})=\frac{\exp\left(-(\mathbf{x}-\mathbf{x}_{0})^{H}\boldsymbol{\Sigma}^{-1}(\mathbf{x}-\mathbf{x}_{0})\right)}{\pi^{M}|\boldsymbol{\Sigma}|}$.

As shown in Fig. 1, once a user is activated, its data transmission can occupy multiple consecutive frames. As such, there exist temporal correlations of the user activities in neighboring frames.

Similar to [28, 29, 30, 6], the dynamic user activity over frames is modeled as a stationary first-order Markov chain with two discrete states. The Markov chains on different devices are assumed to independent and identically distributed (i.i.d.). Let $\lambda_{n,t}\in\{0,1\}$ indicate whether user $n\in\{1,2,\dots,N\}$ is active or not in the $t$th frame. The steady state probability is given by $\text{Pr}(\lambda_{n,t}=1)=p_{a}$, where $p_{a}$ denotes the active probability of each user $n$ in frame $t$. The transition probability matrix is defined as

where $p_{10}=\text{Pr}\{\lambda_{n,t}=0|\lambda_{n,t-1}=1\}$ and other three transition probabilities are defined similarly. By the steady-state assumption, the Markov chain can be completely characterized by two parameters $p_{a}$ and $p_{10}$. The other three transition probabilities can be easily obtained as $p_{01}=p_{a}p_{10}/(1-p_{a})$, $p_{00}=1-p_{01}$, and $p_{11}=1-p_{10}$. Note that $1/p_{10}$ can be considered as the expected number of the consecutive frames occupied by the data transmission of one user, which means that smaller $p_{10}$ leads to stronger temporal correlation of the user activity over frames. Although the temporal correlation can be described more precisely by higher-order Markov chains, we consider the first-order Markov chain here for simplicity.

In order to exploit the temporally correlated user activities, we consider $T$ consecutive frames for signal transmission with $T\geq 2$. In each frame $t\in\{1,2,\dots,T\}$, the transmission is divided into two phases, the pilot phase and data phase. Each user $n$ is pre-assigned a unique pilot sequence $\mathbf{a}_{n}=[a_{n,1},a_{n,2},\dots,a_{n,L}]^{T}\in\mathbb{C}^{L\times 1}$ with length $L$. The BS performs user identification and channel estimation based on the detected pilot sequences in the pilot phase and then decodes the data in the data phase. In this work, we concentrate on the joint user activity detection and channel estimation problem in the pilot phase. Due to the limit of orthogonal resources, the pilot length $L$ is assumed to be much smaller than the number of users $N$, i.e., $L\ll N$, so that the pilot sequences of different users cannot be mutually orthogonal. The pilot sequences in this work are generated from the i.i.d. complex Gaussian distributions with zero mean and variance of $1/L$, and then scaled to have unit power.

The channel in each frame is assumed to satisfy independent Rayleigh distribution. We denote the channel coefficient between user $n$ and the BS in the $t$th frame as $h_{n,t}=\sqrt{\beta_{n}}g_{n,t}$, where $\beta_{n}=P_{n}\gamma_{n}$ is the effective large-scale fading component decided by the transmit power $P_{n}$ and the large-scale attenuation $\gamma_{n}$, $g_{n,t}\in\mathbb{C}$ is the small-scale fading coefficient following i.i.d. circularly symmetric complex Gaussian distributions with zero mean and unit variance, i.e., $g_{n,t}\sim\mathcal{CN}(0,1)$. We adopt the power control strategy in [14], so that the effective large-scale fading component of each user is the same, i.e., $\beta_{n}=\beta,\forall n$. Then the received pilot signals at the BS in continuous $T$ frames can be written as

where $\mathbf{A}=[\mathbf{a}_{1},\mathbf{a}_{2},\dots,\mathbf{a}_{N}]\in\mathbb{C}^{L\times N}$ is the pilot matrix of all users; $\mathbf{\Lambda}=[\boldsymbol{\lambda}_{1},\boldsymbol{\lambda}_{2},\dots,\boldsymbol{\lambda}_{T}]\in\mathbb{C}^{N\times T}$ is the activity matrix with $\boldsymbol{\lambda}_{t}=[\lambda_{1,t},\lambda_{2,t},\dots,\lambda_{N,t}]^{T}\in\mathbb{C}^{N\times 1}$; $\mathbf{H}=[\mathbf{h}_{1},\mathbf{h}_{2},\dots,\mathbf{h}_{T}]\in\mathbb{C}^{N\times T}$ is the channel matrix with $\mathbf{h}_{t}=[h_{1,t},h_{2,t},\dots,h_{N,t}]^{T}\in\mathbb{C}^{N\times 1}$; $\mathbf{X}=[\mathbf{x}_{1},\mathbf{x}_{2},\dots,\mathbf{x}_{T}]\in\mathbb{C}^{N\times T}$ is the effective channel matrix defined as the Hadamard product of $\mathbf{\Lambda}$ and $\mathbf{H}$ with each entry $x_{n,t}=\lambda_{n,t}h_{n,t}$; $\mathbf{W}=[\mathbf{w}_{1},\mathbf{w}_{2},\dots,\mathbf{w}_{T}]\in\mathbb{C}^{L\times T}$ is the additive white complex Gaussian noise (AWGN) matrix where each element has zero mean and variance $\sigma^{2}_{w}$.

Our objective is to jointly detect the active users and estimate their channels by recovering $\mathbf{X}$ from the received signal $\mathbf{Y}$ given the pilot matrix $\mathbf{A}$. This underdetermined linear inverse problem can be referred to as the dynamic compressed sensing (DCS) problem, where the support of the vector $\mathbf{x}_{t}$ slowly changes with $t$. The DCS problem can be solved by optimization-based methods using convex relaxation and the standard convex problem solver, e.g., Dynamic LASSO [5]. It can also be solved by the greed-based algorithms based on OMP and SP by employing the detected active user set in the previous frame as the prior information to assist the detection in the current frame [26, 27]. These two approaches can outperform the traditional CS-based algorithms, but however, are unable to fully utilize the system statistics including the Markov chain-based user activity and the channel distributions. To improve the performance over the optimization-based methods and the greedy-based algorithms, we propose to employ the Bayesian inference method in this work, whose details will be given in the next section.

Following the HyGAMP framework, the GAMP part regards the edges between the nodes $\{x_{n,t}\}$ and $\{p(y_{l,t}|\mathbf{x}_{t})\}$ as the weak edges thanks to their linearizable coupling. Then the messages $\mu^{i}_{(n,t)\rightarrow(l,t)}(x_{n,t})$ and $\upsilon^{i}_{(n,t)\leftarrow(l,t)}(x_{n,t})$ are simplified by using GAMP approximations based on the central limit theorem and Taylor expansion.

The basic GAMP estimation is given in lines 8-17 of Algorithm 1, which can incorporate arbitrary distribution on the input $\mathbf{X}$ and output $\mathbf{Y}$ [32]. Compared with the standard MP algorithm, the number of the update variables in GAMP estimation shrinks from $\mathcal{O}(LN)$ to $\mathcal{O}(L+N)$. In lines 8-11, the GAMP estimation performs the update of the noiseless signal $\mathbf{Z}=\mathbf{A}\mathbf{X}$. The distribution of each element $z_{l,t}$ is obtained as $\mathcal{CN}(z_{l,t};\widehat{p}_{l,t}(i),\tau^{p}_{l,t}(i))$ by the messages from the variable nodes $\{x_{n,t}\}$, where $\widehat{p}_{l,t}(i)$ is the plug-in estimate of $z_{l,t}$ with variance $\tau^{p}_{l,t}(i)$ in the $i$th iteration. By accounting the AWGN channel between $\mathbf{Z}$ and $\mathbf{Y}$, the explicit expression of the MMSE estimator for $z_{l,t}$ and its variance in lines 10-11 are available as

After updating the noiseless signal $z_{n,t}$, the mean and the variance of each element in the residual signal that removes the estimated $\mathbf{X}$ from $\mathbf{Y}$ are given in lines 12-13. Then lines 14-15 give the update of the plug-in estimate $\widehat{r}_{n,t}$ of the true signal $x_{n,t}$, which is modeled as $\widehat{r}_{n,t}(i)=x_{n,t}+n_{n,t}(i)$ with $n_{n,t}(i)\sim\mathcal{CN}(0,\tau^{r}_{n,t}(i))$. To obtain the MMSE estimator of $\mathbf{X}$, we can give the approximate posterior marginal distribution of $x_{n,t}$ as

where $p^{i}_{X}(x_{n,t})=(1-\overrightarrow{p}_{n,t}(i))\delta(x_{n,t})+\overrightarrow{p}_{n,t}(i)\mathcal{CN}(x_{n,t};0,\beta)$ based on the extrinsic message $\mu^{i}_{(n,t)\rightarrow(n,t)}(\lambda_{n,t})=(1-\overrightarrow{p}_{n,t}(i))(1-\lambda_{n,t})+\overrightarrow{p}_{n,t}(i)\lambda_{n,t}$. The detailed derivation of $\overrightarrow{p}_{n,t}(i)$ will be given as (25) in the MP part. By simplifying (13), the approximate posterior marginal distribution $p^{i}(x_{n,t}|\mathbf{Y})$ can be written in the form of a spike and slab probability as

where

Thus, we can obtain the explicit expression of the posterior mean and variance of $x_{n,m}$ in lines 16-17 as

Since $\overrightarrow{p}_{n,t}(i)$ is updated by the MP part at each iteration, the MMSE estimator (18) also changes at each iteration. However, the MMSE estimator on $x_{n,t}$ keeps constant at each iteration in the AMP-based algorithms proposed in [28, 29] for temporally-correlated massive access.

By regarding all the edges in the MP part as strong edges, we utilize the standard MP algorithm to update the activity probability of each user with the extrinsic message $\upsilon^{i}_{(n,t)\leftarrow(n,t)}(\lambda_{n,t})$ from the GAMP part.

First, the extrinsic message $\upsilon^{i}_{(n,t)\leftarrow(n,t)}(\lambda_{n,t})$ is obtained as

where $\upsilon^{i}_{(n,t)\leftarrow(n,t)}(x_{n,t})=\mathcal{CN}(x_{n,t};\widehat{r}_{n,t}(i),\tau^{r}_{n,t}(i))$ and $\overleftarrow{p}_{n,t}(i)$ can be expressed by

With (IV-B), we can implement the statistical dependencies across the frames based on the Markov chain of the activity variables $\{\lambda_{n,t}\}$. We also define $\overleftarrow{q}_{n,t}(i)$ and $\overrightarrow{q}_{n,t}(i)$ as the activity likelihoods of user $n$ contained in the messages $\zeta^{i}_{(n,t)\leftarrow(n,t+1)}(\lambda_{n,t})$ and $\xi^{i}_{(n,t)\rightarrow(n,t)}(\lambda_{n,t})$, respectively. The message $\xi^{i}_{(n,t)\rightarrow(n,t+1)}(\lambda_{n,t})$ for $\forall t=1,\dots,T-1$ is expressed by the product of two Bernoulli distributions as

Then we update the forward message $\xi^{i}_{(n,t)\rightarrow(n,t)}(\lambda_{n,t})$ that represents the useful information conveyed from $\lambda_{n,t-1}$ to $\lambda_{n,t}$. For $t=1$, we always $\xi^{i}_{(n,t)\rightarrow(n,t)}(\lambda_{n,t})=(1-p_{a})(1-\lambda_{n,t})+p_{a}\lambda_{n,t}$, i.e., $\overrightarrow{q}_{n,1}(i)=p_{a},\forall n$, since $p(\lambda_{n,1})$ is only connected with $\lambda_{n,1}$ in the factor graph. For $t\geq 2$, the forward message $\xi^{i}_{(n,t)\rightarrow(n,t)}(\lambda_{n,t})$ can be obtained as