# Asynchronous Grant-Free Random Access: Receiver Design with Partially Uni-
Directional Message Passing and Interference Suppression Analysis

Zhaoji Zhang, Yuhao Chi, Qinghua Guo,
Ying Li, Guanghui Song, and Chongwen Huang Zhaoji Zhang, Yuhao Chi, Ying Li,
and Guanghui Song are with the School of Telecommunications Engineering,
Xidian University, Xi’an 710071, China (email<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>songguanghui@xidian.edu.cn).
Qinghua Guo is with the School of Electrical, Computer and Telecommunications
Engineering, University of Wollongong, Wollongong, NSW 2522, Australia
(e-mail<EMAIL_ADDRESS>Chongwen Huang is with the College of Information
Science and Electronic Engineering, Zhejiang University, Hangzhou 310027,
China, also with the International Joint Innovation Center, Zhejiang
University, Haining 314400, China, and also with the Zhejiang Provincial Key
Laboratory of Information Processing, Communication and Networking (IPCAN),
Hangzhou 310027, China (e-mail: chongwenhuang@zju.edu.cn).

###### Abstract

Massive Machine-Type Communications (mMTC) features a massive number of low-
cost user equipments (UEs) with sparse activity. Tailor-made for these
features, grant-free random access (GF-RA) serves as an efficient access
solution for mMTC. However, most existing GF-RA schemes rely on strict
synchronization, which incurs excessive coordination burden for the low-cost
UEs. In this work, we propose a receiver design for asynchronous GF-RA, and
address the joint user activity detection (UAD) and channel estimation (CE)
problem in the presence of asynchronization-induced inter-symbol interference.
Specifically, the delay profile is exploited at the receiver to distinguish
different UEs. However, a sample correlation problem in this receiver design
impedes the factorization of the joint likelihood function, which complicates
the UAD and CE problem. To address this correlation problem, we design a
partially uni-directional (PUD) factor graph representation for the joint
likelihood function. Building on this PUD factor graph, we further propose a
PUD message passing based sparse Bayesian learning (SBL) algorithm for
asynchronous UAD and CE (PUDMP-SBL-aUADCE). Our theoretical analysis shows
that the PUDMP-SBL-aUADCE algorithm exhibits higher signal-to-interference-
and-noise ratio (SINR) in the asynchronous case than in the synchronous case,
i.e., the proposed receiver design can exploit asynchronization to suppress
multi-user interference. In addition, considering potential timing error from
the low-cost UEs, we investigate the impacts of imperfect delay profile, and
reveal the advantages of adopting the SBL method in this case. Finally,
extensive simulation results are provided to demonstrate the performance of
the PUDMP-SBL-aUADCE algorithm.

###### Index Terms:

Grant-free random access, asynchronization, partially uni-directional message
passing, user activity detection, channel estimation

##  I Introduction

Massive Machine-Type Communications (mMTC) serves as a key application
scenario in the fifth-generation (5G) to support diversified Internet-of-
Things (IoT) applications [1], such as smart cities, smart home, smart
factories, etc. Unlike conventional human-oriented application scenarios, mMTC
is characterized by the following distinctive features, i.e., the massiveness
of potential user equipments (UEs), the sparse activity of each UE, short data
packets from active UEs, as well as the need for low power consumption from
the low-cost devices. Furthermore, these features will become even more
prominent with the evolution of mMTC in future Beyond 5G (B5G) and 6G systems.
On the other hand, in the random access channel (RACH) procedure, each active
UE needs to establish a connection with the base station (BS) before data
transmission. However, these above-mentioned features of mMTC will bring new
challenges to the conventional RACH procedure, such as the shortage of uplink
resources, random access collisions, etc.

To address the above-mentioned challenges, the grant-free random access (GF-
RA) mechanism [2], also termed as grant-free non-orthogonal multiple access
(GF-NOMA), has recently emerged as a promising RACH solution for mMTC.
Specifically, in the GF-RA mechanism, active UEs can directly transmit their
payload data, without applying for the scheduling grant to establish a
connection with the BS. In this way, the GF-RA mechanism can skip the
handshaking procedure of the conventional grant-based random access (GB-RA)
mechanism, which reduces signaling overhead and improves the data-transmission
efficiency. In addition, enabled by different NOMA techniques, the GF-RA
mechanism also facilitates the sharing of access resources among all the
active UEs, which improves the resource utilization efficiency and mitigates
the collision problem. Since only a small amount of control signaling is
exchanged between the BS and active UEs in GF-RA, several critical problems
should be addressed at the BS. Specifically, the BS needs to identify the
active UEs, estimate the channel state information (CSI) from active UEs, and
finally detect their data. These problems are termed as user-activity
detection (UAD), channel estimation (CE), and multi-user data detection (MUD),
respectively. Existing solutions to these problems are briefed as follows.

In some typical GF-RA schemes [3, 4, 5, 6], each active UE transmits its
unique pilot sequence along with its payload data. In this way, the UAD and CE
problems are firstly addressed at the receiver, then the UAD and CE result
will be subsequently used to solve the MUD problem. Furthermore, several novel
schemes [7, 8, 9] have been recently proposed to jointly address the UAD, CE,
and MUD problems at the receiver, so as to further improve the UAD and MUD
performance. Despite that different receiver designs and detection algorithms
have been proposed for GF-RA, most of existing works are built on the
assumption of ideally synchronized signal reception at the BS. In other words,
the data frame and data symbols from different UEs are assumed to reach the BS
at exactly the same time. However, this assumption of ideal synchronization
entails frequent coordination between the BS and UEs, which significantly
increases the coordination burden, and contradicts the low power consumption
requirement of the low-cost UEs. Confronted with this problem, the receiver
design for asynchronous GF-RA has attracted increasing research interest.

###  I-A Literature Review

The orthogonal frequency-division multiplexing (OFDM) technique can employ the
cyclic prefix (CP) to handle the interference caused by different transmission
delay. Therefore, different asynchronous GF-RA schemes [10, 11, 12, 13, 14]
have been proposed under the OFDM-based NOMA framework. In [11, 12, 13], the
CP is used not only to handle the multi-path delay, but also to compensate for
the asynchronization between different users. However, this CP-based strategy
in [11, 12, 13] requires that the CP should be excessively long to accommodate
the large delay among all the UEs. In [14], the CP is only used to handle the
multi-path delay, so that the required CP length can be significantly reduced.
A triangular successive interference cancellation (T-SIC) algorithm was then
proposed in [14] to tackle the interference caused by asynchronization.
However, this T-SIC algorithm requires known user activity and CSI, as well as
sufficiently large power difference among different users to facilitate SIC
decoding. However, these assumptions may become impractical for mMTC with
densely deployed and randomly activated UEs. In addition, a common problem
with the OFDM-based GF-RA schemes [10, 11, 12, 13, 14] is that the OFDM
technique incurs high peak-to-average power ratio (PAPR), which poses
challenges for hardware implementations at the low-cost mMTC devices.

Apart from the OFDM technique, asynchronous GF-RA schemes have also been
extensively studied with other NOMA-enabling techniques, such as the
interleaver-division multiple access (IDMA) [15, 16], code-division multiple
access (CDMA) [17, 18], code-domain NOMA [19, 20, 21], and multiple-input
multiple-output (MIMO) [22, 23, 24, 25, 26]. We further categorize these
asynchronous GF-RA schemes according to different levels of asynchronization,
i.e. asynchronous frame with symbol synchronization [17, 18, 20, 21, 22, 23,
24, 25] and asynchronization within symbol duration [15, 16, 19, 26]. The
details are explained as follows.

The GF-RA schemes in [17, 18, 20, 21, 22, 23, 24, 25] assume that data frames
from different users are not synchronized at the receiver, yet the sampling at
the BS is _ideally_ aligned with the symbol interval, so that symbol-level
synchronization can be achieved. Under this assumption, compressed-sensing
based algorithms have been proposed for UAD and CE [21] or joint delay
learning (DL), UAD and MUD [17, 18], while a learned approximate message
passing (LAMP) algorithm was proposed for joint DL, UAD and CE [20]. In the
MIMO system, a generalized AMP (GAMP) algorithm was proposed for joint UAD and
CE [25], a block coordinate descent algorithm [23] and a LAMP algorithm [24]
were proposed for the joint DL, UAD and CE problem, while a Bilinear GAMP
(BiGAMP) algorithm [25] was further proposed for joint DL, UAD, CE, and MUD.
However, the assumption of symbol synchronization in [17, 18, 20, 21, 22, 23,
24, 25], i.e. sampling perfectly aligned with the symbol interval, can be hard
to realize due to the lack of coordination among the massive UEs.

In contrast to above-mentioned schemes, the asynchronous GF-RA schemes in [15,
16, 19, 26] dropped the infeasible assumption of symbol-level synchronization,
and consider the asynchronization within each symbol duration. A transceiver
structure was proposed in [15] with auxiliary preamble designed for UAD, and a
cross-symbol message passing MUD algorithm was proposed in [16], while a block
AMP algorithm was proposed in [19] for the CE and MUD problem. However, _ideal
knowledge_ on user activity [16, 19] or CSI [16, 15] is still required in
these schemes. Most recently, a turbo AMP (TAMP) algorithm has been designed
in [26] for the joint DL, CE, and MUD problem. However, this algorithm is
designed for the novel _unsourced_ massive access (UMA) scenario, where the BS
is concerned with specific payload data instead of the identity of each active
UE. In other words, the UAD problem is _bypassed_ , but the _exact number of
active UEs_ is still required in [26], while this assumption may become
infeasible for the sourced random access scenario considered in [10, 11, 12,
13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25] and this paper.

###  I-B Motivations and Contributions

According to the literature review, it is concluded that existing asynchronous
GF-RA schemes are still subject to several limitations, such as the infeasible
assumption of symbol-level synchronization or ideal knowledge of the user
activity or CSI, etc. On the other hand, mMTC devices are commonly deployed at
fixed locations with no or low mobility [27, 28, 29], so that the BS can
potentially obtain the delay profile of the UEs111The assumption of available
delay profile is commonly adopted in [14, 16, 19]. UEs can transmit preamble
sequences to facilitate estimation of the delay profile at the BS, which is a
common practice in LTE-based RACH design [30]. Since mMTC devices are commonly
deployed with no or low mobility [27, 28, 29], the transmission delay of each
UE remains constant in a long time. Therefore, acquiring delay profile at the
BS will not incur excessive signaling overhead or frequent coordination for
the low-cost UEs.. Therefore, we are motivated to exploit this delay profile
to distinguish different UEs at the receiver. Furthermore, confronted with the
massiveness of UE and the complicated inter-symbol interference caused by
symbol asynchronization, it is desirable to develop a low-complexity detection
algorithm for the joint UAD and CE problem in asynchronous GF-RA. With the UAD
and CE result, the subsequent MUD problem can be readily addressed with
existing solutions [16].

In this paper, we propose a receiver design for asynchronous GF-RA, where the
delay profile (also called as delay information) is exploited at the BS to
distinguish different UEs. In the presence of asynchronization-induced inter-
symbol interference, we propose a partially uni-directional (PUD) message
passing based sparse Bayesian learning (SBL) algorithm for asynchronous UAD
and CE (PUDMP-SBL-aUADCE). Our main contributions are summarized as follows.

Figure 1: Asynchronous grant-free random access system model.

(i) A receiver design is proposed where the delay profile is exploited to
distinguish different UEs. In order to address the sample correlation problem
in this receiver design, a PUD factor graph is established. Based on this PUD
factor graph, a low-complexity message passing based solution, i.e. the PUDMP-
SBL-aUADCE algorithm is proposed to address the joint UAD and CE problem for
asynchronous GF-RA.

(ii) The message passing procedure is analyzed, which proves that the PUDMP-
SBL-aUADCE algorithm exhibits higher signal-to-interference-and-noise ratio
(SINR) in the asynchronous case than in the synchronous case. In other words,
the proposed receiver can effectively exploit asynchronization to suppress
multi-user interference (MUI), which improves the UAD and CE accuracy.

(iii) Considering the potential timing error resulting from the low-cost
oscillators or the mobility of the devices, we investigate the impacts of
imperfect delay profile at the BS, and further explain the advantages of
adopting the SBL method in the proposed PUDMP-SBL-aUADCE algorithm.

The rest of this paper is organized as follows. The system model and receiver
design are introduced in Section II. The proposed PUDMP-SBL-aUADCE algorithm
is explained in Section III, together with the complexity analysis, SINR
analysis, as well as the discussion on the impacts of imperfect delay
information. Simulation results are provided in Section IV, and the paper is
finally concluded in Section V.

_Notations:_ Scalars and vectors are written in italic letters and boldface
lower-case letters, respectively, while matrices are in boldface upper-case
letters. All the vectors are column vectors, and $(\cdot)^{T}$ is the
transpose operation. $\mathbb{E}[\cdot]$ and $\mathbb{V}[\cdot]$ take the
expectation and variance of a random variable, respectively.
$X\sim\mathcal{CN}(\mu,v)$ and $X\sim\mathcal{N}(\mu,v)$ represent that a
random variable $X$ follows a complex Gaussian distribution and a real
Gaussian distribution with mean $\mu$ and variance $v$, respectively.
$Gam(\gamma;\epsilon,\eta)$ and $\mathcal{CN}(x|\mu,v)$ denote the probability
density function (pdf) of a Gamma distribution and a complex Gaussian
distribution, respectively. The symbol $\propto$ means “proportional to”.
$\delta(x)$ is the Dirac function with $\delta(x)=0$ for $\forall x\neq 0$,
and $\int_{-\infty}^{\infty}\delta(x)dx=1$. In this paper, the terms UE,
device, and user are used interchangeably.

##  II System Model

Figure 2: (a) Signal waveform $s(t)$ with (b) its auto-correlation
$\rho_{\tau}$, and (c) the MUI factor $\rho^{2}_{\tau}+\rho^{2}_{T_{s}-\tau}$
of a rectangular waveform.

###  II-A Transmitter Side

As shown in Fig. 1, we consider an asynchronous cellular system, where one BS
serves $K$ potential UEs. The BS is equipped with $M$ antennas, while each UE
is equipped with one antenna. We further denote
$\mathcal{K}=\\{1,2,\ldots,K\\}$ and $\mathcal{M}=\\{1,2,\ldots,M\\}$ as the
UE-index set and antenna-index set, respectively. In each round of random
access, we assume that each UE is independently activated with a small
probability $p_{a}$, and the active UEs transmit their payload data in a
grant-free manner. Specifically, if the $k$-th UE is activated, it will
transmit its user-specific pilot sequence
$\bm{p}_{k}=[p^{1}_{k},\ldots,p_{k}^{L_{p}}]^{T}$, followed by a guard
interval (GI) and its payload data sequence
$\bm{d}_{k}=[d^{1}_{k},\ldots,d_{k}^{L_{d}}]^{T}$, where $L_{p}$ and $L_{d}$
denote the pilot length and the data length, respectively. The set of pilot-
symbol indexes is denoted as $\mathcal{L}_{p}=\\{1,2,\ldots,L_{p}\\}$.

Before carrier modulation, each pilot symbol and data symbol is passed through
a shaping filter [16, 26], which is characterized by a pulse-shaping signal
waveform $s(t)$. We note that the proposed method in this paper does not rely
on any particular pulse shape, but for simplicity, we assume a rectangular
waveform $s(t)$ in Fig. 2(a), which has unit energy during the symbol period
$T_{s}$, i.e. $\int_{0}^{T_{s}}s^{2}(t)dt=1$, and becomes zero outside the
interval $[0,T_{s}]$. The auto-correlation $\rho_{\tau}$ and the MUI factor
$\rho^{2}_{\tau}+\rho^{2}_{T_{s}-\tau}$ of this rectangular waveform $s(t)$
are plotted in Fig. 2(b) and Fig. 2(c) respectively, and they will be later
explained in the following sections. The waveform for other commonly-used
shaping filters can be found in [16], and therefore not presented here for
brevity. In this way, the baseband waveform signal for the $l$-th pilot symbol
$p_{k}^{l}$ of the $k$-th UE is obtained as
$p_{k}^{l}s\big{(}t-(l-1)T_{s}\big{)}$. Then, the baseband waveform signal is
modulated, and transmitted to the BS.

Denote $\tau_{k}$ as the transmission delay of the $k$-th UE, which is mainly
determined by the distance from the $k$-th UE to the BS. We do not perform
synchronization to compensate for different transmission delay, which relieves
the coordination burden for the low-cost UEs. Furthermore, we will later show
in Section III-E that the proposed receiver design can effectively exploit the
asynchronization to suppress MUI. It is assumed that the delay information
$\\{\tau_{k},k\in\mathcal{K}\\}$ can be obtained by the BS, since mMTC devices
are mostly deployed at fixed locations with no or low mobility [27, 28, 29].
Still, we will also investigate the impacts of imperfect delay information in
Section III-F, which may be potentially caused by device mobility or the
timing error of the low-cost local oscillators. In addition, the maximum
relative delay between any two UEs is defined as
$\Delta^{\max}_{\tau}=\max\limits_{\forall
k,k^{\prime}\in\mathcal{K}}|\tau_{k}-\tau_{k^{\prime}}|$. Then, the length of
the guard interval (GI) should be configured to be larger than
$\Delta^{\max}_{\tau}$, so that the received pilot signal will not overlap
with the received data signal.

###  II-B Receiver Side

The delay profile $\\{\tau_{k},k\in\mathcal{K}\\}$ is exploited at the
receiver to distinguish different UEs. Specifically, the matched filtering
operation and sampling operation are performed at the BS on each antenna, and
each matched filter is aligned with the delay $\tau_{k}$ of one UE
$k\in\mathcal{K}$.

Take the $m$-th antenna as an example, and denote $y^{m}(t)$ as the received
baseband waveform signal, which is written as [16, 26]

$y^{m}(t)=\sum\limits^{K}_{k=1}\sum\limits^{L_{p}}_{l=1}\alpha_{k}g^{m}_{k}p^{l}_{k}s\Big{(}t-(l-1)T_{s}-\tau_{k}\Big{)}+n^{m}(t),$
(1)

where $g^{m}_{k}=\bar{g}^{m}_{k}\exp^{-j2\pi f_{c}\tau_{k}}$ is the baseband
equivalent channel gain [26], $\bar{g}^{m}_{k}\sim\mathcal{CN}(0,1)$ is the
channel coefficient from the $k$-th UE to the $m$-th antenna, $f_{c}$
represents the carrier modulation frequency, $n^{m}(t)$ denotes the additive
white Gaussian noise (AWGN), $\alpha_{k}$ is the activity indicator, i.e.
$\alpha_{k}=1$ if the $k$-th UE is activated, and $\alpha_{k}=0$ otherwise. In
addition, we assume a block fading channel, i.e., the channel gains
$\\{g^{m}_{k},m\in\mathcal{M},k\in\mathcal{K}\\}$ remain constant during each
round of random access.

For notation simplicity, the UEs are indexed by the value of the transmission
delay, i.e. $\tau_{k^{\prime}}\leq\tau_{k}$ for $\forall k^{\prime}<k$, and
the maximum relative delay is $\Delta^{\max}_{\tau}=\tau_{K}-\tau_{1}$. Since
we mainly focus on the asynchronization within each symbol duration, it is
assumed that the maximum relative delay is shorter than one symbol period,
i.e. $\Delta^{\max}_{\tau}<T_{s}$. This assumption does not sacrifice any
generality, and we take a two-UE case with $\tau_{1}=0.2T_{s}$ and
$\tau_{2}=1.5T_{s}$ as an example for explanation. In this case, according to
(1), the $l$-th pilot symbol of UE 1 is overlapped with the $(l-2)$-th and
$(l-1)$-th pilot symbols of UE 2 at the receiver side. Simply increasing the
pilot index of UE 2 by one, we have the effective delay
$\tau^{\text{eff}}_{2}=\tau_{2}-T_{s}=0.5T_{s}$ for UE 2. In this way, the
assumption that $\Delta^{\max}_{\tau}<T_{s}$ can always hold for any UE number
$K$.

The received pilot waveform signal $y^{m}(t)$ is passed through $K$ user-
oriented matched filters in parallel, where the $k$-th matched filter is
aligned with the delay $\tau_{k}$ of the $k$-th UE. Then, denote
$y^{m}_{k}(l)$ as the $l$-th sampled output of the $k$-th matched filter on
the $m$-th antenna, and we have

$\begin{split}&y^{m}_{k}(l)\overset{(a)}{=}\int^{lT_{s}+\tau_{k}}_{(l-1)T_{s}+\tau_{k}}y^{m}(t)s\Big{(}t-(l-1)T_{s}-\tau_{k}\Big{)}dt,\\\
&\overset{(b)}{=}\alpha_{k}g^{m}_{k}p^{l}_{k}+\sum\limits_{k^{\prime}<k}\alpha_{k^{\prime}}{g}^{m}_{k^{\prime}}\Big{(}p^{l}_{k^{\prime}}\rho_{\tau_{k}-\tau_{k^{\prime}}}+p^{l+1}_{k^{\prime}}\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}\Big{)}\\\
&\ \ \
+\sum\limits_{k^{\prime}>k}\alpha_{k^{\prime}}{g}^{m}_{k^{\prime}}\Big{(}p^{l-1}_{k^{\prime}}\rho_{\tau_{k^{\prime}}-\tau_{k}}+p^{l}_{k^{\prime}}\rho_{T_{s}-(\tau_{k^{\prime}}-\tau_{k})}\Big{)}+n^{m}_{k}(l),\\\
&\overset{(c)}{=}\sum\limits^{K}_{k^{\prime}=1}\alpha_{k^{\prime}}g^{m}_{k^{\prime}}\bar{p}^{l}_{k,k^{\prime}}+n^{m}_{k}(l),\\\
\end{split}$ (2)

where equation ($b$) of (2) is obtained by substituting (1) into equation
($a$), and $n^{m}_{k}(l)$ is the noise sample with distribution
$n^{m}_{k}(l)\sim\mathcal{CN}(0,\sigma^{2}_{n})$. Define the _auto-
correlation_ $\rho_{\tau}$ of $s(t)$ as
$\rho_{\tau}\overset{\Delta}{=}\int_{0}^{T_{s}}s(t)s(t-\tau)dt$, which has the
property $\rho_{-\tau}=\rho_{\tau}$ by definition. It is shown that
$\rho_{\tau}$ is a function of the delay $\tau$, and $\rho_{\tau}$ for the
rectangular waveform is plotted in Fig. 2(b). Furthermore, we denote
$\bar{p}^{l}_{k,k^{\prime}}$ as the $l$-th _effective_ pilot symbol of the
$k^{\prime}$-th UE at the $k$-th matched filter, so that the notations can be
unified as in equation ($c$) of (2), and $\bar{p}^{l}_{k,k^{\prime}}$ is
written as

$\bar{p}^{l}_{k,k^{\prime}}{=}\left\\{\begin{array}[]{l}p^{l}_{k},\ \ \ \ \ \
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \text{if}\
k^{\prime}=k,\\\
p^{l}_{k^{\prime}}\rho_{\tau_{k}-\tau_{k^{\prime}}}+p^{l+1}_{k^{\prime}}\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})},\
\ \text{if}\ k^{\prime}<k,\\\
p^{l-1}_{k^{\prime}}\rho_{\tau_{k^{\prime}}-\tau_{k}}+p^{l}_{k^{\prime}}\rho_{T_{s}-(\tau_{k^{\prime}}-\tau_{k})},\,\,\
\text{if}\ k^{\prime}>k.\end{array}\right.$ (3)

Figure 3: Multi-user inter-symbol interference caused by asynchronization.

As shown in equation ($b$) of (2) and Fig. 3, the $l$-th pilot symbol
$p^{l}_{k}$ of UE $k$ is interfered by two adjacent pilot symbols from another
UE $k^{\prime}$ due to the symbol asynchronization. If all the UEs share the
same delay, $y^{m}_{k}(l)$ in (2) will reduce to the symbol-synchronized case,
i.e. $y^{m}_{k}(l)=\sum^{K}_{k=1}\alpha_{k}g^{m}_{k}p^{l}_{k}+n^{m}_{k}(l)$,
since we have $\rho_{0}=1$ and $\rho_{T_{s}}=0$.

After the matched filtering and sampling operations on each antenna, we can
get the expression of the discrete received pilot signal in (2), which
facilitates subsequent UAD and CE.

##  III Partially Uni-Directional Message Passing Based Sparse Bayesian
Learning for Asynchronous UAD and CE

###  III-A Problem Formulation

To facilitate joint UAD and CE at the BS, we first define the _effective
channel gain_ $h_{k}^{m}$ as $h_{k}^{m}\overset{\Delta}{=}\alpha_{k}g_{k}^{m}$
from the $k$-th UE to the $m$-th antenna. Furthermore, we organize all the
effective channel gains in a matrix
$\bm{H}=[\bm{h}_{1},\ldots,\bm{h}_{K}]^{T}$ with size $K\times M$, where
$\bm{h}_{k}=[h^{1}_{k},\ldots,h_{k}^{m},\ldots,h^{M}_{k}]^{T}$ is the
effective channel-gain vector of UE $k$. We can also write
$\bm{H}=[\bm{h}^{1},\ldots,\bm{h}^{m},\ldots,\bm{h}^{M}]$ with
$\bm{h}^{m}=[h^{m}_{1},\ldots,h^{m}_{K}]^{T}$. It is noted that $\alpha_{k}=0$
for inactive UEs, and the matrix $\bm{H}$ therefore exhibits the _row-wise_
sparsity. If the $k$-th UE is inactive, we have
$\bm{h}^{T}_{k}=\bm{0}_{1\times M}$ for the $k$-th row of $\bm{H}$. If the
$k$-th UE is active, we have $\bm{h}^{T}_{k}\sim\mathcal{CN}(\bm{0}_{1\times
M},\bm{I}_{M})$, where $\bm{I}_{M}$ is an identity matrix with size $M\times
M$, indicating that each element of $\bm{h}_{k}$ is independently and
identically distributed (_i.i.d_). According to the row-wise sparsity of
$\bm{H}$, we follow the sparse Bayesian learning (SBL) method, and adopt a
hierarchical prior probability for $\bm{H}$, i.e.

$p(\bm{H})=\int_{\bm{\gamma}}p(\bm{H}|\bm{\gamma})p(\bm{\gamma})d\bm{\gamma},$
(4)

where $p(\bm{H}|\bm{\gamma})$ is the conditional prior pdf of $\bm{H}$,
$p(\bm{\gamma})$ is the prior pdf of
$\bm{\gamma}=[\gamma_{1},\ldots,\gamma_{k},\ldots,\gamma_{K}]^{T}$, and
$\gamma_{k}$ is the activity-related hyper-parameter for the $k$-th UE.

We further organize all the received pilot signals on the $m$-th antenna in a
matrix
$\bm{Y}^{m}=[\bm{y}^{m}_{1},\ldots,\bm{y}^{m}_{k},\ldots,\bm{y}^{m}_{K}]$ with
size $L_{p}\times K$, where
$\bm{y}^{m}_{k}=[{y}^{m}_{k}(1),\ldots,{y}^{m}_{k}(l),\ldots,{y}^{m}_{k}(L_{p})]^{T}$
is the sampled output from the $k$-th matched filter on the $m$-th antenna,
and $y^{m}_{k}(l)$ is given in (2). Then, the collection of $\bm{Y}^{m}$ on
all the antennas is denoted as $\\{\bm{Y}^{m},m\in\mathcal{M}\\}$.
Furthermore, we denote $\sigma^{2}_{n}$ as the variance of the noise sample
$n^{m}_{k}(l)$ in (2), and $\lambda=1/\sigma^{2}_{n}$ as the noise precision.
Both $\sigma^{2}_{n}$ and $\lambda$ are assumed unknown to the BS, while we
assume that the prior pdf $p(\lambda)$ follows $p(\lambda)\propto 1/\lambda$.
In this way, the joint a posteriori pdf of $\bm{H}$, $\bm{\gamma}$, and
$\lambda$ is written as

$\begin{split}&p\Big{(}\bm{H},\bm{\gamma},\lambda|\\{\bm{Y}^{m},m\in\mathcal{M}\\}\Big{)}\\\
&\propto
p\Big{(}\\{\bm{Y}^{m},m\in\mathcal{M}\\}|\bm{H},\lambda\Big{)}p(\bm{H}|\bm{\gamma})p(\bm{\gamma})p(\lambda),\\\
&\overset{(d)}{=}p(\lambda)\prod_{m=1}^{M}p\Big{(}\bm{Y}^{m}|\bm{h}^{m},\lambda\Big{)}\prod_{k=1}^{K}\prod_{m=1}^{M}p(h^{m}_{k}|\gamma_{k})\prod_{k=1}^{K}p(\gamma_{k}),\end{split}$
(5)

where equation ($d$) of (5) is obtained by the fact that the effective channel
gain $\bm{h}_{k}$ is independent among different UE, and that the received
pilot signal $\bm{Y}^{m}$ on different antennas are also mutually independent.
It is assumed that
$p(h^{m}_{k}|\gamma_{k})=\mathcal{CN}(h^{m}_{k};0,\gamma^{-1}_{k})$, so that
we have $\gamma_{k}=1$ if the $k$-th UE is active, and $\gamma_{k}=+\infty$
otherwise. Therefore, $\gamma_{k}$ is dubbed the activity-related hyper-
parameter. Following the SBL method [5], we assume
$p(\gamma_{k})=Gam(\gamma_{k};\epsilon,\eta)$, so that the UAD problem can be
solved by estimating $\gamma_{k},k\in\mathcal{K}$.

The joint UAD and CE problem is equivalent to finding $\bm{H}$, i.e. detecting
the nonzero channel gain vector $\bm{h}_{k}$ and further estimating its value.
However, finding the maximum a posteriori solution of $\bm{H}$ that maximizes
$p\big{(}\bm{H}|\\{\bm{Y}^{m},m\in\mathcal{M}\\}\big{)}=\int_{\bm{\gamma},\lambda}p\big{(}\bm{H},\bm{\gamma},\lambda|\\{\bm{Y}^{m},m\in\mathcal{M}\\}\big{)}d\bm{\gamma}d\lambda$
incurs prohibitively high computational complexity. As an alternative, we
establish a factor graph representation for the factorization in (5), and
further propose a low-complexity message passing based solution. However, we
observe a _sample correlation_ problem in the receiver design, which impedes
the factorization of the joint likelihood function
$p\big{(}\bm{Y}^{m}|\bm{h}^{m},\lambda\big{)}$ in (5), and therefore
complicates the UAD and CE problem. To address this problem, a novel PUD
factor graph is established, based on which the PUDMP-SBL-aUADCE algorithm is
derived. The details are explained as follows.

###  III-B Sample Correlation Problem and PUD Factor Graph

We consider the term $p\big{(}\bm{Y}^{m}|\bm{h}^{m},\lambda\big{)}$ in (5),
which is the joint likelihood function of all the sampled outputs
$\\{\bm{y}^{m}_{k},k\in\mathcal{K}\\}$ from $K$ different matched filters on
the $m$-th antenna. According to (2), the $k$-th column $\bm{y}^{m}_{k}$ of
$\bm{Y}^{m}$ is written as

$\bm{y}^{m}_{k}=\bar{\bm{P}}_{k}\bm{h}^{m}+\bm{n}^{m}_{k},$ (6)

where the element $\bar{p}^{l}_{k,k^{\prime}}$ on the $l$-th row and
$k^{\prime}$-th column of $\bar{\bm{P}}_{k}$ is given in (3). It is noted that
the noise samples in $\bm{n}^{m}_{k}$, i.e.
$\\{{n}^{m}_{k}(l),l\in\mathcal{L}_{p}\\}$ are _mutually independent_ , since
they are sampled from non-overlapping time intervals. Hence, we have

$p(\bm{y}^{m}_{k}|\bm{h}^{m},\lambda)=\prod_{l=1}^{L_{p}}\mathcal{CN}\Big{(}{y}^{m}_{k}(l);(\bar{\bm{p}}^{l}_{k})^{T}\bm{h}^{m},\lambda^{-1}\Big{)},$
(7)

where $(\bar{\bm{p}}^{l}_{k})^{T}$ is the $l$-th row of $\bar{\bm{P}}_{k}$,
${y}^{m}_{k}(l)$ is the $l$-th element of $\bm{y}^{m}_{k}$. However, we cannot
further factorize the joint likelihood function
$p\Big{(}\bm{Y}^{m}|\bm{h}^{m},\lambda\Big{)}$ as
$\prod_{k=1}^{K}p(\bm{y}^{m}_{k}|\bm{h}^{m},\lambda)$ because the sampled
output from different matched filters, i.e. different columns of $\bm{Y}^{m}$
are correlated. Specifically, consider the example in Fig. 3, where the $l$-th
sampled output of the $k$-th matched filter $y^{m}_{k}(l)$ is given in (2).
The noise sample $n^{m}_{k}(l)$ in $y^{m}_{k}(l)$ is expressed as

$n^{m}_{k}(l)=\int^{lT_{s}+\tau_{k}}_{(l-1)T_{s}+\tau_{k}}n^{m}(t)s\Big{(}t-(l-1)T_{s}-\tau_{k}\Big{)},$
(8)

and we can obtain a similar expression for the noise sample $n^{m}_{k-1}(l)$
in $y^{m}_{k-1}(l)$ from the ($k-1$)-th matched filter. Since the Gaussian
noise $n(t)$ is white and circularly symmetric, the correlation coefficient
$\bar{\rho}[n^{m}_{k}(l),n^{m}_{k-1}(l)]$ between $n^{m}_{k}(l)$ and
$n^{m}_{k-1}(l)$, and the correlation coefficient
$\bar{\rho}[n^{m}_{k}(l),n^{m}_{k-1}(l+1)]$ between $n^{m}_{k}(l)$ and
$n^{m}_{k-1}(l+1)$ are derived as

$\begin{split}&\bar{\rho}[n^{m}_{k}(l),n^{m}_{k-1}(l)]\\\
&=\frac{\mathbb{E}\Big{[}\big{(}n^{m}_{k}(l)\\!-\\!\mathbb{E}[n^{m}_{k}(l)]\big{)}\big{(}n^{*m}_{k\\!-\\!1}(l)\\!-\\!\mathbb{E}[n^{*m}_{k\\!-\\!1}(l)]\big{)}\Big{]}}{\sqrt{\mathbb{V}\big{(}n^{m}_{k}(l)\big{)}\mathbb{V}\big{(}n^{m}_{k\\!-\\!1}(l)\big{)}}},\\\
&=\frac{\mathbb{E}[n^{m}_{k}(l)n^{*m}_{k-1}(l)]}{\sigma_{n}^{2}}=\rho_{\tau_{k}-\tau_{k-1}},\\\
&\bar{\rho}[n^{m}_{k}(l),n^{m}_{k-1}(l+1)]=\rho_{T_{s}-(\tau_{k}-\tau_{k-1})},\end{split}$
(9)

where $n^{*m}_{k-1}(l)$ is the conjugate of $n^{m}_{k-1}(l)$. We can easily
generalize (9) as
$\bar{\rho}[n^{m}_{k}(l),n^{m}_{k^{\prime}}(l)]=\rho_{\tau_{k}-\tau_{k^{\prime}}}$
and
$\bar{\rho}[n^{m}_{k}(l),n^{m}_{k^{\prime}}(l+1)]=\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}$
for $\forall k^{\prime}<k$. It is shown in (9) that the noise samples from
different matched filters are correlated. This correlation problem impedes
further factorization of the joint likelihood function, i.e.
$p\big{(}\bm{Y}^{m}|\bm{h}^{m},\lambda\big{)}\neq\prod_{k=1}^{K}p(\bm{y}^{m}_{k}|\bm{h}^{m},\lambda)$,
which leads to a challenge in designing low-complexity message passing
algorithms.

In order to facilitate message passing based solution of the joint UAD and CE
problem, we need to establish a factor graph representation for the
factorization in (5). Messages passed on the factor graph are assumed to be
independent. However, the sampled output from different matched filters are
actually correlated as explained above, and this correlation problem will
drastically undermine the convergence of the message passing algorithm. To
address this correlation problem, we design a partially uni-directional (PUD)
factor graph representation for the joint likelihood function
$p\big{(}\bm{Y}^{m}|\bm{h}^{m},\lambda\big{)}$ as in Fig. 4.

Figure 4: Partially uni-directional factor graph for the $m$-th antenna.

The PUD factor graph in Fig. 4 is established for the channel gain vector
$\bm{h}^{m}=[h^{m}_{1},\ldots,h^{m}_{k},\ldots,h^{m}_{K}]^{T}$ and received
signal $\bm{Y}^{m}$ on the $m$-th antenna. Specifically, the variable node
(VN) $h^{m}_{k}$ represents the channel gain from the $k$-th UE to the $m$-th
antenna. The VN $z^{m,k}_{l}$ is defined as the noiseless received signal
$z^{m,k}_{l}\overset{\Delta}{=}(\bar{\bm{p}}^{l}_{k})^{T}\bm{h}^{m}$. Since
the noise samples are mutually independent within the same matched filter but
correlated between different matched filters, we set $K$ noise precision VNs
$\\{\lambda_{k},k\in\mathcal{K}\\}$, one for each matched filter. The function
node (FN) $f^{m,k}_{\delta(l)}$ represents the equality constraint between two
types of VNs, i.e.
$f^{m,k}_{\delta(l)}(z^{m,k}_{l},\bm{h}^{m})=\delta\big{(}z^{m,k}_{l}-(\bar{\bm{p}}^{l}_{k})^{T}\bm{h}^{m}\big{)}$.
The FN $f_{y(l)}^{m,k}$ represents the pdf of the received signal
${y}^{m}_{k}(l)$, i.e.
$f_{y(l)}^{m,k}(z^{m,k}_{l},\lambda_{k})=\mathcal{CN}\big{(}{y}^{m}_{k}(l);z^{m,k}_{l},\lambda_{k}^{-1}\big{)}$.
In addition, since the noise samples between different antennas are also
independent, each noise precision VN $\lambda_{k}$ is also connected to the
FNs $\big{\\{}f^{m^{\prime},k}_{y(l)},l\in\mathcal{L}_{p}\big{\\}}$ located on
the PUD factor graph of the $m^{\prime}$-th antenna, where
$m^{\prime}\in\mathcal{M}/m$.

Figure 5: Factor graph representation of the posterior probability
factorization.

The principle of partially uni-directional message passing is explained as
follows. As shown in Fig. 4, the VNs $\\{h^{m}_{k},k\in\mathcal{K}\\}$ and FNs
$\\{f^{m,k}_{\delta(l)},l\in\mathcal{L}_{p},k\in\mathcal{K}\\}$ are fully
connected. However, only the messages passed from the FNs
$\\{f^{m,k}_{\delta(l)},l\in\mathcal{L}_{p}\\}$ are employed to update the
belief of $h^{m}_{k}$. In other words, only the messages obtained from the
$k$-th matched filter are used to update the channel estimate of UE $k$. On
the other hand, the messages passed along the dashed edges in Fig. 4 are not
used by UE $k$ to update $h^{m}_{k}$, since these messages are obtained from
the other matched filters $\\{k^{\prime}\in\mathcal{K}/k\\}$, and therefore
correlated with the messages from the $k$-th matched filter. Apart from the
messages passed between FNs
$\\{f^{m,k}_{\delta(l)},l\in\mathcal{L}_{p},k\in\mathcal{K}\\}$ and VNs
$\\{h^{m}_{k},k\in\mathcal{K}\\}$, the remaining message passing procedures
between the VNs and FNs in Fig. 4 are still bi-directional, which explains why
this factor graph is _partially_ uni-directional. The detailed message passing
procedure will be explained in Section III-C.

Based on the PUD factor graph in Fig. 4, we establish a complete factor graph
in Fig. 5 for the factorization in equation ($d$) of (5). This factor graph is
composed of $M$ parallel PUD factor graphs, and additional sets of VNs
$\\{\gamma_{k},k\in\mathcal{K}\\}$, $\epsilon$, and FNs
$\\{f_{h}^{m,k},m\in\mathcal{M},k\in\mathcal{K}\\}$,
$\\{f_{\gamma}^{k},k\in\mathcal{K}\\}$, $f_{\lambda}$. Here,
$\\{\gamma_{k},k\in\mathcal{K}\\}$ and $\epsilon$ are related to the prior
variance of the VNs $\\{h^{m}_{k},k\in\mathcal{K},M\in\mathcal{M}\\}$. FN
$f_{h}^{m,k}$ and $f_{\lambda}$ represent the pdf
$f_{h}^{m,k}(h^{m}_{k},\gamma_{k})=p(h^{m}_{k}|\gamma_{k})$ in (5) and
$p(\lambda)\propto 1/\lambda$, while each FN $f_{\gamma}^{k}$ represents the
Gamma-distribution pdf
$f_{\gamma}^{k}(\gamma_{k},\epsilon)=p(\gamma_{k}|\epsilon)=Gam(\gamma_{k};\epsilon,\eta)$.
For reading convenience, FNs and related functions on the factor graph are
summarized in TABLE I.

######  Remark 1.

The pdf $p(\gamma_{k})$ in equation ($d$) of (5) is represented with a FN
$f_{\gamma}^{k}$ in the factor graph, while this FN $f_{\gamma}^{k}$ is
connected to a learnable shape parameter VN $\epsilon$. In the message passing
procedure described below, the shape parameter $\epsilon$ is updated
iteratively. According to [32], this shape parameter $\epsilon$ functions as a
selective amplifier for the VNs $\\{\gamma_{k},k\in\mathcal{K}\\}$. Therefore,
updating the shape parameter VN $\epsilon$ benefits the robustness of the
proposed PUDMP-SBL-aUADCE algorithm, especially under the case with non-ideal
delay information. Related details will be explained in Section III-F.

TABLE I: FNs in the factor graph and related functions FN | Function | PDF
---|---|---
$f_{\lambda}$ | $p(\lambda)$ | $\propto\lambda^{-1}$
$f^{m,k}_{y(l)}$ | $f_{y(l)}^{m,k}(z^{m,k}_{l},\lambda_{k})$ | $\mathcal{CN}\big{(}{y}^{m}_{k}(l);z^{m,k}_{l},{\lambda^{-1}_{k}}\big{)}$
$f^{m,k}_{\delta(l)}$ | $f^{m,k}_{\delta(l)}(z^{m,k}_{l},\bm{h}^{m})$ | $\delta\big{(}z^{m,k}_{l}-(\bar{\bm{p}}^{l}_{k})^{T}\bm{h}^{m}\big{)}$
$f_{h}^{m,k}$ | $p(h^{m}_{k}|\gamma_{k})$ | $\mathcal{CN}(h^{m}_{k};0,\gamma^{-1}_{k})$
$f_{\gamma}^{k}$ | $p(\gamma_{k}|\epsilon)$ | $Gam(\gamma_{k};\epsilon,\eta)$
|  |

###  III-C Message Passing Procedure and Message Update Rules

In this subsection, we explain the _forward_ (from the left to the right)
message passing procedure and the _backward_ (from the right to the left)
message passing procedure in Fig. 5, as well as related message update rules
at different VNs and FNs.

####  III-C1 Messages Update at FN $f^{m,k}_{\delta(l)}$

We first consider the forward message passing procedure at FN
$f^{m,k}_{\delta(l)}$ in the $t$-th iteration, and the forward message
$I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}$ passed from FN
$f^{m,k}_{\delta(l)}$ to VN $h^{m}_{k}$ is derived as

$\begin{split}I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})&\\!=\\!\int
f^{m,k}_{\delta(l)}(z^{m,k}_{l},\bm{h}^{m})I^{t}_{z^{m,k}_{l}\to
f^{m,k}_{\delta(l)}}(z^{m,k}_{l})\\\
&\\!\times\\!\\!\prod\limits_{k^{\prime}\neq
k}\\!\\!I^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}(h^{m}_{k^{\prime}})dz^{m,k}_{l}d\\{h^{m}_{k^{\prime}},k^{\prime}\in\mathcal{K}/k\\},\\\
&\\!=\\!\mathcal{CN}(h^{m}_{k};\mu^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}},v^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}),\end{split}$ (10)

where

$\begin{split}v^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}&=\Big{(}{\hat{\lambda}_{k}}^{-1}+\sum\limits_{k^{\prime}\neq
k}|\bar{p}^{l}_{k,k^{\prime}}|^{2}v^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}\Big{)}\Big{/}{|\bar{p}^{l}_{k,k}|^{2}},\\\
\mu^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}&=\Big{(}y^{m}_{k}(l)-\sum\limits_{k^{\prime}\neq
k}\bar{p}^{l}_{k,k^{\prime}}\mu^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}\Big{)}\Big{/}{\bar{p}^{l}_{k,k}},\end{split}$ (11)

where $I^{t}_{z^{m,k}_{l}\to
f^{m,k}_{\delta(l)}}(z^{m,k}_{l})=\mathcal{CN}(z^{m,k}_{l};y^{m}_{k}(l),{\hat{\lambda}_{k}}^{-1})$
is the forward message passed from VN $z^{m,k}_{l}$ to FN
$f^{m,k}_{\delta(l)}$, which will be later derived in (29).
$I^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}(h^{m}_{k^{\prime}})=\mathcal{CN}(h^{m}_{k^{\prime}};\mu^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}},v^{t-1}_{h^{m}_{k^{\prime}}\to f^{m,k}_{\delta(l)}})$ is
the backward message passed from VN $h^{m}_{k^{\prime}}$ to FN
$f^{m,k}_{\delta(l)}$ in the $(t-1)$-th iteration, which will be later derived
in (18).

Similarly, we can also derive the backward message
$I^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}(z^{m,k}_{l})$ passed from FN
$f^{m,k}_{\delta(l)}$ to VN $z^{m,k}_{l}$ in the $t$-th iteration, i.e.

$\begin{split}&I^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}(z^{m,k}_{l})\\\
&=\int
f^{m,k}_{\delta(l)}(z^{m,k}_{l},\bm{h}^{m})\prod\limits_{k^{\prime}\in\mathcal{K}}I^{t}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}(h^{m}_{k})d\\{h^{m}_{k},k\in\mathcal{K}\\},\\\
&=\mathcal{CN}(z^{m,k}_{l};\mu^{t}_{f^{m,k}_{\delta(l)}\to
z^{m,k}_{l}},v^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}),\end{split}$ (12)

where

$\begin{split}v^{t}_{f^{m,k}_{\delta(l)}\to
z^{m,k}_{l}}&=\sum\limits_{k^{\prime}\in\mathcal{K}}|\bar{p}^{l}_{k,k^{\prime}}|^{2}v^{t}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}},\\\ \mu^{t}_{f^{m,k}_{\delta(l)}\to
z^{m,k}_{l}}&=\sum\limits_{k^{\prime}\in\mathcal{K}}\bar{p}^{l}_{k,k^{\prime}}\mu^{t}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}.\end{split}$ (13)

####  III-C2 Messages Update at VN $h^{m}_{k}$

For each VN $h^{m}_{k}$, its belief $b^{t}(h^{m}_{k})$ in the $t$-th iteration
is obtained as

$\begin{split}b^{t}(h^{m}_{k})&\propto I^{t}_{f^{m,k}_{h}\to
h^{m}_{k}}(h^{m}_{k})\prod_{l\in\mathcal{L}_{p}}I^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k}),\\\
&\propto\mathcal{CN}(h^{m}_{k};\mu^{t}_{h^{m}_{k}},v^{t}_{h^{m}_{k}}),\end{split}$
(14)

with

$\begin{split}v^{t}_{h^{m}_{k}}&=\Big{(}\sum\limits_{l\in\mathcal{L}_{p}}\frac{1}{v^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}}+\hat{\gamma}_{k}\Big{)}^{-1},\\\
\mu^{t}_{h^{m}_{k}}&=v^{t}_{h^{m}_{k}}\Bigg{(}\sum\limits_{l\in\mathcal{L}_{p}}\frac{\mu^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}}{v^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}}\Bigg{)},\end{split}$
(15)

where $I^{t}_{f^{m,k}_{h}\to
h^{m}_{k}}(h^{m}_{k})=\mathcal{CN}(h^{m}_{k};0,\hat{\gamma}^{-1}_{k})$ is the
message passed from FN $f^{m,k}_{h}$ to VN $h^{m}_{k}$, which will be derived
in (23). The term $\prod_{l\in\mathcal{L}_{p}}I^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})$ implies that the belief of $h^{m}_{k}$ is updated
according to the principle of _partially uni-directional_ message passing in
Fig. 4. Then, we derive the following backward messages passed from VN
$h^{m}_{k}$ to different FNs $f^{m,k}_{\delta(l)}$ and
$f^{m,k^{\prime}}_{\delta(l)}$.

Denote $I^{t}_{h^{m}_{k}\to f^{m,k}_{\delta(l)}}(h^{m}_{k})$ as the backward
message passed from the VN $h^{m}_{k}$ of the $k$-th UE to FN
$f^{m,k}_{\delta(l)}$ in the $k$-th matched filter, and it is derived as

$I^{t}_{h^{m}_{k}\to
f^{m,k}_{\delta(l)}}(h^{m}_{k})\propto\frac{b^{t}(h^{m}_{k})}{I^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})}.$ (16)

Denote $I^{t}_{h^{m}_{k}\to f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})$ as the
backward message passed from the VN $h^{m}_{k}$ of the $k$-th UE to the FN
$f^{m,k^{\prime}}_{\delta(l)}$ in the $k^{\prime}$-th matched filter, and it
is derived as

$I^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})\\!\propto\\!\left\\{\begin{array}[]{l}\\!\\!\frac{b^{t}(h^{m}_{k})}{I^{t}_{f^{m,k}_{\delta(l-1)}\\!\to\\!h^{m}_{k}}\\!\\!(h^{m}_{k})I^{t}_{f^{m,k}_{\delta(l)}\\!\to\\!h^{m}_{k}}\\!\\!(h^{m}_{k})},\
\text{for }k^{\prime}<k,\\\
\\!\\!\frac{b^{t}(h^{m}_{k})}{I^{t}_{f^{m,k}_{\delta(l)}\\!\to\\!h^{m}_{k}}\\!\\!(h^{m}_{k})I^{t}_{f^{m,k}_{\delta(l+1)}\\!\to\\!h^{m}_{k}}\\!\\!(h^{m}_{k})},\
\text{for }k^{\prime}>k.\end{array}\right.$ (17)

We take the case with $k^{\prime}<k$ to explain the derivation in (17).
According to (11), the forward messages $I^{t}_{f^{m,k}_{\delta(l-1)}\to
h^{m}_{k}}(h^{m}_{k})$ and $I^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})$ involve the sampled output $y^{m}_{k}(l-1)$ and
$y^{m}_{k}(l)$, respectively. It is shown in Fig. 3 that both $y^{m}_{k}(l-1)$
and $y^{m}_{k}(l)$ are correlated with $y^{m}_{k^{\prime}}(l)$ for
$k^{\prime}<k$. Therefore, the forward messages
$I^{t}_{f^{m,k}_{\delta(l-1)}\to h^{m}_{k}}(h^{m}_{k})$ and
$I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ are also correlated with
$y^{m}_{k^{\prime}}(l)$. When $y^{m}_{k^{\prime}}(l)$ is used to calculate
$I^{t+1}_{f^{m,k^{\prime}}_{\delta(l)}\to
h^{m}_{k^{\prime}}}(h^{m}_{k^{\prime}})$ at FN $f^{m,k^{\prime}}_{\delta(l)}$
as in (10), the incoming message $I^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})$ should be independent from
$y^{m}_{k^{\prime}}(l)$ to avoid the correlation problem. Therefore,
$I^{t}_{f^{m,k}_{\delta(l-1)}\to h^{m}_{k}}(h^{m}_{k})$ and
$I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ should be excluded from
$b^{t}(h_{k}^{m})$ to obtain the backward message $I^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})$ in (17). Similarly,
$I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ and
$I^{t}_{f^{m,k}_{\delta(l+1)}\to h^{m}_{k}}(h^{m}_{k})$ should be excluded
from $b^{t}(h_{k}^{m})$ for $k^{\prime}>k$. According to (16) and (17), the
backward message update at VN $h^{m}_{k}$ is tedious for different UE indexes.
Sacrificing a little useful information, we obtain a more unified message
update rule for $\forall k^{\prime}\in\mathcal{K}$, i.e.,

$\begin{split}I^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})&\propto\frac{b^{t}(h^{m}_{k})}{\prod\limits_{\zeta=l-1}^{l+1}I^{t}_{f^{m,k}_{\delta(\zeta)}\to
h^{m}_{k}}(h^{m}_{k})},\\\ &=\mathcal{CN}(h^{m}_{k};\mu^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}},v^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}),\end{split}$ (18)

with

$\begin{split}v^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}&=\Big{(}\frac{1}{v^{t}_{h^{m}_{k}}}-\sum\limits_{\zeta=l-1}^{l+1}\frac{1}{v^{t}_{f^{m,k}_{\delta(\zeta)}\to
h^{m}_{k}}}\Big{)}^{-1},\\\ \mu^{t}_{h^{m}_{k}\to
f^{m,k^{\prime}}_{\delta(l)}}&=v^{t}_{h^{m}_{k}}\Bigg{(}\frac{\mu^{t}_{h^{m}_{k}}}{v^{t}_{h^{m}_{k}}}-\sum\limits_{\zeta=l-1}^{l+1}\frac{\mu^{t}_{f^{m,k}_{\delta(\zeta)}\to
h^{m}_{k}}}{v^{t}_{f^{m,k}_{\delta(\zeta)}\to h^{m}_{k}}}\Bigg{)}.\end{split}$
(19)

In other words, we exclude three messages $I^{t}_{f^{m,k}_{\delta(l-1)}\to
h^{m}_{k}}(h^{m}_{k})$, $I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$,
and $I^{t}_{f^{m,k}_{\delta(l+1)}\to h^{m}_{k}}(h^{m}_{k})$ from the belief
$b^{t}(h^{m}_{k})$ to obtain the backward message from VN $h^{m}_{k}$ to an
arbitrary FN $f^{m,k^{\prime}}_{\delta(l)}$ with $\forall
k^{\prime}\in\mathcal{K}$. Replacing (16) and (17) with (18) may lose some
useful uncorrelated information, while the loss is negligible if $L_{p}\gg 1$.

####  III-C3 Message Update at FN $f_{h}^{m,k}$, VN $\gamma_{k}$, FN
$f_{\gamma}^{k}$, and VN $\epsilon$

The forward message $I^{t}_{f_{h}^{m,k}\to\gamma_{k}}(\gamma_{k})$ passed from
FN $f_{h}^{m,k}$ to VN $\gamma_{k}$ is derived according to the mean-field
(MF) rule, i.e.

$\begin{split}I^{t}_{f_{h}^{m,k}\to\gamma_{k}}(\gamma_{k})&=\exp\Bigg{\\{}\int
b^{t}(h^{m}_{k})\ln f_{h}^{m,k}(h^{m}_{k},\gamma_{k})dh^{m}_{k}\Bigg{\\}},\\\
&\propto\gamma_{k}\exp\Big{\\{}-\gamma_{k}\Big{(}|\mu^{t}_{h^{m}_{k}}|^{2}+v^{t}_{h^{m}_{k}}\Big{)}\Big{\\}}.\end{split}$
(20)

The belief $b^{t}(\gamma_{k})$ of $\gamma_{k}$ is obtained by combining all
the incoming messages to VN $\gamma_{k}$,

$\begin{split}b^{t}(\gamma_{k})&\propto
I^{t-1}_{f^{k}_{\gamma}\to\gamma_{k}}(\gamma_{k})\prod\limits^{M}_{m=1}I^{t}_{f_{h}^{m,k}\to\gamma_{k}}(\gamma_{k}),\\\
&\propto{\gamma_{k}}^{\hat{\epsilon}+M-1}\exp\Bigg{\\{}\\!-\\!\gamma_{k}\Bigg{[}\eta+\sum\limits_{m=1}^{M}\Big{(}|\mu^{t}_{h^{m}_{k}}|^{2}\\!+\\!v^{t}_{h^{m}_{k}}\Big{)}\Bigg{]}\Bigg{\\}},\end{split}$
(21)

where
$I^{t-1}_{f^{k}_{\gamma}\to\gamma_{k}}(\gamma_{k})=Gam(\gamma_{k};\hat{\epsilon},\eta)$
is the backward message from FN $f_{\gamma}^{k}$ to VN $\gamma_{k}$ in the
previous iteration, which will be derived in (26). Then, the hyper-parameter
$\hat{\gamma}_{k}$ is estimated as

$\hat{\gamma}_{k}\overset{\Delta}{=}\mathbb{E}_{b^{t}(\gamma_{k})}\big{[}\gamma_{k}\big{]}=\frac{\hat{\epsilon}+M}{\eta+\sum\limits_{m=1}^{M}\Big{(}|\mu^{t}_{h^{m}_{k}}|^{2}+v^{t}_{h^{m}_{k}}\Big{)}},$
(22)

and $\eta$ is commonly set as $0$ for the SBL method [32]

With the belief $b^{t}(\gamma_{k})$ of $\gamma_{k}$, we can obtain the
backward message $I^{t}_{f_{h}^{m,k}\to h^{m}_{k}}(h^{m}_{k})$ passed from FN
$f_{h}^{m,k}$ to VN $h^{m}_{k}$, which is derived according to the following
MF rule,

$\begin{split}I^{t}_{f_{h}^{m,k}\to h^{m}_{k}}(h^{m}_{k})&=\exp\Bigg{\\{}\int
b^{t}(\gamma_{k})\ln
f_{h}^{m,k}(h^{m}_{k},\gamma_{k})d\gamma_{k}\Bigg{\\}},\\\
&\propto\mathcal{CN}(h^{m}_{k};0,\hat{\gamma}^{-1}_{k}).\end{split}$ (23)

We further employ the MF rule to derive the forward message
$I^{t}_{f_{\gamma}^{k}\to\epsilon}(\epsilon)$ passed from FN $f_{\gamma}^{k}$
to VN $\epsilon$, i.e.

$\begin{split}I^{t}_{f_{\gamma}^{k}\to\epsilon}(\epsilon)&=\exp\Bigg{\\{}\int
b^{t}(\gamma_{k})\ln
f_{\gamma}^{k}(\gamma_{k},\epsilon)d\gamma_{k}\Bigg{\\}},\\\
&\propto\frac{\eta^{\epsilon}}{\Gamma(\epsilon)}\exp\Bigg{\\{}(\epsilon-1)\mathbb{E}_{b^{t}(\gamma_{k})}\big{[}\ln\gamma_{k}\big{]}-\eta\hat{\gamma}_{k}\Bigg{\\}}.\end{split}$
(24)

Combine all incoming messages to VN $\epsilon$, the belief $b^{t}(\epsilon)$
is

$\begin{split}b^{t}(\epsilon)&\propto\prod\limits_{k=1}^{K}I^{t}_{f_{\gamma}^{k}\to\epsilon}(\epsilon),\\\
&\propto\Big{(}\frac{\eta^{\epsilon}}{\Gamma(\epsilon)}\Big{)}^{K}\exp\Bigg{\\{}(\epsilon\\!-\\!1)\sum\limits_{k=1}^{K}\Big{(}\mathbb{E}_{b^{t}(\gamma_{k})}\big{[}\ln\gamma_{k}\big{]}\\!-\\!\eta\hat{\gamma}_{k}\Big{)}\Bigg{\\}},\end{split}$
(25)

and the backward message $I^{t}_{f_{\gamma}^{k}\to\gamma_{k}}(\gamma_{k})$
passed from FN $f_{\gamma}^{k}$ to VN $\gamma_{k}$ is further derived as

$\begin{split}I^{t}_{f_{\gamma}^{k}\to\gamma_{k}}(\gamma_{k})&=\exp\Bigg{\\{}\int
b^{t}(\epsilon)\ln f_{\gamma}^{k}(\gamma_{k},\epsilon)d\epsilon\Bigg{\\}},\\\
&\propto Gam(\gamma_{k};\hat{\epsilon},\eta),\end{split}$ (26)

where
$\hat{\epsilon}\overset{\Delta}{=}\mathbb{E}_{b^{t}(\epsilon)}\big{[}\epsilon\big{]}$.
However, it is difficult to directly calculate $\hat{\epsilon}$ by definition.
Here, we use another empirical yet effective solution [32] to update
$\hat{\epsilon}$ in a more straightforward way, i.e.,

$\hat{\epsilon}=\frac{1}{2}\sqrt{\ln\Bigg{(}\frac{1}{K}\sum\limits_{k=1}^{K}\hat{\gamma}_{k}\Bigg{)}-\frac{1}{K}\sum\limits_{k=1}^{K}\ln\hat{\gamma}_{k}}.$
(27)

####  III-C4 Message Update at FN $f^{m,k}_{y(l)}$ and VN $z^{m,k}_{l}$

The backward message $I^{t}_{f^{m,k}_{y(l)}\to z^{m,k}_{l}}$ passed from FN
$f^{m,k}_{y(l)}$ to VN $z^{m,k}_{l}$ is derived as

$\begin{split}I^{t}_{f^{m,k}_{y(l)}\to
z^{m,k}_{l}}(z^{m,k}_{l})&\\!=\\!\exp\Bigg{\\{}\\!\int
b^{t}\big{(}\lambda_{k}\big{)}\ln
f_{y(l)}^{m,k}\big{(}z^{m,k}_{l},\lambda_{k}\big{)}d\lambda_{k}\\!\Bigg{\\}},\\\
&\propto\mathcal{CN}\Big{(}z^{m,k}_{l};y^{m}_{k}(l),\hat{\lambda}^{-1}_{k}\Big{)},\end{split}$
(28)

where $b^{t}(\lambda_{k})$ and $\hat{\lambda}^{-1}_{k}$ will be later
explained in (LABEL:B_Lambda) and (34), respectively.

For each VN $z^{m,k}_{l}$, the forward message $I^{t}_{z^{m,k}_{l}\to
f^{m,k}_{\delta(l)}}(z^{m,k}_{l})$ is simply obtained as

$I^{t}_{z^{m,k}_{l}\to
f^{m,k}_{\delta(l)}}(z^{m,k}_{l})=I^{t}_{f^{m,k}_{y(l)}\to
z^{m,k}_{l}}(z^{m,k}_{l}).$ (29)

In this way, the belief $b^{t}(z^{m,k}_{l})$ of $z^{m,k}_{l}$ is obtained as

$\begin{split}b^{t}(z^{m,k}_{l})&\propto I^{t}_{f^{m,k}_{y(l)}\to
z^{m,k}_{l}}(z^{m,k}_{l})I^{t}_{f^{m,k}_{\delta(l)}\to
z^{m,k}_{l}}(z^{m,k}_{l}),\\\
&\propto\mathcal{CN}(z^{m,k}_{l};\mu^{t}_{z^{m,k}_{l}},v^{t}_{z^{m,k}_{l}}),\end{split}$
(30)

with

$\begin{split}v^{t}_{z^{m,k}_{l}}&=\Bigg{(}\frac{1}{v^{t}_{f^{m,k}_{\delta(l)}\to
z^{m,k}_{l}}}+\hat{\lambda}_{k}\Bigg{)}^{-1},\\\
\mu^{t}_{z^{m,k}_{l}}&=v^{t}_{z^{m,k}_{l}}\Bigg{(}y^{m}_{k}(l)\hat{\lambda}_{k}+\frac{\mu^{t}_{f^{m,k}_{\delta(l)}\to
z^{m,k}_{l}}}{v^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}}\Bigg{)}.\end{split}$
(31)

With the belief $b^{t}(z^{m,k}_{l})$, we can now calculate the forward message
$I^{t}_{f^{m,k}_{y(l)}\to\lambda_{k}}(\lambda_{k})$ passed from FN
$f^{m,k}_{y(l)}$ to VN $\lambda_{k}$. According to the MF rule, we have

$\begin{split}I^{t}_{f^{m,k}_{y(l)}\to\lambda_{k}}(\lambda_{k})&\\!=\\!\exp\Bigg{\\{}\\!\int
b^{t}\big{(}z^{m,k}_{l}\big{)}\ln
f_{y(l)}^{m,k}\big{(}z^{m,k}_{l},\lambda_{k}\big{)}dz^{m,k}_{l}\\!\Bigg{\\}},\\\
&\propto\mathcal{CN}\Big{(}z^{m,k}_{l};y^{m}_{k}(l),\hat{\lambda}^{-1}_{k}\Big{)}.\end{split}$
(32)

Then, the belief $b^{t}(\lambda_{k})$ of each VN $\lambda_{k}$ is derived as

$\begin{split}&b^{t}(\lambda_{k})\propto
f(\lambda)\prod\limits^{M}_{m=1}\prod\limits^{L_{p}}_{l=1}I^{t}_{f^{m,k}_{y(l)}\to\lambda_{k}}(\lambda_{k}),\\\
&\propto\lambda_{k}^{ML_{p}-1}\exp\\!\Bigg{\\{}\\!\\!-\\!\lambda_{k}\\!\sum\limits^{M}_{m=1}\sum\limits^{L_{p}}_{l=1}\Big{[}\big{|}y^{m}_{k}(l)\\!-\\!\mu^{t}_{z^{m,k}_{l}}\big{|}^{2}+v^{t}_{z^{m,k}_{l}}\Big{]}\Bigg{\\}}.\end{split}$
(33)

In this way, the noise precision $\hat{\lambda}_{k}$ is estimated as

$\hat{\lambda}_{k}\overset{\Delta}{=}\mathbb{E}_{b^{t}(\lambda_{k})}\Big{[}\lambda_{k}\Big{]}=\frac{ML_{p}}{\sum\limits^{M}_{m=1}\sum\limits^{L_{p}}_{l=1}\Big{[}\big{|}y^{m}_{k}(l)\\!-\\!\mu^{t}_{z^{m,k}_{l}}\big{|}^{2}+v^{t}_{z^{m,k}_{l}}\Big{]}}.$
(34)

###  III-D Algorithm Summary and Complexity Analysis

Input: $\\{\bm{Y}^{m},m\in\mathcal{M}\\}$,
$\\{\bm{\bar{P}}_{k},k\in\mathcal{K}\\}$

Output: UAD $\\{\hat{\alpha}_{k},k\in\mathcal{K}\\}$ and CE
$\\{\hat{h}^{m}_{k},m\in\mathcal{M},k\in\mathcal{K}\\}$

Initialize:

;

$\text{Set}\ \hat{\epsilon}=10^{-3},\ \hat{\lambda}_{k}=1,\
\hat{\gamma}_{k}=1\ \text{for}\ \forall k\text{, and set}\ t=0.$

$\text{Set}\ v^{t}_{h^{m}_{k}\to f^{m,k^{\prime}}_{\delta(l)}}=1,\
\mu^{t}_{h^{m}_{k}\to f^{m,k^{\prime}}_{\delta(l)}}=0\ \text{for}\ \forall k,\
\forall k^{\prime},\ \forall l,\ \forall m$.

for  _$t=1:N_{it}$_ do

1\. Update $I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ with (10) and
(11).

2\. Update $b^{t}(h^{m}_{k})$ with (14) and (15).

3\. Update $\hat{\gamma}_{k}$ with (22).

4\. Update $\hat{\epsilon}$ with (27).

5\. Update $I^{t}_{h^{m}_{k}\to f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})$ with
(18) and (19).

6\. Update $I^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}(z^{m,k}_{l})$ with (12)
and (13).

7\. Update $b^{t}(z^{m,k}_{l})$ with (30) and (31).

8\. Update $\hat{\lambda}_{k}$ with (34).

UAD and CE Decision: If $\hat{\gamma}_{k}<\gamma_{th}$, $\hat{\alpha}_{k}=1$
and $\hat{h}^{m}_{k}=\mu^{N_{it}}_{h^{m}_{k}}$ for $\forall m\in\mathcal{M}$.
Otherwise, $\hat{\alpha}_{k}=0$ and $\hat{h}^{m}_{k}=0$ for $\forall
m\in\mathcal{M}$.

Algorithm 1 PUDMP-SBL-aUADCE Algorithm

According to the descriptions above, the proposed PUDMP-SBL-aUADCE algorithm
is summarized as in Algorithm 1. The input information includes the received
pilot signal $\bm{Y}^{m}$ on all the $M$ antennas, and the equivalent pilot
matrix $\bm{\bar{P}}_{k}$ of each matched filter $k\in\mathcal{K}$. Here, the
element $\bar{p}^{l}_{k,k^{\prime}}$ on the $l$-th row and $k^{\prime}$-th
column of matrix $\bar{\bm{P}}_{k}$ is given in (3). After the initialization
step, the message passing procedure of the PUDMP-SBL-aUADCE algorithm is
performed with $N_{it}$ iterations, producing the activity-detection result
$\\{\hat{\alpha}_{k},k\in\mathcal{K}\\}$, where $\hat{\alpha}_{k}=1$ if
$\hat{\gamma}_{k}<\gamma_{th}$, and $\hat{\alpha}_{k}=0$ otherwise. The UAD
threshold is empirically set as $\gamma_{th}=10$ for the simulations in
Section IV. We further obtain the CE result for the UEs that are detected as
active, i.e. $\hat{h}^{m}_{k}=\mu^{N_{it}}_{h^{m}_{k}}$ if
$\hat{\alpha}_{k}=1$, where $\mu^{N_{it}}_{h^{m}_{k}}$ is the mean of the
belief distribution $b^{t}(h^{m}_{k})$ in the $N_{it}$-th iteration, which is
given by (15).

We further analyze the computational complexity of the proposed PUDMP-SBL-
aUADCE algorithm, which is assessed by the number of multiplication/division
operations. According to the partially uni-directional message passing
principle, the message $I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$
in step 1 of Algorithm 1 is independently updated for different UE $k$.
Considering each UE index $k$, pilot-symbol index $l$, and antenna index $m$,
updating $I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ totally incurs
$\mathcal{O}(K^{2}L_{p}M)$ operations in each iteration. Updating
$b^{t}(h^{m}_{k})$ and $\hat{\gamma}_{k}$ in step 2 and step 3 incurs
$\mathcal{O}(KL_{p}M)$ operations and $\mathcal{O}(KM)$ operations,
respectively. According to (27), the update of $\hat{\epsilon}$ in step 4 only
needs $K+1$ logarithmic operations and $1$ additional operation to take the
square root. Then, $\mathcal{O}(KL_{p}M)$ operations are required to update
$I^{t}_{h^{m}_{k}\to f^{m,k^{\prime}}_{\delta(l)}}(h^{m}_{k})$ in step 5,
while $\mathcal{O}(K^{2}L_{p}M)$ operations are required for
$I^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}(z^{m,k}_{l})$ in step 6. Finally,
updating both $b^{t}(z^{m,k}_{l})$ in step 7 and $\hat{\lambda}_{k}$ in step 8
will incur $\mathcal{O}(KL_{p}M)$ operations. To sum up, the proposed PUDMP-
SBL-aUADCE algorithm totally incurs $\mathcal{O}(N_{it}K^{2}L_{p}M)$
multiplication/division operations.

Despite that the overall computational complexity scales with the square of
the UE number $K$, it is noted that the computationally-intensive messages
$I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ and
$I^{t}_{f^{m,k}_{\delta(l)}\to z^{m,k}_{l}}(z^{m,k}_{l})$ are _independently_
updated for different UE $k$, which facilitates parallel implementation.
Therefore, we can allocate the computation task to $K$ distributed processors,
which is one major advantage of the proposed message passing algorithm. In
this way, the computational complexity at each local processor still scales
_linearly_ with the UE number $K$, the pilot length $L_{p}$, and the antenna
number $M$.

###  III-E SINR Analysis for Synchronous and Asynchronous Case

The synchronous case is equivalent to setting the same delay $\tau_{k}=\tau$
for $\forall k\in\mathcal{K}$. In this case, the proposed PUDMP-SBL-aUADCE
algorithm will reduce to the MP-BSBL algorithm in [5]. To make this paper
self-contained, we explain how to modify the PUDMP-SBL-aUADCE algorithm in
Algorithm 1 for synchronous GF-RA. Firstly, we adopt the same matched filter
to all the UEs, so that the received signal in (6) will be re-written as
$\bm{y}^{m}_{\text{syn}}=\bm{P}\bm{h}^{m}+\bm{n}^{m}$, where $\bm{P}$ is the
pilot matrix with element $p^{l}_{k}$ on its $l$-th row and $k$-th column. In
this way, the factor graph in Fig. 4 will contain $K$ VNs
$\\{h^{m}_{k},k\in\mathcal{K}\\}$, and only $L_{p}$ FNs
$\\{f^{m}_{\delta(l)},l\in\mathcal{L}_{p}\\}$, since we have
$f^{m,k}_{\delta(l)}=f^{m}_{\delta(l)}$ for $\forall k\in\mathcal{K}$ in the
synchronous case. Accordingly, there will be only $L_{p}$ VNs
$\\{z^{m}_{l},l\in\mathcal{L}_{p}\\}$, $L_{p}$ FNs
$\\{f^{m}_{\delta(l)},l\in\mathcal{L}_{p}\\}$, and one VN $\lambda$ in Fig. 4
and Fig. 5 in the synchronous case. Then, in step 2 of the message passing
algorithm, we will update the belief $b^{t}_{\text{syn}}(h^{m}_{k})$ of VN
$h^{m}_{k}$ with all the incoming messages $\\{I^{t}_{f^{m}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k}),l\in\mathcal{L}_{p}\\}$ and $I^{t}_{f^{m,k}_{h}\to
h^{m}_{k}}(h^{m}_{k})$. Furthermore, the message $I^{t}_{h^{m}_{k}\to
f^{m}_{\delta(l)}}(h^{m}_{k})$ in step 5 will be updated according to (16),
i.e. $I^{t}_{h^{m}_{k}\to f^{m}_{\delta(l)}}(h^{m}_{k})\propto
b^{t}_{\text{syn}}(h^{m}_{k})/I^{t}_{f^{m}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})$. Apart from step 2 and step 5, the remaining steps of
Algorithm 1 will remain _exactly_ the same for the synchronous case.

In the following, we will analyze the SINR at FN $f^{m}_{\delta(l)}$ for the
synchronous case and at FN $f^{m,k}_{\delta(l)}$ for the asynchronous case.
The result is presented in the following proposition.

######  Proposition 1.

Assuming that the synchronous case and the asynchronous case share the same
estimated noise precision $\hat{\lambda}_{k}$ and backward-passing variances
$\\{v^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}},k^{\prime}\in\mathcal{K}/k\\}$, the forward message from
FN $f^{m,k}_{\delta(l)}$ to VN $h^{m}_{k}$ has higher SINR in the asynchronous
case than in the synchronous case.

_Proof_ : According to the update rule of $I^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})$ in (10) and (11), the incoming messages
$I^{t-1}_{h^{m}_{k^{\prime}}\to f^{m,k}_{\delta(l)}}(h^{m}_{k^{\prime}})$ from
the other UEs $k^{\prime}\in\mathcal{K}/k$ are treated as interference at FN
$f^{m,k}_{\delta(l)}$. Therefore, under the asynchronous case, the SINR for
the message $I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ can be
represented as

$\text{SINR}^{\text{asy}}_{k}{=}\frac{|\bar{p}^{l}_{k,k}|^{2}}{\hat{\lambda}^{-1}_{k}+\sum\limits_{k^{\prime}\neq
k}v^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}|\bar{p}^{l}_{k,k^{\prime}}|^{2}}.$ (35)

Similarly, under the synchronous case, the SINR for the message
$I^{t}_{f^{m}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ can be expressed as

$\text{SINR}^{\text{syn}}_{k}{=}\frac{|p^{l}_{k}|^{2}}{\hat{\lambda}^{-1}+\sum\limits_{k^{\prime}\neq
k}v^{t-1}_{h^{m}_{k^{\prime}}\to f^{m}_{\delta(l)}}|p^{l}_{k^{\prime}}|^{2}},$
(36)

and we have $\bar{p}^{l}_{k,k}=p^{l}_{k}$ in (35) according to (3). Assuming
that $\hat{\lambda}_{k}=\hat{\lambda}$ and $v^{t-1}_{h^{m}_{k^{\prime}}\to
f^{m,k}_{\delta(l)}}=v^{t-1}_{h^{m}_{k^{\prime}}\to f^{m}_{\delta(l)}}$ in
(35) and (36), the difference between $\text{SINR}^{\text{asy}}_{k}$ and
$\text{SINR}^{\text{syn}}_{k}$ is then only determined by the difference
between $|\bar{p}^{l}_{k,k^{\prime}}|^{2}$ and $|p^{l}_{k^{\prime}}|^{2}$.
Without loss of any generality, take the $k^{\prime}$-th UE with
$k^{\prime}<k$ as an example, where we have
$\bar{p}^{l}_{k,k^{\prime}}=p^{l}_{k^{\prime}}\rho_{\tau_{k}-\tau_{k^{\prime}}}+p^{l+1}_{k^{\prime}}\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}$
according to (3). Assuming that each pilot symbol has unit energy, we
therefore have $\mathbb{E}\big{(}|p^{l}_{k^{\prime}}|^{2}\big{)}=1$ for
$\forall k^{\prime}\in\mathcal{K},l\in\mathcal{L}_{p}$. Then, the expectation
on $|p^{l}_{k,k^{\prime}}|^{2}$ is written as

$\begin{split}\mathbb{E}\big{(}|p^{l}_{k,k^{\prime}}|^{2}\big{)}&=\rho^{2}_{\tau_{k}-\tau_{k^{\prime}}}+\rho^{2}_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}\\\
&+2\rho_{\tau_{k}-\tau_{k^{\prime}}}\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}\mathbb{R}\Big{(}(p^{l}_{k^{\prime}})^{*}p^{l+1}_{k^{\prime}}\Big{)},\end{split}$
(37)

where $\mathbb{R}(\cdot)$ takes the real part of a complex number. Since each
pilot symbol has unit energy, it is easy to show that
$\mathbb{E}\big{(}|p^{l}_{k,k^{\prime}}|^{2}\big{)}$ in (37) is maximized when
$p^{l}_{k^{\prime}}=p^{l+1}_{k^{\prime}}$, which is an impractical assumption
for pilot design. Still, conditioned on this impractical assumption of
$p^{l}_{k^{\prime}}=p^{l+1}_{k^{\prime}}$, we have

$\begin{split}&\max\Bigg{(}\mathbb{E}\big{(}|p^{l}_{k,k^{\prime}}|^{2}\big{)}\Bigg{)}=\mathbb{E}\Big{(}|p^{l}_{k,k^{\prime}}|^{2}\Big{|}p^{l}_{k^{\prime}}=p^{l+1}_{k^{\prime}}\Big{)}\\\
&=|\rho_{\tau_{k}-\tau_{k^{\prime}}}+\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}|^{2}\leq
1=\mathbb{E}\big{(}|p^{l}_{k^{\prime}}|^{2}\big{)},\end{split}$ (38)

where
$|\rho_{\tau_{k}-\tau_{k^{\prime}}}+\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}|^{2}=1$
holds only for the rectangular waveform $s(t)$, which has the property
$\rho_{\tau}+\rho_{T_{s}-\tau}=\rho_{0}$. To sum up,
$\mathbb{E}\big{(}|p^{l}_{k,k^{\prime}}|^{2}\big{)}\leq\mathbb{E}\big{(}|p^{l}_{k^{\prime}}|^{2}\big{)}$,
and the equality holds only under the impractical assumption that
$p^{l}_{k^{\prime}}=p^{l+1}_{k^{\prime}}$.

The conclusion obtained from (38) is applicable to the general pilot design.
We further consider the _i.i.d_ Gaussian pilot symbols, i.e.
$p^{l}_{k}\sim\mathcal{CN}(0,1)$ for $\forall k\in\mathcal{K}$ and $\forall
l\in\mathcal{L}_{p}$, which is a common choice for pilot design in message
passing based detection algorithms [3, 9, 4]. In this case, we have
$\mathbb{R}\Big{(}(p^{l}_{k^{\prime}})^{*}p^{l+1}_{k^{\prime}}\Big{)}=0$ for
(37), and Proposition 1 can be proven with the Cauchy-Schwartz inequality,

$\begin{split}&\mathbb{E}\Big{(}|\bar{p}^{l}_{k,k^{\prime}}|^{2}\big{|}\text{Gaussian
pilot}\Big{)}=\rho^{2}_{\tau_{k}-\tau_{k^{\prime}}}+\rho^{2}_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})},\\\
&<\int_{0}^{T_{s}}s^{2}\Big{(}t-(\tau_{k}-\tau_{k^{\prime}})\Big{)}dt\int_{0}^{T_{s}}s^{2}\Big{(}t\Big{)}dt\\\
&\ \ \
+\int_{0}^{T_{s}}s^{2}\Big{(}t+T_{s}-(\tau_{k}-\tau_{k^{\prime}})\Big{)}dt\int_{0}^{T_{s}}s^{2}\Big{(}t\Big{)}dt=1.\end{split}$
(39)

Therefore, we have
$\mathbb{E}\Big{(}|\bar{p}^{l}_{k,k^{\prime}}|^{2}\big{|}\text{Gaussian
pilot}\Big{)}<\mathbb{E}\big{(}|p^{l}_{k^{\prime}}|^{2}\big{)}$. In this way,
both conclusions from (38) and (39) corroborate the claim of
$\mathbb{E}\Big{(}\text{SINR}^{\text{asy}}_{k}\Big{)}>\mathbb{E}\Big{(}\text{SINR}^{\text{syn}}_{k}\Big{)}$
in Proposition 1. $\hfill\blacksquare$

For Gaussian pilots, according to (39), the factor
$\rho^{2}_{\tau_{k}-\tau_{k^{\prime}}}+\rho^{2}_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}$
determines the multi-user interference (MUI) from the $k^{\prime}$-th UE to
the $k$-th UE during the update of $I^{t}_{f^{m,k}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})$. Given the signal waveform $s(t)$, this MUI factor is a
function of the relative delay $|\tau_{k}-\tau_{k^{\prime}}|$ between the
$k^{\prime}$-th UE and the $k$-th UE, and we illustrate this function in Fig.
2(c). For two arbitrary UEs $k$ and $k^{\prime}$, we assume that $\tau_{k}$
and $\tau_{k^{\prime}}$ are uniformly and randomly distributed on the interval
$[0,T_{s}]$. Then, the average of the MUI factor is further plotted as the
dashed line in Fig. 2(c). In contrast to the asynchronous case, the
synchronous case has a larger MUI factor, i.e.
$\mathbb{E}\big{(}|p^{l}_{k^{\prime}}|^{2}\big{)}=1$. As a result, the message
$I^{t}_{f^{m,k}_{\delta(l)}\to h^{m}_{k}}(h^{m}_{k})$ under the asynchronous
case will be more reliable than the message $I^{t}_{f^{m}_{\delta(l)}\to
h^{m}_{k}}(h^{m}_{k})$ under the synchronous case. Furthermore, according to
(14) and (15), the belief of the VN $h^{m}_{k}$ will also be more reliable
under the asynchronous case. In summary, the proposed method exploits
asynchronization between different UEs to suppress MUI, leading to improved
detection accuracy.

###  III-F Asynchronous Case with Imperfect Delay Information

So far, we have derived and analyzed the PUDMP-SBL-aUADCE algorithm based on
the assumption that the BS has _perfect_ knowledge of the delay information
$\\{\tau_{k},k\in\mathcal{K}\\}$. In this subsection, we consider a more
practical scenario with imperfect delay information, which may be potentially
caused by device mobility or timing error of the low-cost oscillators at the
devices. Specifically, assume that the $k$-th UE is activated, and denote
$\tilde{\tau}_{k}$ as its _actual_ transmission delay in this round of GF-RA.
On the other hand, the BS only knows the imperfect delay information
$\tau_{k}$, and the delay error is denoted as
$\tilde{\tau}^{\text{err}}_{k}\overset{\Delta}{=}\tilde{\tau}_{k}-\tau_{k}$.
In the following, we will analyze the impacts of the imperfect delay
information, as well as the benefits of adopting the SBL method in this case.

Similar to (1), denote $\tilde{y}^{m}(t)$ as the received baseband waveform
signal with actual transmission delay
$\\{\tilde{\tau}_{k},k\in\mathcal{K}\\}$, and we have

$\tilde{y}^{m}(t)=\sum\limits^{K}_{k=1}\sum\limits^{L_{p}}_{l=1}\alpha_{k}g^{m}_{k}p^{l}_{k}s\Big{(}t-(l-1)T_{s}-\tilde{\tau}_{k}\Big{)}+\tilde{n}(t).$
(40)

Then, the BS will perform the matched filtering and sampling operations on
$\tilde{y}^{m}(t)$, using the imperfect delay information
$\\{\tau_{k},k\in\mathcal{K}\\}$. In this case, similar to (2), the $l$-th
sampled output $\tilde{y}^{m}_{k}(l)$ of the $k$-th matched filter on the
$m$-th antenna is written as

$\begin{split}&\tilde{y}^{m}_{k}(l)=\int^{lT_{s}+\tau_{k}}_{(l-1)T_{s}+\tau_{k}}\tilde{y}^{m}(t)s\Big{(}t-(l-1)T_{s}-\tau_{k}\Big{)}dt,\\\
&\overset{(e)}{=}\rho_{|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l}_{k}+\rho_{T_{s}-|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}+\tilde{n}^{m}_{k}(l)\\\
&\ \ \
+\sum\limits_{\\{k^{\prime},2T_{s}\geq\tilde{\tau}_{k^{\prime}}-\tau_{k}>T_{s}\\}}\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!h^{m}_{k^{\prime}}\Big{(}p^{l-2}_{k^{\prime}}\rho_{2T_{s}+\tau_{k}-\tilde{\tau}_{k^{\prime}}}+p^{l-1}_{k^{\prime}}\rho_{\tilde{\tau}_{k^{\prime}}-\tau_{k}-T_{s}}\Big{)}\\\
&\ \ \
+\sum\limits_{\\{k^{\prime},T_{s}\geq\tilde{\tau}_{k^{\prime}}-\tau_{k}>0\\}}\\!\\!\\!\\!\\!\\!\\!\\!h^{m}_{k^{\prime}}\Big{(}p^{l-1}_{k^{\prime}}\rho_{T_{s}+\tau_{k}-\tilde{\tau}_{k^{\prime}}}+p^{l}_{k^{\prime}}\rho_{\tilde{\tau}_{k^{\prime}}-\tau_{k}}\Big{)}\\\
&\ \ \
+\sum\limits_{\\{k^{\prime},T_{s}>\tau_{k}-\tilde{\tau}_{k^{\prime}}\geq
0\\}}\\!\\!\\!\\!\\!\\!\\!\\!h^{m}_{k^{\prime}}\Big{(}p^{l}_{k^{\prime}}\rho_{\tau_{k}-\tilde{\tau}_{k^{\prime}}}+p^{l+1}_{k^{\prime}}\rho_{T_{s}+\tilde{\tau}_{k^{\prime}}-\tau_{k}}\Big{)}\\\
&\ \ \
+\sum\limits_{\\{k^{\prime},2T_{s}>\tau_{k}-\tilde{\tau}_{k^{\prime}}\geq
T_{s}\\}}\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!h^{m}_{k^{\prime}}\Big{(}p^{l+1}_{k^{\prime}}\rho_{\tau_{k}-\tilde{\tau}_{k^{\prime}}-T_{s}}+p^{l+2}_{k^{\prime}}\rho_{2T_{s}+\tilde{\tau}_{k^{\prime}}-\tau_{k}}\Big{)},\\\
&\overset{(f)}{=}\rho_{|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l}_{k}+\rho_{T_{s}-|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}+\tilde{n}^{m}_{k}(l)+\tilde{\Theta}^{m}_{k}(l),\\\
&\overset{(g)}{=}\rho_{|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l}_{k}+\sum\limits_{k^{\prime}\neq
k}h^{m}_{k^{\prime}}\bar{p}^{l}_{k,k^{\prime}}+\tilde{w}^{m}_{k}(l),\\\
\end{split}$ (41)

where $h^{m}_{k}=\alpha_{k}g^{m}_{k}$, $\tilde{n}^{m}_{k}(l)$ is the AWGN,
equation ($e$) in (41) is derived based on the assumption that
$|\tilde{\tau}^{\text{err}}_{k}|\leq\frac{1}{2}T_{s}$ for $\forall
k\in\mathcal{K}$, and the sign function is $\text{sign}(x)=-1,0,1$ for $x<0$,
$x=0$, and $x>0$, respectively. With some straightforward derivation, we can
readily extend equation ($e$) in (41) to a more general case with
$|\tilde{\tau}^{\text{err}}_{k}|>\frac{1}{2}T_{s}$. However, we will show that
$|\tilde{\tau}^{\text{err}}_{k}|>\frac{1}{2}T_{s}$ indicates a serious
mismatch problem between the known delay information $\tau_{k}$ and the actual
delay $\tilde{\tau}_{k}$, which causes detection failure for the $k$-th UE.

To start with, we consider the first term on the right hand side (RHS) of
equation ($e$) in (41). In the first term, the same unknown coefficient
$\rho_{|\tilde{\tau}^{\text{err}}_{k}|}$ affects the VNs
$\\{h^{m}_{k},m\in\mathcal{M}\\}$ of the $k$-th UE on all the $M$ antennas.
Note that if the $k$-th UE is active, we have $h^{m}_{k}\sim\mathcal{CN}(0,1)$
and
$\rho_{|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}\sim\mathcal{CN}(0,\rho^{2}_{|\tilde{\tau}^{\text{err}}_{k}|})$
for $\forall m\in\mathcal{M}$. On the other hand, we adopt the conditional
prior distribution
$p(h^{m}_{k}|\gamma_{k})=\mathcal{CN}(h^{m}_{k}|0,\gamma^{-1}_{k})$ with
$p(\gamma_{k})=Gam(\gamma_{k};\epsilon,\eta)$ in the SBL method. Both
parameters $\gamma_{k}$ and $\epsilon$ are iteratively updated in the proposed
PUDMP-SBL-aUADCE algorithm. In other words, the PUDMP-SBL-aUADCE algorithm can
implicitly learn the unknown power factor
$\rho^{2}_{|\tilde{\tau}^{\text{err}}_{k}|}$ caused by delay error via
estimating parameter $\gamma_{k}$. In addition, as pointed out in Remark 1,
the shape parameter $\epsilon$ functions as a selective amplifier for the VNs
$\\{\gamma_{k},k\in\mathcal{K}\\}$ [32], and updating $\epsilon$ benefits the
robustness of the PUDMP-SBL-aUADCE algorithm against unknown delay error
$\\{\tilde{\tau}^{\text{err}}_{k},k\in\mathcal{K}\\}$.

Then, we consider the second term
$\rho_{T_{s}-|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}$
on the RHS of equation ($e$) in (41). Note that even if the factor
$\rho_{|\tilde{\tau}^{\text{err}}_{k}|}$ is accurately learned from the PUDMP-
SBL-aUADCE algorithm, it is still difficult to detect
$\text{sign}(\tilde{\tau}^{\text{err}}_{k})$ and
$p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}$. Therefore, a practical
yet effective solution is to treat the second term
$\rho_{T_{s}-|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}$
as a noise term. Furthermore, if
$|\tilde{\tau}^{\text{err}}_{k}|>\frac{1}{2}T_{s}$, the power of this noise
term
$\rho_{T_{s}-|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}$
will exceed the power of the first term
$\rho_{|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l}_{k}$, which may lead to
failed estimation of $\rho_{|\tilde{\tau}^{\text{err}}_{k}|}$ and $h^{m}_{k}$.
In other words, if $|\tilde{\tau}^{\text{err}}_{k}|>\frac{1}{2}T_{s}$, the
matched filter aligned with imperfect delay $\tau_{k}$ will fail to match the
target pilot symbol $p^{l}_{k}$, which seriously undermines the accuracy of
the PUDMP-SBL-aUADCE algorithm.

Finally, we consider the last four accumulative-addition terms on the RHS of
equation ($e$). For notation convenience, the sum of these four terms is
denoted as $\tilde{\Theta}^{m}_{k}(l)$ in equation ($f$) of (41), and we will
have $\tilde{\Theta}^{m}_{k}(l)=\sum\limits_{k^{\prime}\neq
k}h^{m}_{k^{\prime}}\bar{p}^{l}_{k,k^{\prime}}$ if
$\tilde{\tau}^{\text{err}}_{k^{\prime}}=0$ for $\forall
k^{\prime}\in\mathcal{K}/k$. However, due to the existence of the random delay
error $\\{\tilde{\tau}^{\text{err}}_{k},k\in\mathcal{K}\\}$, it is difficult
to determine $\tilde{\Theta}^{m}_{k}(l)$. Again, we treat the difference
between $\tilde{\Theta}^{m}_{k}(l)$ and $\sum\limits_{k^{\prime}\neq
k}h^{m}_{k^{\prime}}\bar{p}^{l}_{k,k^{\prime}}$ as an equivalent noise. In
this way, we can write $\tilde{y}^{m}_{k}(l)$ as in equation ($g$) of (41),
where $\tilde{w}^{m}_{k}(l)$ collects all the _equivalent_ noise terms, i.e.
$\tilde{w}^{m}_{k}(l)=\rho_{T_{s}-|\tilde{\tau}^{\text{err}}_{k}|}h^{m}_{k}p^{l-\text{sign}(\tilde{\tau}^{\text{err}}_{k})}_{k}+\tilde{n}^{m}_{k}(l)+\big{(}\tilde{\Theta}^{m}_{k}(l)-\sum\limits_{k^{\prime}\neq
k}h^{m}_{k^{\prime}}\bar{p}^{l}_{k,k^{\prime}}\big{)}$. In this way, equation
($g$) is employed by the proposed PUDMP-SBL-aUADCE algorithm for UAD and CE,
where the noise power of $\tilde{w}^{m}_{k}(l)$ is estimated by updating
$\hat{\lambda}_{k}$ in (34).

##  IV Simulation Results

Figure 6: Convergence performance on (a) AER and (b) CE-MSE with $K=200$,
$p_{a}=0.1$, and $\text{SNR}=10\text{dB}$. Figure 7: Impacts of pilot length
$L_{p}$ on (a) AER and (b) CE-MSE performance with $K=200$, $M=10$,
$p_{a}=0.2$, and $N_{it}=80$.

In this section, extensive simulations are conducted to examine the
performances of the PUDMP-SBL-aUADCE algorithm. Specifically, we consider the
UAD and CE accuracy, and define the activity-detection error rate (AER) and
channel-estimation mean square error (CE-MSE) as performance metrics, i.e.,
$\text{AER}={\sum^{K}_{k=1}|\alpha_{k}-\hat{\alpha}_{k}|}/{K}$ and $\text{CE-
MSE}={\sum^{K}_{k=1}\sum^{M}_{m=1}|h^{m}_{k}-\hat{h}^{m}_{k}|^{2}}/{KM}$.
Through the following simulations, we investigate the convergence performance,
and the impacts of different system parameters, i.e. the pilot length $L_{p}$,
the antenna number $M$, the activation probability $p_{a}$, the delay error
$\tilde{\tau}^{\text{err}}$, and the signal-to-noise ratio (SNR) which is
defined as $\text{SNR}=10\log_{10}\frac{1}{\sigma^{2}_{n}}$. For comparison,
we also simulate the performances of the unitary AMP based SBL (UAMP-SBL)
algorithm [31, 32], the block orthogonal matching pursuit (BOMP) algorithm
[33], the orthogonal AMP (OAMP) algorithm [34], the genie-aided minimum mean
square error (GA-MMSE) algorithm, as well as the MP-BSBL algorithm [5] which
can be considered as a synchronous version of the PUDMP-SBL-aUADCE algorithm.
Here, the GA-MMSE algorithm is aided with ideal user activity information, and
therefore it provides the CE-MSE lower bound. We will further show that the
proposed PUDMP-SBL-aUADCE algorithm can closely approach this lower bound
within a wide range of SNR.

Figure 8: Impacts of antenna number $M$ on (a) AER and (b) CE-MSE performance
with $K=200$, $Lp=20$, $p_{a}=0.1$, and $N_{it}=80$. Figure 9: Impacts of
activation probability $p_{a}$ on (a) AER and (b) CE-MSE performance with
$K=200$, $M=10$, $\text{SNR}=10\text{dB}$, and $N_{it}=80$.

###  IV-A Convergence Performance

To start with, we simulate the convergence performance of the PUDMP-SBL-aUADCE
algorithm with different iteration number $N_{it}$. The AER and CE-MSE results
are illustrated in Fig. 6(a) and Fig. 6(b), respectively. It is shown that
with the same pilot length $L_{p}=30$ and antenna number $M=10$, the PUDMP-
SBL-aUADCE algorithm in the asynchronous case could significantly outperform
the MP-BSBL algorithm in the synchronous case, especially in terms of the AER
performance. This performance gain from asynchronization is consistent with
the conclusion from Proposition 1. Consider only the asynchronous case and fix
$M=10$, it is shown that increasing $L_{p}$ from $20$ to $30$ could improve
the convergence speed of the PUDMP-SBL-aUADCE algorithm, as well as the UAD
and CE accuracy. Fixing $L_{p}=20$ and increasing $M$ from $10$ to $30$ for
the PUDMP-SBL-aUADCE algorithm, we can also observe a performance gain for the
AER, while the CE-MSE performance will not get improved with larger $M$. This
observation can be explained from the fact that the channel gains
$\\{h^{m}_{k},m\in\mathcal{M}\\}$ for each active UE $K$ are _i.i.d._ on
different antennas. In this case, as long as the UAD accuracy is sufficiently
high, increasing the antenna number $M$ will not lead to higher CE accuracy.

###  IV-B Impacts of Pilot Length

As shown in Fig. 6, the pilot length has remarkable impacts on the convergence
and accuracy of the PUDMP-SBL-aUADCE algorithm, and simulation results with
different pilot length are plotted in Fig. 7. We can also observe the
asynchronization gain from Fig. 7, i.e. with the same
$\text{SNR}=10\text{dB}$, the PUDMP-SBL-aUADCE algorithm in the asynchronous
case still outperforms the MP-BSBL algorithm in the synchronous case. For
example, fix $\text{SNR}=10\text{dB}$ and $L_{p}=40$, the AER of the PUDMP-
SBL-aUADCE algorithm is lower than that of the MP-BSBL algorithm by _more than
two magnitude orders_. In addition, it is shown that increasing the pilot
length could also effectively improve the CE-MSE performances.

Figure 10: Impacts of $\sigma_{\tau}$ on (a) AER and (b) CE-MSE performance
with $K=200$, $M=10$, $\text{SNR}=10\text{dB}$, and $N_{it}=80$.

###  IV-C Impacts of Antenna Number

The detection performances with different antenna number $M$ are shown in Fig.
8. Again, we can still observe the asynchronization gain when SNR is fixed as
10 dB. Furthermore, it is shown that increasing the antenna number $M$ can
always improve the UAD detection accuracy. This observation is explained by
the fact that with more antennas, more messages
$\\{\mu^{t}_{h^{m}_{k}},v^{t}_{h^{m}_{k}},m\in\mathcal{M}\\}$ will be employed
in (22) to update the activity-related hyper-parameter $\hat{\gamma}_{k}$ for
each UE $k$. However, the CE-MSE performance converges at $M=16$, and further
increasing $M$ will not improve the CE accuracy. Explanation for this
observation can be found in Section IV-A.

###  IV-D Impacts of Activation Probability

The impacts of the activation probability $p_{a}$ on the PUDMP-SBL-aUADCE
algorithm are illustrated in Fig. 9. It is shown that, when the random access
scenario is crowded with more active UEs, both the AER and the CE-MSE
performances will inevitably deteriorate due to the increased MUI. As
explained in Proposition 1, the PUDMP-SBL-aUADCE algorithm can effectively
mitigate the MUI in the asynchronous case. Therefore, given the same pilot
length $L_{p}$, the PUDMP-SBL-aUADCE algorithm can always outperform the MP-
BSBL algorithm within a wide range of $p_{a}$. Especially, in the low $p_{a}$
regime, the PUDMP-SBL-aUADCE algorithm can improve the AER of the MP-BSBL
algorithm by _more than two magnitude orders_ , which makes the PUDMP-SBL-
aUADCE algorithm particularly favorable to mMTC applications that are featured
by sparsely activated devices.

###  IV-E Impacts of Imperfect Delay Information

We consider a more practical scenario with imperfect delay information
potentially caused by device mobility or timing error. Recall that the delay
error $\tilde{\tau}^{\text{err}}_{k}$ is defined as the difference between the
actual transmission delay $\tilde{\tau}_{k}$ and the imperfect delay
information $\tau_{k}$ at the BS, i.e.
$\tilde{\tau}^{\text{err}}_{k}=\tilde{\tau}_{k}-\tau_{k}$. Without loss of any
generality, we assume that
$\tilde{\tau}^{\text{err}}_{k}\sim\mathcal{N}\big{(}0,(\sigma_{\tau}T_{s})^{2}\big{)}$
for $\forall k\in\mathcal{K}$. Note that the assumption of Gaussian-
distributed $\tilde{\tau}^{\text{err}}_{k}$ is adopted simply for simulation
convenience, so that the severity of the delay error can be controlled with
one parameter $\sigma_{\tau}$. The performances of the PUDMP-SBL-aUADCE
algorithm under different $\sigma_{\tau}$ are illustrated in Fig. 10. It is
shown that, even in the presence of random delay error
$\tilde{\tau}^{\text{err}}_{k}$, the PUDMP-SBL-aUADCE algorithm still exhibits
higher UAD and CE accuracy than the MP-BSBL algorithm. In addition, it is
shown that both algorithms will fail to work when $\sigma_{\tau}$ approaches
$0.2$. This observation can be explained by the fact that a larger
$\sigma_{\tau}$ is more likely to incur the mismatch problem with
$|\tilde{\tau}^{\text{err}}_{k}|>\frac{1}{2}T_{s}$. As explained in Section
III-F, if $|\tilde{\tau}^{\text{err}}_{k}|>\frac{1}{2}T_{s}$, the matched
filter aligned with imperfect delay $\tau_{k}$ will fail to match the target
pilot symbol, which eventually leads to failed UAD and CE.

Figure 11: CE-MSE performance under different SNR with $K=200$, $M=32$,
$L_{p}=30$, and $p_{a}=0.1$. Solid lines represent the performances of
different algorithms in the asynchronous case, while dashed lines represent
the performances of different algorithms in the synchronous case.

###  IV-F Performance Comparison under Different SNR

For comparison convenience, all the sampled outputs on different antennas are
re-organized in the linear mixing model,

$\bar{\bm{Y}}=\bar{\bm{P}}\bm{H}+\bar{\bm{N}},$ (42)

where
$\bar{\bm{Y}}=[\bm{Y}_{1}^{T},\ldots,\bm{Y}_{k}^{T},\ldots,\bm{Y}_{K}^{T}]^{T}$
has the size $KL_{p}\times M$,
$\bm{Y}_{k}=[\bm{y}^{1}_{k},\ldots,\bm{y}^{m}_{k},\ldots,\bm{y}^{M}_{k}]$ has
the size $L_{p}\times M$, and $\bm{Y}_{k}$ collects the sampled outputs from
the $k$-th matched filter on all the $M$ antennas, while $\bm{y}^{m}_{k}$ is
given in (6).
$\bar{\bm{P}}=[\bar{\bm{P}}^{T}_{1},\ldots,\bar{\bm{P}}^{T}_{k},\ldots,\bar{\bm{P}}^{T}_{K}]^{T}$,
and the element $\bar{p}^{l}_{k,k^{\prime}}$ on the $l$-th row and
$k^{\prime}$-th column of $\bar{\bm{P}}_{k}$ is given in (3). $\bm{H}$ is the
effective channel matrix explained in Section III-A. Each column
$\bar{\bm{n}}^{m},m\in\mathcal{M}$ of $\bar{\bm{N}}$ is distributed as
$\bar{\bm{n}}^{m}\sim\mathcal{CN}(\bm{0},\bar{\bm{\Sigma}})$ with the noise
covariance matrix $\bar{\bm{\Sigma}}$ written as
$\bar{\bm{\Sigma}}=\bm{\Sigma}\sigma^{2}_{n}$. Denote $\bm{\Sigma}_{i,j}$ as
the element on the $i$-th row and $j$-th column of $\bm{\Sigma}$. According to
(9), $\bm{\Sigma}$ is a symmetric matrix with
$\bm{\Sigma}_{(k-1)L_{p}+l,(k^{\prime}-1)L_{p}+l}=\rho_{\tau_{k}-\tau_{k^{\prime}}}$
and
$\bm{\Sigma}_{(k-1)L_{p}+l,(k^{\prime}-1)L_{p}+l+1}=\rho_{T_{s}-(\tau_{k}-\tau_{k^{\prime}})}$
for $\forall k^{\prime}\leq k$. In this way, we can adopt the GA-MMSE
algorithm, a modified version of the UAMP-SBL algorithm [32], the BOMP
algorithm [33], and the OAMP algorithm [34] for comparison. The active-UE
number $K_{a}$ is assumed known to the BOMP algorithm, which provides an
_ideal stopping criterion_ for the matching pursuit procedure in the BOMP
algorithm. As the UAMP-SBL algorithm developed in [32] is for white noise, a
noise-whitening operation is incorporated by left-multiplying
$\bm{\Sigma}^{-1/2}$ on both sides of (42), which is called as UAMP-SBL-
Whitening. In addition, the GA-MMSE algorithm is aided with ideal user
activity information, so that it provides the CE-MSE lower bound. The
simulation results are illustrated in Fig. 11.

It is shown in Fig. 11 that, in the asynchronous case, both the GA-MMSE
algorithm and the BOMP algorithm outperform their counterparts in the
synchronous case. Similarly, the PUDMP-SBL-aUADCE algorithm also outperforms
the MP-BSBL algorithm. The UAMP-SBL-Whitening algorithm also outperforms the
algorithms in the synchronous case. However, the computational complexity of
this noise-whitening operation is $\mathcal{O}(K^{3}L_{p}^{3})$. In contrast,
the proposed receiver design employs the PUD factor graph with related message
update rules to address this sample correlation problem. It is shown that the
proposed PUDMP-SBL-aUADCE algorithm can closely approach the CE-MSE lower
bound within a wide range of SNR. Similarly, the correlation feature of the
noise samples in $\bar{\bm{N}}$ is also exploited by the Linear MMSE (LMMSE)
module of the OAMP algorithm, so that the OAMP algorithm also exhibits the
bound-approaching performances. However, this LMMSE module incurs
$\mathcal{O}(K^{3}L^{3}_{p})$ computations for the OAMP algorithm, which is
prohibitively high for mMTC applications with a massive UE number $K$. Recall
from Section III-D that the computations of the PUDMP-SBL-aUADCE algorithm is
inherently feasible for parallel implementations. Therefore, the PUDMP-SBL-
aUADCE algorithm serves as a bound-approaching yet practical solution for the
joint UAD and CE problem in asynchronous GF-RA.

##  V Conclusions

In this paper, a PUDMP-SBL-aUADCE algorithm was proposed to address the joint
UAD and CE problem for asynchronous GF-RA. Theoretical analysis shows that the
PUDMP-SBL-aUADCE algorithm provides higher SINR in the asynchronous case than
in the synchronous case. Considering the potential timing error from the low-
cost UEs, the impacts of imperfect delay profile were investigated, revealing
the advantages of adopting the SBL method in this case. Finally, extensive
simulation results were provided to demonstrate the superior performances of
the PUDMP-SBL-aUADCE algorithm.

## References

  * [1] W. Zhan, C. Xu, X. Sun, and J. Zou, “Toward optimal connection management for massive machine-type communications in 5G system,” _IEEE Internet Things J.,_ vol. 8, no. 17, pp. 13237-13250, Sept. 2021.
  * [2] J. Choi, J. Ding, N.-P. Le, and Z. Ding, “Grant-free random access in machine-type communication: Approaches and challenges,” _IEEE Wireless Commun.,_ vol. 29, no. 1, pp. 151-158, Feb. 2022.
  * [3] L. Liu and W. Yu, “Massive connectivity with massive MIMO-Part I: Device activity detection and channel estimation,” _IEEE Trans. Signal Process.,_ vol. 66, no. 11, pp. 2933-2946, Mar. 2018.
  * [4] Z. Zhang _et al._ , “User activity detection and channel estimation for grant-free random access in LEO satellite-enabled Internet of Things,” _IEEE Internet Things J.,_ vol. 7, no. 9, pp. 8811-8825, Sept. 2020.
  * [5] Y. Zhang, Q. Guo, Z. Wang, J. Xi, and N. Wu, “Block sparse Bayesian learning based joint user activity detection and channel estimation for grant-free NOMA systems,” _IEEE Trans. Veh. Technol.,_ vol. 67, no. 10, pp. 9631-9640, Oct. 2018.
  * [6] Z. Zhang, Y. Li, C. Huang, Q. Guo, C. Yuen, and Y. L. Guan, “DNN-aided block sparse Bayesian learning for user activity detection and channel estimation in grant-free non-orthogonal random access,” _IEEE Trans. Veh. Technol.,_ vol. 68, no. 12, pp. 12000-12012, Dec. 2019.
  * [7] X. Bian, Y. Mao, and J. Zhang, “Joint activity detection, channel estimation, and data decoding for grant-free massive random access,” _IEEE Internet Things J.,_ Early Access, 2023.
  * [8] S. Jiang, X. Yuan, X. Wang, C. Xu, and W. Yu, “Joint user identification, channel estimation, and signal detection for grant-free NOMA,” _IEEE Trans. Wireless Commun.,_ vol. 19, no. 10, pp. 6960-6976, Oct. 2020.
  * [9] Z. Zhang, Q. Guo, Y. Li, M. Jin, and C. Huang, “Variational Bayesian inference clustering based joint user activity and data detection for grant-free random access in mMTC,” _IEEE Internet Things J.,_ Early Access, 2023.
  * [10] S. Kim, J. Kim, and D. Hong, “A new non-orthogonal transceiver for asynchronous grant-free transmission systems,” _IEEE Trans. Wireless Commun.,_ vol. 20, no. 3, pp. 1889-1902, Mar. 2021.
  * [11] G. Sun, Y. Li, X. Yi, W. Wang, X. Gao, and L. Wang, “OFDMA based massive grant-free transmission in the presence of timing offset,” in _Proc. 13th Int. Conf. Wireless Commun. Signal Process.,_ Changsha, China, 2021, pp. 1-6.
  * [12] K. Andreev, S. S. Kowshik, A. Frolov, and Y. Polyanskiy, “Low complexity energy efficient random access scheme for the asynchronous fading MAC,” in _Proc. IEEE 90th Veh. Technol. Conf. (VTC-Fall),_ Honolulu, HI, USA, 2019, pp. 1-5.
  * [13] X. Chen, L. Liu, D. Guo, and G. W. Wornell, “Asynchronous massive access and neighbor discovery using OFDMA,” _IEEE Trans. Inf. Theory,_ vol. 69, no. 4, pp. 2364-2384, Apr. 2023.
  * [14] H. Haci, H. Zhu, and J. Wang, “Performance of non-orthogonal multiple access with a novel asynchronous interference cancellation technique,” _IEEE Trans. Commun.,_ vol. 65, no. 3, pp. 1319-1335, Mar. 2017.
  * [15] S. Kim, H. Kim, H. Noh, Y. Kim and, D. Hong, “Novel transceiver architecture for an asynchronous grant-free IDMA system,” _IEEE Trans. Wireless Commun.,_ vol. 18, no. 9, pp. 4491-4504, Sept. 2019.
  * [16] J. Liu, Y. Li, G. Song, and Y. Sun, “Detection and analysis of symbol-asynchronous uplink NOMA with equal transmission power,” _IEEE Wireless Commun. Lett.,_ vol. 8, no. 4, pp. 1069-1072, Aug. 2019.
  * [17] R. Zhang, X. Wu, C. Li, J. Yan, and S. Zhang, “Active user detection and data recovery for quasi-asynchronous grant-free random access,” in _Proc. 13th Int. Conf. Wireless Commun. Signal Process.,_ Changsha, China, 2021, pp. 1-5.
  * [18] H. F. Schepker, C. Bockelmann, and A. Dekorsy, “Coping with CDMA asynchronicity in compressive sensing multi-user detection,” in _Proc. IEEE 77th Veh. Technol. Con. (VTC-Spring),_ Dresden, Germany, 2013, pp. 1-5.
  * [19] X. Lin, L. Kuang, Z. Ni, C. Jiang, and S. Wu, “Approximate message passing-based detection for asynchronous NOMA,” _IEEE Commun. Lett.,_ vol. 24, no. 3, pp. 534-538, Mar. 2020.
  * [20] W. Zhu, M. Tao, X. Yuan, and Y. Guan, “Asynchronous massive connectivity with deep-learned approximate message passing,” in _Proc. IEEE Int. Conf. Commun.,_ Dublin, Ireland, Jun. 2020, pp. 1-6.
  * [21] W. Zhang, J. Li, X. Zhang, and S. Zhou, “A joint user activity detection and channel estimation scheme for packet-asynchronous grant-free access,” _IEEE Wireless Commun. Lett.,_ vol. 11, no. 2, pp. 338-342, Feb. 2022.
  * [22] T. Ding, X. Yuan, and S. C. Liew, “Sparsity learning-based multiuser detection in grant-free massive-device multiple access,” _IEEE Trans. Wireless Commun.,_ vol. 18, no. 7, pp. 3569-3582, Jul. 2019.
  * [23] L. Liu and Y. Liu, “An efficient algorithm for device detection and channel estimation in asynchronous IoT systems,” in _Proc. IEEE ICASSP,_ Toronto, ON, Canada, Jun. 2021, pp. 4815-4819.
  * [24] W. Zhu, M. Tao, X. Yuan, and Y. Guan, “Deep-learned approximate message passing for asynchronous massive connectivity,” _IEEE Trans. Wireless Commun.,_ vol. 20, no. 8, pp. 5434-5448, Aug. 2021.
  * [25] R. B. D. Renna and R. C. de Lamare, “Dynamic message scheduling based on activity-aware residual belief propagation for asynchronous mMTC,” _IEEE Wireless Commun. Lett.,_ vol. 10, no. 6, pp. 1290-1294, Jun. 2021.
  * [26] S. Jiang, C. Xu, X. Yuan, Z. Han, Z. Wang, and X. Wang, “Bayesian receiver design for asynchronous massive connectivity,” _IEEE Trans. Wireless Commun.,_ vol. 22, no. 1, pp. 301-316, Jan. 2023.
  * [27] _Machine to Machine (M2M) Communication Study Report (Draft)_ , IEEE P802.16 PPC, Nov. 2010.
  * [28] 3GPP TR 22.368, _Service Requirements for Machine-Type Communications (MTC),_ v.11.1.0.
  * [29] K. S. Ko _et al._ , “A novel random access for fixed-location machine-to-machine communications in OFDMA based systems,” _IEEE Commun. Lett.,_ vol. 16, no. 9, pp. 1428-1431, Sep. 2012.
  * [30] Z. Wang and V. W. S. Wong, “Optimal access class barring for stationary machine type communication devices with timing advance information,” _IEEE Trans. Wireless Commun.,_ vol. 14, no. 10, pp. 5374-5387, Oct. 2015.
  * [31] Q. Guo and J. Xi, “Approximate message passing with unitary transformation,” Apr. 2015, _arXiv:1504.04799._ [Online]. Available: http://arxiv.org/abs/1504.04799
  * [32] M. Luo, Q. Guo, M. Jin, Y. C. Eldar, D. Huang, and X. Meng, “Unitary approximate message passing for sparse Bayesian learning,” _IEEE Trans. Signal Process.,_ vol. 69, pp. 6023-6039, 2021.
  * [33] Y. C. Eldar, P. Kuppinger, and H. Bolcskei, “Block-sparse signals: Uncertainty relations and efficient recovery,” _IEEE Trans. Signal Process.,_ vol. 58, no. 6, pp. 3042-3054, Jun. 2010.
  * [34] Y. Cheng, L. Liu, and L. Ping, “Orthogonal AMP for massive access in channels with spatial and temporal correlations,” _IEEE J. Sel. Areas Commun.,_ vol. 39, no. 3, pp. 726-740, Mar. 2021.