# Cascaded Channel Estimation for
Intelligent Reflecting Surface Assisted
Multiuser MISO Systems

Huayan Guo, , and Vincent K. N. Lau  The work of H. Guo and V. K. N. Lau was
supported by the Hong Kong Research Grants Council Collaborative Research Fund
2020/21 Collaborative Research Project Grant (CRPG) C6012-20GF. The authors
are with the Department of Electronics and Computer Engineering, The Hong Kong
University of Science and Technology, Hong Kong 999077, China (e-mail:
{eeguohuayan, eeknlau}@ust.hk).

###### Abstract

This paper investigates the uplink cascaded channel estimation for
intelligent-reflecting-surface (IRS)-assisted multi-user multiple-input-
single-output systems. We focus on a sub-6 GHz scenario where the channel
propagation is not sparse and the number of IRS elements can be larger than
the number of BS antennas. A novel channel estimation protocol without the
need of on-off amplitude control to avoid the reflection power loss is
proposed. In addition, the pilot overhead is substantially reduced by
exploiting the common-link structure to decompose the cascaded channel
coefficients by the multiplication of the common-link variables and the user-
specific variables. However, these two types of variables are highly coupled,
which makes them difficult to estimate. To address this issue, we formulate an
optimization-based joint channel estimation problem, which only utilizes the
covariance of the cascaded channel. Then, we design a low-complexity
alternating optimization algorithm with efficient initialization for the non-
convex optimization problem, which achieves a local optimum solution. To
further enhance the estimation accuracy, we propose a new formulation to
optimize the training phase shifting configuration for the proposed protocol,
and then solve it using the successive convex approximation algorithm.
Comprehensive simulations verify that the proposed algorithm has supreme
performance compared to various state-of-the-art baseline schemes.

###### Index Terms:

Intelligent reflecting surface (IRS), reconfigurable intelligent surface
(RIS), channel estimation, multiple-input multiple-output (MIMO).

##  I Introduction

Intelligent reflecting surface (IRS) is an artificial planar array consisting
of numerous reconfigurable passive elements with the capability of
manipulating the impinging electromagnetic signals and offering anomalous
reflections [1, 2, 3, 4]. Many recent studies have indicated that IRS is a
promising solution to build a programmable wireless environment via steering
the incident signal in fully customizable ways to enhance the spectral and
energy efficiency of legacy systems [5, 6, 7, 8]. Most contributions in this
area focus on joint active and passive procoding design with various
objectives and constraints [9, 10, 11, 12, 13]. The potential gains claimed by
these works highly depend on the availability of accurate channel state
information (CSI). However, channel estimation is a challenging task for the
IRS-assisted system because there are no sensing elements or radio frequency
chains, and thus there is no baseband processing capability in the IRS.

Figure 1: The IRS-assisted multiuser MISO communication system.

(a) Protocol for SU-MISO system

(b) Directly extended protocol for MU-MISO system

Figure 2: The conventional uplink channel estimation protocols without
exploiting the common-link structure.

Some early-attempted works [14, 15, 16] estimate the uplink cascaded IRS
channel for single user (SU) multiple-input-single-output (MISO) systems using
the protocol shown in Figure 2(a). In these works, the cascaded channel is
equivalently represented as a traditional $M\times N$ MIMO channel, where $M$
and $N$ is the base station (BS) array size and the IRS size, respectively,
and the sensing matrix for channel reconstruction consists of the phase
shifting vectors in consecutive training timeslots. Many works [17, 18, 19,
20, 21, 22] directly extend the SU protocol to multi-user (MU) MISO systems,
in which $K$ users transmit orthogonal pilot sequences in each training
timeslot, as shown in Figure 2(b). Based on this protocol, the on-off IRS
state (amplitude) control strategy is proposed in [17, 18, 19] to better
decompose the MU cascaded channel coefficients for easier channel estimation
of the cascaded channel for each user.

It is pointed out by [23] that direct application of the SU protocol on an MU-
MISO system fails to exploit the structural property and results in
substantially larger pilot overheads. Intuitively, all the cascaded channels
share a common BS-IRS link and it is possible to reduce the pilot overhead
since the number of independent variables is $MN+NK$ instead of $MNK$. One
algorithm is proposed in [24] with the idea of sequentially estimating the BS-
IRS channel and the IRS-user channels. However, it requires that the BS can
work at full-duplex mode. For the cascaded channel estimation, a new channel
estimation protocol is proposed in [23]. Specifically, the cascaded channel of
one reference user is firstly estimated based on the SU protocol, and then
other users’ channels are estimated by only estimating the ratios of their
channel coefficients to the reference channel, which can be referred to as the
relative channels. The overall training overhead is reduced from $NK$ to
$K+N+\lceil\frac{N}{M}\rceil(K-1)$. However, there is an error propagation
issue associated with this scheme since a low-accuracy estimation on the
reference channel may jeopardize the estimations of the relative channels.
Moreover, some IRS elements need to be switched off while estimating the
relative channels for coefficients decomposition [23].

For IRS design, the “off” state means no reflection (i.e., perfectly absorbing
the incident signals), and hence it is difficult [25, 26, 27] and also
attracts additional implementation costs since this state is unnecessary for
data transmission after the channel estimation. In addition, switching off the
IRS elements causes reflection power loss, which will lower the receive
signal-to-noise ratio (SNR). Some recent works attempt to overcome this issue
using “always-ON” training schemes. In [20, 21, 22], cascaded channel
estimation algorithms based on tensor decomposition are proposed for MU-MISO
systems, without requiring selected IRS elements to be off using the protocol
in Figure 2(b). In particular, the training phase shifts are optimized to
minimize the mean squared error (MSE), and it has been verified that the
discrete-Fourier-transform (DFT)-based training phase shifting configuration
is optimal in this scenario. However, the pilot overhead is $NK$ since the
protocol in Figure 2(b) does not utilize the common-link property. In [28], an
always-ON training scheme is proposed, which extends the protocol in [23] to
the double-IRS aided system. However, the number of BS antennas needs to be
equal to or larger than the number of IRS elements (i.e., $M\geq N$) to
guarantee a full-rank measurement matrix to estimate the relative channels
between the reference user and the other users.111 According to the property
${\rm{rank}}({\bf A}\otimes{\bf B})={\rm{rank}}({\bf A}){\rm{rank}}({\bf B})$,
the rank of the measurement matrix in equation (40) of [28] cannot be larger
than $(K-1)M$ while the targeted rank is $(K-1)N$.  This assumption is quite
restrictive as the number of elements of the IRS ($N$) is usually larger than
the number of antennas at the BS.

Another critical problem is the feasibility issue when the channel statistical
prior information is utilized to improve the channel estimation accuracy,
although this is a common idea for conventional MIMO channel estimation. In
[21] and [22], statistical knowledge of the individual BS-IRS link and IRS-
user links is required. However, these messages are not available in practice
since none of the existing algorithms, to the best of the authors’ knowledge,
can reconstruct the individual channel coefficients when $M<N$ due to the
ambiguity issue, which introduces a random scaling on the estimated individual
channels [17, 29]. In [23], the linear minimum mean squared error (LMMSE)
estimator is adopted to estimate the relative channels. However, the LMMSE
estimator requires the covariance of the relative channels, which is also
difficult since in general the distribution of the relative channels is heavy-
tailed.

In this paper, we focus on the uplink cascaded IRS channel estimation for an
MU-MISO system in sub-6 GHz bands. In these frequency bands, the IRS design
will be easier, but the channel propagation is not sparse.222If the IRS
operates at millimeter wave frequencies, some sparse channel estimation
algorithms have been further proposed based on the protocol in Figure 2(b)
with less training timeslots by exploiting the potential channel sparsity [29,
30, 31, 32, 33, 34].  In addition, the number of IRS elements in the system
can be larger than the number of BS antennas (i.e., $N>M$). We propose a
holistic solution to address the aforementioned issues. In particular, a novel
always-ON training protocol is designed; meanwhile the common-link structure
is utilized to reduce the pilot overhead. Furthermore, an optimization-based
cascaded channel estimation framework, which is flexible to utilize more
practical channel statistical prior information, is proposed. The following
summarizes our key contributions.

  * •

Always-ON Channel Estimation Protocol Exploiting the Common-Link Structure: We
propose a novel channel estimation protocol without the need for on-off
amplitude control to avoid the reflection power loss. Meanwhile, the common-
link structure is exploited and the pilot overhead is reduced to
$K+N+\lceil\frac{N}{M}\rceil(K-1)$.

In addition, the proposed protocol is applicable with any number of elements
at the IRS (also including $N\leq M$). Further, it does not need a “reference
user”, and as such, the estimation performance is enhanced owing to the
multiuser diversity.

  * •

Optimization-Based Cascaded Channel Estimation Framework: Since there is no
on-off amplitude control, the cascaded channel coefficients are highly
coupled. In order to exploit the common-link structure, we decompose the
cascaded channel coefficients by the multiplication of the common-link
variables and the user-specific variables, and then an optimization-based
joint channel estimation problem is formulated based on the maximal a
posterior probability (MAP) rule. The proposed optimization-based approach is
flexible to incorporate different kinds of channel statistical prior setups.
Specifically, we utilize the combined statistical information of the cascaded
channels, which is a weaker requirement compared to statistical knowledge of
the individual BS-IRS and IRS-user channels. Then, a low-complexity
alternating optimization algorithm is proposed to achieve a local optimum
solution. Simulation results demonstrated that the optimization solution with
proposed protocol achieves a more than $15$ dB gain compared to the benchmark.

  * •

Training Phase Shifting Optimization for the Proposed Protocol: The phase
shifting configuration can substantially enhance the channel estimation
performance of the cascaded IRS channel because the phase shifting vectors are
important components in the measurement matrix for channel reconstruction.
However, traditional solutions [20, 21, 22] of phase shifting optimization for
SU cascaded IRS channel estimation cannot be directly applied to the MU case
due to that the cascaded channel coefficients are highly coupled when the
common-link structure is exploited. We propose a new formulation to optimize
the phase shifting configuration, which maximizes the average reflection gain
of the IRS. Simulation results further verify the proposed configuration
achieves a more than $3$ dB gain compared to the state-of-the-art baselines.

##  II System Model

###  II-A System Model of MU-MISO IRS Systems

This paper investigates the uplink channel estimation in a narrow-band IRS-
aided MU-MISO communication system that consists of one BS with $M$ antennas,
one IRS with $N$ elements, and $K$ single-antenna users,333 We adopt the
single-antenna-user setup here for ease of presentation. The signal model can
be directly extended to the setup when users have multiple antennas by
transmitting orthogonal uplink pilot sequences in different antennas.  as
illustrated in Figure 1. Let ${\bf h}_{{\rm d},k}\in{\mathbb{C}}^{M\times 1}$
denote the BS-user channel (a.k.a., the direct channel) for user $k$, ${\bf
G}\in{\mathbb{C}}^{N\times M}$ denote the common BS-IRS channel, and ${\bf
h}_{{\rm r},k}\in{\mathbb{C}}^{N\times 1}$ denote the IRS-user channel for
user $k$. We assume quasi-static block fading for all the channels such that
the channel coefficients remain constant within one channel coherence
interval, and they are independent and identically distributed (i.i.d.)
between coherence intervals. Note that the quasi-static model considers the
worse case scenarios where the temporal correlations between blocks are not
exploited. In practice, the pilot overheads can be further reduced if one
exploits temporal correlations of the channel blocks [35, 36], but this is
outside the scope of the paper.

The received baseband signal at the BS is given by

$\displaystyle{\bf y}_{t}$ $\displaystyle=\sum_{k=1}^{K}\left({\bf h}_{{\rm
d},k}+{\bf G}^{\rm T}{\bf\Theta}_{t}{\bf h}_{{\rm r},k}\right)x_{k,t}+{{\bf
z}_{t}},$ (1)

where $t$ is the time index, $x_{k,t}$ is the transmit pilot symbol from user
$k$, ${\bf z}_{t}\sim{\cal{CN}}({\bm{0}},\sigma_{0}^{2}{\bf I}_{M})$ is the
additive white Gaussian noise (AWGN), and
${\bf\Theta}_{t}\in{\mathbb{C}}^{N\times N}$ is the IRS reflection coefficient
matrix. It is known that ${\bf\Theta}_{t}$ is a diagonal matrix such that
${\bf\Theta}_{t}={\rm diag}({\bm{\theta}}_{t})$, where
${\bm{\theta}}_{t}=[e^{\jmath\varphi_{t,1}},e^{\jmath\varphi_{t,2}},\cdots,e^{\jmath\varphi_{t,N}}]^{\rm
T}$ is the phase shifting vector from the IRS.444 In practice, the reflection
efficiency cannot be 1, which is known as the reflection loss. However, for
ease of presentation, this loss can be absorbed into the path loss of ${\bf
G}$ since it is a constant value.

###  II-B IRS Cascaded Channel Model

Denote the cascaded channel related to the BS’s $m$-th antenna and the $k$-th
user by

$\displaystyle{\bf h}_{{\rm I},k,m}={\rm{diag}}({\bf h}_{{\rm r},k}){\bf
g}_{m},$ (2)

where ${\bf g}_{m}$ is the $m$-th column in ${\bf G}=[{\bf g}_{1},\cdots,{\bf
g}_{M}]$. The cascaded channel over all BS antennas is given by ${\bf H}_{{\rm
I},k}=[{\bf h}_{{\rm I},k,1}^{\rm T},\cdots,{\bf h}_{{\rm I},k,M}^{\rm
T}]^{\rm T}$, and we have

$\displaystyle{\bf H}_{{\rm I},k}={\bf G}^{\rm T}{\rm{diag}}({\bf h}_{{\rm
r},k}).$ (3)

Substituting (3) into (1), the received signal is given by

$\displaystyle{\bf y}_{t}$ $\displaystyle=\sum_{k=1}^{K}\left({\bf h}_{{\rm
d},k}+{\bf H}_{{\rm I},k}{\bm{\theta}}_{t}\right)x_{k,t}+{\bf z}_{t}.$ (4)

In [17, 18, 19] and [20, 21, 22], ${\bf H}_{{\rm I},k}$ for all $k$ are
estimated without exploiting the implicit common link structure behind the
${\bf H}_{{\rm I},k}$ for all $k$. As a result, there are $MNK$ variables to
be estimated and this poses a heavy penalty on the required pilot overheads in
the MU-MISO system. On the other hand, from (3), we can see that the cascaded
channels $\left\\{{\bf H}_{{\rm I},1},{\bf H}_{{\rm I},2},\cdots,{\bf H}_{{\rm
I},K}\right\\}$ all share a common BS-IRS link $\bf G$. Specifically, ${\bf
H}_{{\rm I},k}$ is the multiplication of the common $\bf G$ and the user-
specific ${\rm{diag}}({\bf h}_{{\rm r},k})$. In other words, the cascaded
channels are not independent variables. In fact, the common link structure
$\bf G$ should be exploited in the channel estimation. As a result, the total
number of independent variables is reduced to $MN+NK$.

We assume ${\mathbb{E}}\left[{\bf h}_{{\rm I},k,m}\right]={\bf 0}$ for all $k$
and $m$.555 The proposed algorithm in this paper is applicable for the case
when an LoS link exists and ${\mathbb{E}}\left[{\bf h}_{{\rm
I},k,m}\right]\neq{\bf 0}$. Simply substitute ${\mathbb{E}}\left[{\bf h}_{{\rm
I},k,m}\right]$ into the prior distribution model in (22). Note that since the
rank-1 LoS link usually is very strong and easier to be estimated, it will
dominate the power of the channel coefficients, and the NMSE will be better
than the NLoS scenario investigated in this paper.  The covariance of the
cascaded channel ${\bf h}_{{\rm I},k,m}$ is given by

$\displaystyle{\bf C}_{m}^{(k)}={\mathbb{E}}\left[{\bf h}_{{\rm I},k,m}{\bf
h}_{{\rm I},k,m}^{\rm H}\right].$ (5)

In this paper, we focus on the case when ${\bf C}_{m}^{(k)}$ is a full rank
matrix for all $m$ and $k$, which is generally true in sub-6 GHz bands. We
design a channel estimation algorithm and phase shifting configuration scheme
by utilizing knowledge of ${\bf C}_{m}^{(k)}$. Note that in [21] and [22], the
channel estimation algorithms for the cascaded channels require knowledge of
the covariance of the IRS-BS link $\bf G$ as well as the covariance of the
IRS-user links ${\bf h}_{{\rm r},k}$. Also note that knowledge of the
covariance of the cascaded channel ${\bf C}_{m}^{(k)}$ is a weaker requirement
compared to knowledge of the individual covariances $\bf G$ and ${\bf h}_{{\rm
r},k}$.

##  III Proposed Channel Estimation Protocol

###  III-A Overview of the Selected On-Off Channel Estimation Protocol

The selected on-off channel estimation protocol in [23] is illustrated in
Figure 3, and consists of three stages. In stage I, the BS-user channels are
estimated by switching off all the IRS elements. In stage II, a reference user
is selected, which is indexed by user $1$, and its cascaded channel ${\bf
H}_{{\rm I},1}$ is estimated using the algorithm for SU-MISO cases [15]. In
stage III, the other $K-1$ users’ cascaded channels are estimated by
exploiting the common-link property.

Figure 3: Selected on-off channel estimation protocol in [23].

We focus on the estimation in stage III. Substituting ${\bf H}_{{\rm
I},1}={\bf G}^{\rm T}{\rm{diag}}({\bf h}_{{\rm r},1})$ into (1), the received
signal is given by

$\displaystyle{\bf y}_{t}$ $\displaystyle=\sum_{k=1}^{K}{\bf h}_{{\rm
d},k}x_{k,t}+{\bf H}_{{\rm I},1}{\bm{\theta}}_{t}x_{1,t}$ (6)
$\displaystyle\qquad+\sum_{k=2}^{K}{\bf H}_{{\rm I},1}{\rm
diag}({\bm{\theta}}_{t}){\bf h}_{{\rm u},k}x_{k,t}+{{\bf z}_{t}},$

where ${\bf h}_{{\rm u},k}={\rm diag}({\bf h}_{{\rm r},1})^{-1}{\bf h}_{{\rm
r},k}$ for all $k=2,3,\cdots,K$, which are the user-specific variables to be
estimated in this stage after exploiting ${\bf H}_{{\rm I},1}$ as the common-
link variable. In [23], to estimate ${\bf h}_{{\rm u},k}$, only the $k$-th
user sends $x_{k,t}=1$ and all the other users are inactive such that
$x_{j,t}=0$ for all $j\neq k$. The received signal is given by

$\displaystyle{\bf y}_{t}$ $\displaystyle={\bf h}_{{\rm d},k}+{\bf H}_{{{\rm
I},1}}{\rm diag}({\bm{\theta}}_{t}){\bf h}_{{\rm u},k}+{\bf z}_{t}.$ (7)

By only switching on the first $M$ IRS elements with
$[\theta_{t,1},\theta_{t,2},\cdots,\theta_{t,M}]^{\rm T}={\bf 1}$, the first
$M$ coefficients in ${\bf h}_{{\rm u},k}$ can be estimated by

$\displaystyle\begin{bmatrix}{h}_{{\rm u},k,1}\\\ \vdots\\\ {h}_{{\rm
u},k,M}\end{bmatrix}={\begin{bmatrix}H_{{{\rm I},1},11}&\dots&H_{{{\rm
I},1},1M}\\\ \vdots&\ddots&\vdots\\\ H_{{{\rm I},1},M1}&\dots&H_{{{\rm
I},1},MM}\end{bmatrix}}^{-1}({\bf y}_{i}-{\bf h}_{{\rm d},k}).$ (8)

Note that one may adopt the LMMSE estimator to achieve better performance if
the covariance of ${\bf h}_{{\rm u},k}={\rm diag}({\bf h}_{{\rm
r},1})^{-1}{\bf h}_{{\rm r},k}$ is available [23, Section V]. In the next
timeslot, the next $M$ IRS elements are switched on with
$[\theta_{t,M+1},\theta_{t,M+2},\cdots,\theta_{t,2M}]^{\rm T}={\bf 1}$ while
the other elements are switched off to estimate the next $M$ coefficients in
${\bf h}_{{\rm u},k}$. The estimation continues in this way until all the
coefficients in ${\bf h}_{{\rm u},k}$ are estimated, which finally costs
$I=\lceil\frac{N}{M}\rceil$ timeslots. The overall pilot overhead of the
protocol in [23] is $K+N+\lceil\frac{N}{M}\rceil(K-1)$.

###  III-B Always-ON Channel Estimation Protocol

We propose a novel always-ON channel estimation protocol without switching off
selected IRS elements. The proposed protocol consists of two stages, as
illustrated in Figure 4. In particular, stage I contains $L_{1}+1$ timeslots,
where $L_{1}=\lceil\frac{N}{M}\rceil$, and each timeslot contains $K$ samples.
Stage II contains $L_{2}=N-L_{1}$ timeslots, and each timeslot contains only
one sample. Therefore, we have $K(L_{1}+1)+L_{2}$ received samples in total.
As defined in (1), the $t$-th received sample is given by

$\displaystyle{\bf y}_{t}$ $\displaystyle=\sum_{k=1}^{K}\left({\bf h}_{{\rm
d},k}+{\bf G}^{\rm T}{\rm diag}({\bm{\theta}}_{t}){\bf h}_{{\rm
r},k}\right)x_{k,t}+{{\bf z}_{t}},$ (9)

where $t=1,2,\cdots,K(L_{1}+1)+L_{2}$.

Figure 4: Proposed always-ON protocol.

####  III-B1 Received Signal Samples in Each Training Timeslot

To facilitate analysis, we introduce new notations on the received signal
samples within one training timeslot indexed by $\ell=0,1,\cdots,N$.

  * •

Stage I (Users send orthogonal pilot sequences ${\bf X}$): Stage I consists of
training timeslots $\ell=0,1,\cdots,L_{1}$. Denote by ${\bf
X}\in{\mathbb{C}}^{K\times K}$, where ${\bf X}^{\rm H}{\bf X}=K{\bf I}_{K}$,
the orthogonal pilot sequences consisting of unit-modulus elements. In the
$\ell$-th timeslot, $K$ users transmit ${\bf X}$ by $K$ samples, and the IRS
is configured by the phase shifting vector ${\bm{\theta}}_{\ell}$. The $K$
received samples in the $\ell$-th timeslot are ${\bf y}_{K\ell+1},{\bf
y}_{K\ell+2},\cdots,{\bf y}_{K\ell+K}$. Let ${{\bf Y}}_{\ell}=[{\bf
y}_{K\ell+1},{\bf y}_{K\ell+2},\cdots,{\bf y}_{K\ell+K}]$. We have

$\displaystyle{{\bf Y}}_{\ell}$ $\displaystyle=\left({\bf H}_{\rm d}+{\bf
G}^{\rm T}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm r}\right){\bm{X}}+{{\bf
Z}}_{\ell},$ (10)

where

$\displaystyle{\bf H}_{\rm d}$ $\displaystyle=[{\bf h}_{{\rm d},1},\cdots,{\bf
h}_{{\rm d},K}]\in{\mathbb{C}}^{M\times K},$ (11) $\displaystyle{\bf H}_{\rm
r}$ $\displaystyle=[{\bf h}_{{\rm r},1},\cdots,{\bf h}_{{\rm
r},K}]\in{\mathbb{C}}^{N\times K}$ (12)

are the stacked channel coefficient matrices, $\ell$ is the timeslot index,
and ${{\bf Z}}_{\ell}=[{\bf z}_{K\ell+1},{\bf z}_{K\ell+2},\cdots,{\bf
z}_{K\ell+K}]$ denotes the noise.

  * •

Stage II (Users send pilot $\bar{\bf x}$, which is the first column of $\bf
X$): Stage II consists of training timeslots $\ell=L_{1}+1,L_{1}+2,\cdots,N$.
We denote the first column of $\bf X$ by

$\bar{\bf x}=[{\bar{x}}_{1},\cdots,{\bar{x}}_{K}]^{\rm T}.$ (13)

In the $\ell$-th timeslot, the users transmit $\bar{\bf x}$, while the IRS is
configured by ${\bm{\theta}}_{\ell}$. The received signal in stage II is
denoted by

$\displaystyle\bar{\bf y}_{\ell}$ $\displaystyle=\left({\bf H}_{\rm d}+{\bf
G}^{\rm T}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm r}\right)\bar{\bf
x}+{\bf z}_{\ell+(K-1)L_{1}+K},$ (14)

where $\bar{\bf y}_{\ell}={\bf y}_{\ell+(K-1)L_{1}+K}$, and
$\ell=L_{1}+1,L_{1}+2,\cdots,N$.

$\bar{\bf y}_{\ell}$ may extend to all the $N+1$ timeslots in the protocol by

$\displaystyle\bar{\bf y}_{\ell}=\begin{cases}{\bf y}_{1+K\ell},\;&{\rm
for}\;\ell=0,\cdots,L_{1},\\\ {\bf y}_{\ell+(K-1)L_{1}+K},&{\rm
for}\;\ell=L_{1}+1,\cdots,N.\end{cases}$ (15)

The pilot overhead in stage I is $(\lceil\frac{N}{M}\rceil+1)K$, and the
overhead in stage II is $N-\lceil\frac{N}{M}\rceil$. The overall pilot
overhead of the proposed protocol is $K+N+\lceil\frac{N}{M}\rceil(K-1)$. Note
that the pilot overhead is significantly reduced by about $M$ times compared
to the protocol in [17, 18, 19] and [20, 21, 22], which require $NK$ pilots.

####  III-B2 Signal Pre-processing

In the proposed protocol, we set ${\bm{\theta}}_{0}=-{\bm{\theta}}_{1}$ to
decouple the estimation on the BS-user channel666 The BS-user channel can be
estimated based on ${{\bf Y}}_{0}$ and ${{\bf Y}}_{1}$ by the linear minimum
mean squared error estimator, which is similar to the methods in existing
works [21, 22, 23] (See Appendix B). and on the cascaded IRS channel. To
facilitate the estimation of the cascaded IRS channel, signal pre-processing
to remove the BS-user channel from the received signals is performed.

  * •

_Pre-processing on ${{\bf Y}}_{\ell}$ for $\ell=1,2,\cdots,L_{1}$_: For
$\ell=1$, we have

$\displaystyle{\bf R}_{1}$ $\displaystyle=\frac{1}{2}\left({{\bf Y}}_{1}-{{\bf
Y}}_{0}\right)$ (16) $\displaystyle={\bf G}^{\rm T}{\rm
diag}({\bm{\theta}}_{1}){\bf H}_{\rm r}{\bf X}+\tilde{\bf Z}_{1},$

where the elements in $\tilde{\bf Z}_{1}$ follow i.i.d.
${\cal{CN}}({0},\frac{1}{2}\sigma_{0}^{2})$. For $\ell=2,3,\cdots,L_{1}$, we
have

$\displaystyle{\bf R}_{\ell}$ $\displaystyle={{\bf
Y}}_{\ell}-\frac{1}{2}\left({{\bf Y}}_{0}+{{\bf Y}}_{1}\right)$ (17)
$\displaystyle={\bf G}^{\rm T}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm
r}{\bf X}+\tilde{\bf Z}_{\ell},$

where the elements in $\tilde{\bf Z}_{\ell}$ follow
${\cal{CN}}({0},\frac{3}{2}\sigma_{0}^{2})$. Note that the difference between
(16) and (17) is that they have different noise variances.

  * •

_Pre-processing on $\bar{\bf y}_{\ell}$ for $\ell=1,2,\cdots,N$_: In the same
manner, the BS-user channel is removed from $\bar{\bf y}_{\ell}$ for
$\ell=1,2,\cdots,N$:

$\displaystyle\bar{\bf r}_{\ell}$ $\displaystyle=\bar{\bf
y}_{\ell}-\frac{1}{2}\left(\bar{\bf y}_{0}+\bar{\bf y}_{1}\right)$ (18)
$\displaystyle={\bf G}^{\rm T}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm
r}\bar{\bf x}+\bar{\bf z}_{\ell}$ $\displaystyle={\bf G}^{\rm T}{\rm
diag}({\bf H}_{\rm r}\bar{\bf x}){\bm{\theta}}_{\ell}+\bar{\bf z}_{\ell},$

where the elements in $\bar{\bf z}_{\ell}$ follow
${\cal{CN}}({0},\frac{3}{2}\sigma_{0}^{2})$.

######  Lemma 1 (Effectiveness of the proposed protocol).

The cascaded channels ${\bf H}_{{\rm I},k}={\bf G}^{\rm T}{\rm{diag}}({\bf
h}_{{\rm r},k})$ ($k=1,2,\cdots,K$) for all $K$ users can be perfectly
recovered with probability one by adopting an orthogonal phase shifting
configuration matrix
${\bm{\Phi}}=[{\bm{\theta}}_{1},{\bm{\theta}}_{2},\cdots,{\bm{\theta}}_{N}]$
whose elements are all non-zero, if there is no noise, and ${\bf G}={\bf
F}_{\rm R}{\ddot{\bf G}}{\bf F}_{\rm B}^{\rm T}$ and ${\bf h}_{{\rm r},k}={\bf
F}_{\rm R}{\ddot{\bf h}}_{{\rm r},k}$ where ${\bf F}_{\rm
B}\in{\mathbb{C}}^{M\times M}$ and ${\bf F}_{\rm R}\in{\mathbb{C}}^{N\times
N}$ is respectively the angular domain basis for the BS antenna array and the
IRS, ${\ddot{\bf G}}=[{\ddot{\bf g}}_{1},{\ddot{\bf g}}_{2},\cdots,{\ddot{\bf
g}}_{M}]$, and ${\ddot{\bf g}}_{m}$ and ${\ddot{\bf h}}_{{\rm r},k}$ for all
$m$ and $k$ are pairwise independent following zero-mean multivariate normal
distributions.

###### Proof:

See Appendix A. ∎

Compared with [23] and [28], the proposed protocol has two main differences,
which provide the opportunity to keep all the IRS elements ON and to reduce
the pilot overhead by utilizing the common-link structure at the same time.
Here, we try to explain the intuition by supposing we adopt a similar
estimation algorithm to those in [23] and [28], i.e., first estimate the
reference channel and then estimate the relative channels. Firstly, instead of
requiring a specific reference user, we design a virtual reference channel
${\bf H}_{\rm v}={\bf G}^{\rm T}{\rm diag}({\bf H}_{\rm r}\bar{\bf x})$, which
is fair for all users and can be reconstructed using the $N$ observations in
(18). Secondly, in our protocol, the relative channels can be estimated by
using ${\bf R}_{\ell}$ in (16) and (17) resulting a measurement matrix $[({\bf
H}_{\rm v}{\rm diag}({\bm{\theta}}_{1}))^{\rm T},\cdots,({\bf H}_{\rm v}{\rm
diag}({\bm{\theta}}_{L_{1}}))^{\rm T}]^{\rm T}$. One can see that there are
$L_{1}$ different ${\bm{\theta}}_{\ell}$ in the measurement matrix instead of
a fixed one as in [28]. When $L_{1}\geq\lceil\frac{N}{M}\rceil$, it is
possible that the rank of the measurement matrix becomes $N$ with a proper
phase shifting configuration to obtain a reasonable estimation.777 It is seen
that the selected on-off protocol in [23] can be treated as a special case of
our protocol by setting $\bar{\bf x}=[1,0,\cdots,0]^{\rm T}$ and selecting
different parts of the IRS elements to be ON and OFF to obtain $L_{1}$
different ${\bm{\theta}}_{\ell}$, as introduced in Section III-A.  Note that
the above two-step channel estimation algorithm is only for explaining the
intuition of the proposed approach. In the next section, we will propose an
optimization-based cascaded channel estimation algorithm that may achieve more
reliable performance.

##  IV Optimization-Based MU-Cascaded IRS Channel Estimation

In this section, we propose an optimization-based channel estimation on the
cascaded IRS channel ${\bf h}_{{\rm I},k,m}$ for all $m$ and $k$ based on the
pre-processed observations $\\{{\bf R}_{\ell},\bar{\bf r}_{\ell}\\}$ in (16),
(17) and (18). We consider a general decomposition on the cascaded channel,
which is friendly in utilizing the channel prior knowledge and the common-link
structure across the multiple users. Specifically, we adopt the MAP approach
to estimate the cascaded channel ${\bf h}_{{\rm I},k,m}$ for all $m$ and $k$
given the pre-processed observations $\\{{\bf R}_{\ell},\bar{\bf
r}_{\ell}\\}$. An alternating optimization algorithm with efficient
initialization is further proposed to achieve a local optimum of the MAP
problem.

###  IV-A MAP Problem Formulation

As shown in (2), the cascaded channel ${\bf h}_{{\rm I},k,m}$ can be
decomposed by the common BS-IRS channel ${\bf g}_{m}$ and the IRS-user channel
${\bf h}_{{\rm r},k}$ as follows:

$\displaystyle{\bf h}_{{\rm I},k,m}={\rm{diag}}({\bf h}_{{\rm r},k}){\bf
g}_{m}.$ (19)

However, the main challenge to estimate the individual ${\bf h}_{{\rm r},k}$
and ${\bf g}_{m}$ is that in the MAP formulation, prior distributions of ${\bf
g}_{m}$ and ${\bf h}_{{\rm r},k}$ will be needed, but it is difficult to
obtain individual covariances of ${\bf g}_{m}$ and ${\bf h}_{{\rm r},k}$ based
on the covariance of the cascaded channel.888 One possible way to estimate the
covariance of the cascaded channel ${\bf C}_{m}^{(k)}$ is using the the
maximum likelihood estimator $\hat{\bf
C}_{m}^{(k)}=\frac{1}{J}\sum_{j=1}^{J}\left[\hat{\bf h}_{{\rm
I},k,m}(j)\hat{\bf h}_{{\rm I},k,m}(j)^{\rm H}\right]$, where $\\{\hat{\bf
h}_{{\rm I},k,m}(j)\\}$ are the estimated historical cascaded channels in the
past $j=1,2,\cdots,J$ transmission frames. Note that similar covariances are
also required by the LMMSE estimators for the selected on-off channel
estimation protocol in [23] (See equations (72) and (86) in [23]).

To address this issue, we consider a more general auxiliary variable set for
the cascaded channel decomposition:

$\displaystyle{\cal A}=\left\\{\left.\left\\{{\bf H}_{\rm g},{\bf H}_{\rm
u}\right\\}\right|{\bf h}_{{\rm I},k,m}={\rm{diag}}({\bf h}_{{\rm u},k}){\bf
h}_{{\rm g},m},\forall k,\forall m\right\\},$ (20)

where ${\bf H}_{\rm g}=[{\bf h}_{{\rm g},1},\cdots,{\bf h}_{{\rm g},M}]$ is
the common-link variable whose $m$-th column is ${\bf h}_{{\rm g},m}$, and
${\bf H}_{\rm u}=[{\bf h}_{{\rm u},1},\cdots,{\bf h}_{{\rm u},K}]$ is the
user-specific variable whose $k$-th column is ${\bf h}_{{\rm u},k}$. One may
verify that $\\{{\bf G},{\bf H}_{\rm r}\\}\in{\cal A}$ according to (19).
Based on this, we formulate an optimization problem on ${\bf H}_{\rm g}$ and
${\bf H}_{\rm u}$ using the MAP approach, which is given by

$\displaystyle{\mathcal{P}}{(\text{A})}\quad\max_{{\bf H}_{\rm g},{\bf H}_{\rm
u}}\;f_{{\rm A}}({\bf H}_{\rm g},{\bf H}_{\rm u}),$

where the objective function is given by999 One may substitute different prior
setups to $p({\bf h}_{{\rm I},k,m})$ in (22). In addition, if channel prior
knowledge is unavailable, one may simply remove $p({\bf h}_{{\rm I},k,m})$
from (21), and the optimization becomes the maximum likelihood approach.

$\displaystyle f_{{\rm A}}({\bf H}_{\rm g},{\bf H}_{\rm
u})=\sum_{m=1}^{M}\sum_{k=1}^{K}\ln p({\bf h}_{{\rm I},k,m})$ (21)
$\displaystyle\quad+\sum_{\ell=1}^{L_{1}}\ln p({\tilde{\bf R}}_{\ell}|{\bf
H}_{\rm g},{\bf H}_{\rm u})+\sum_{\ell=L_{1}+1}^{N}\ln p(\bar{\bf
r}_{\ell}|{\bf H}_{\rm g},{\bf H}_{\rm u}),$

and $\tilde{\bf R}_{\ell}={{\bf R}}_{\ell}{\bf X}^{-1}$. Note that $f_{{\rm
A}}({\bf H}_{\rm g},{\bf H}_{\rm u})$ in (21) only requires the prior
distribution of ${\bf h}_{{\rm I},k,m}$, which is given by

$\displaystyle p({\bf h}_{{\rm I},k,m})\propto e^{-{\bf h}_{{\rm g},m}^{\rm
H}{\rm{diag}}({\bf h}_{{\rm u},k})^{\rm H}{{\bf
C}_{m}^{(k)}}^{-1}{\rm{diag}}({\bf h}_{{\rm u},k}){\bf h}_{{\rm g},m}},$ (22)

where $\propto$ denotes equality up to a scaling that is independent of the
variables (i.e., ${\bf h}_{{\rm I},k,m}$ for (22)). The likelihood functions
$p({\tilde{\bf R}}_{\ell}|{\bf H}_{\rm g},{\bf H}_{\rm u})$ and $p(\bar{\bf
r}_{\ell}|{\bf H}_{\rm g},{\bf H}_{\rm u})$ are given by

$\displaystyle p(\tilde{{\bf R}}_{\ell}|{\bf H}_{\rm g},{\bf H}_{\rm u})$
$\displaystyle\propto e^{-\sigma_{\ell}^{-2}\left\|\tilde{{\bf R}}_{\ell}-{\bf
H}_{\rm g}^{\rm T}{\rm diag}\left({\bm{\theta}}_{\ell}\right){\bf H}_{{\rm
u}}\right\|^{2}_{\rm F}},$ (23) $\displaystyle p(\bar{\bf r}_{\ell}|{\bf
H}_{\rm g},{\bf H}_{\rm u})$ $\displaystyle\propto
e^{-\bar{\sigma}_{\ell}^{-2}\left\|\bar{\bf r}_{\ell}-{\bf H}_{\rm g}^{\rm
T}{\rm diag}\left({\bm{\theta}}_{\ell}\right){\bf H}_{{\rm u}}\bar{\bf
x}\right\|^{2}_{2}},$ (24)

where $\sigma_{1}^{2}=\frac{1}{2K}\sigma_{0}^{2}$,
$\sigma_{\ell}^{2}=\frac{3}{2K}\sigma_{0}^{2}$ for $\ell=2,3,\cdots,L_{1}$,
and $\bar{\sigma}_{\ell}^{2}=\frac{3}{2}\sigma_{0}^{2}$ for
$\ell=L_{1}+1,L_{1}+2,\cdots,N$. Finally, after dropping all the irrelevant
constant terms, the objective function is equivalently written as

$\displaystyle f_{{\rm A}}({\bf H}_{\rm g},{\bf H}_{\rm
u})=-\sum_{\ell=1}^{L_{1}}{\frac{1}{\sigma_{\ell}^{2}}\left\|\tilde{{\bf
R}}_{\ell}-{\bf H}_{\rm g}^{\rm T}{\rm
diag}\left({\bm{\theta}}_{\ell}\right){\bf H}_{{\rm u}}\right\|^{2}_{\rm F}}$
(25)
$\displaystyle\quad-\sum_{\ell=L_{1}+1}^{N}{\frac{1}{\bar{\sigma}_{\ell}^{2}}\left\|\bar{\bf
r}_{\ell}-{\bf H}_{\rm g}^{\rm T}{\rm
diag}\left({\bm{\theta}}_{\ell}\right){\bf H}_{{\rm u}}\bar{\bf
x}\right\|^{2}_{2}}$ $\displaystyle\quad-\sum_{m=1}^{M}\sum_{k=1}^{K}{\bf
h}_{{\rm g},m}^{\rm H}{\rm{diag}}({\bf h}_{{\rm u},k})^{\rm H}{{\bf
C}_{m}^{(k)}}^{-1}{\rm{diag}}({\bf h}_{{\rm u},k}){\bf h}_{{\rm g},m}.$

Note that ${\mathcal{P}}{(\text{A})}$ does not have a unique solution but all
the solutions are equivalent for the purpose of estimation of the cascaded
channel ${\bf h}_{{\rm I},k,m}$ for all $m$ and $k$.

######  Lemma 2 (Equivalence of the solution of ${\mathcal{P}}{(\text{A})}$).

Let ${\bf H}_{\rm g}^{\star}$ and ${\bf H}_{\rm u}^{\star}$ be an optimal
solution of ${\mathcal{P}}{(\text{A})}$, then ${\rm diag}\left({\bf
a}\right){\bf H}_{\rm g}^{\star}$ and ${\rm diag}\left({\bf a}\right)^{-1}{\bf
H}_{\rm u}^{\star}$ is also an optimal solution of ${\mathcal{P}}{(\text{A})}$
for any coefficient ${\bf a}=[a_{1},a_{2},\cdots,a_{N}]^{\rm T}$ with
$|a_{n}|\neq 0$ for all $n=1,2,\cdots,N$.

###### Proof:

See Appendix C. ∎

As a result, there is ambiguity in estimating individual channels from solving
${\mathcal{P}}{(\text{A})}$. Nevertheless, the cascaded channel is unique
regardless of the coefficient ${\bf a}$.

###  IV-B Channel Estimation Algorithm based on Alternative Optimization

Solving ${\mathcal{P}}{(\text{A})}$ is difficult due to the optimization
variables being coupled in the likelihood functions (23) and (24).
Fortunately, we will show that ${\mathcal{P}}{(\text{A})}$ is actually bi-
convex (see Lemma 3 and Lemma 4), which can be solved by alternative
optimization. In particular, we decompose ${\mathcal{P}}{(\text{A})}$ into two
convex sub-problems, and the optimal solutions for these sub-problems will be
derived accordingly.

####  IV-B1 Optimize ${\bf H}_{\rm u}$

We investigate the optimization of ${\bf H}_{\rm u}$ while ${\bf H}_{\rm g}$
are fixed. After dropping all irrelevant terms, the sub-problem is given by

$\displaystyle{\mathcal{P}}{({\text{A}}_{{\rm u}})}\quad\min_{{\bf H}_{\rm
u}}\;f_{{\rm u}}({\bf H}_{\rm u}),$

where

$\displaystyle f_{{\rm u}}({\bf H}_{\rm
u})=\sum_{\ell=1}^{L_{1}}{\frac{1}{\sigma_{\ell}^{2}}\left\|\tilde{{\bf
R}}_{\ell}-{\bf D}_{\ell}{\bf H}_{{\rm u}}\right\|^{2}_{\rm F}}$ (26)
$\displaystyle\quad+\sum_{\ell=L_{1}+1}^{N}{\frac{1}{\bar{\sigma}_{\ell}^{2}}\left\|\bar{\bf
r}_{\ell}-{\bf D}_{\ell}{\bf H}_{{\rm u}}\bar{\bf
x}\right\|^{2}_{2}}+\sum_{k=1}^{K}{\bf h}_{{\rm u},k}^{\rm H}{\bf C}_{{\rm
u},k}{\bf h}_{{\rm u},k},$

and

$\displaystyle{\bf D}_{\ell}$ $\displaystyle={\bf H}_{\rm g}^{\rm T}{\rm
diag}\left({\bm{\theta}}_{\ell}\right),\;{\text{for}}\;\ell=1,2,\cdots,N,$
(27) $\displaystyle{\bf C}_{{\rm u},k}$
$\displaystyle=\sum_{m=1}^{M}{\rm{diag}}({\bf h}_{{\rm g},m})^{\rm H}{{\bf
C}_{m}^{(k)}}^{-1}{\rm{diag}}({\bf h}_{{\rm g},m}),\quad\forall k.$ (28)

######  Lemma 3 (Convexity of ${\mathcal{P}}{({\text{A}}_{{\rm u}})}$).

For any fixed ${\bf H}_{\rm g}$, the objective function of
${\mathcal{P}}{({\text{A}}_{{\rm u}})}$ is a convex quadratic function of the
vectorization of ${\bf H}_{\rm u}$, which is denoted by ${\rm{vec}}({\bf
H}_{\rm u})$.

###### Proof:

See Appendix D. ∎

Based on Lemma 3, the optimal solution of ${\mathcal{P}}{({\text{A}}_{{\rm
u}})}$ is the root of the first order derivative of $f_{{\rm
u}}({\rm{vec}}({\bf H}_{\rm u}))$, which is given by

$\displaystyle{\rm{vec}}({\bf H}_{\rm u}^{\star})$
$\displaystyle={\bm{\Lambda}}_{{\rm u}}^{-1}{\bm{\nu}}_{{\rm u}},$ (29)

where

$\displaystyle{\bm{\Lambda}}_{{\rm u}}$ $\displaystyle={\bf C}_{{\rm u}}+{\bf
I}_{K}\otimes\left(\sum_{\ell=1}^{L_{1}}\frac{1}{\sigma_{\ell}^{2}}{\bf
D}_{\ell}^{\rm H}{\bf D}_{\ell}\right)$ $\displaystyle\quad+\left(\bar{\bf
x}^{\ast}\bar{\bf x}^{\rm
T}\right)\otimes\left(\sum_{\ell=L_{1}+1}^{N}\frac{1}{\bar{\sigma}_{\ell}^{2}}|\bar{x}_{k}|^{2}{\bf
D}_{\ell}^{\rm H}{\bf D}_{\ell}\right),$ (30) $\displaystyle{\bf C}_{{\rm u}}$
$\displaystyle={\rm{blkdiag}}({\bf C}_{{\rm u},1},{\bf C}_{{\rm
u},2},\cdots,{\bf C}_{{\rm u},K}),$ (31)

and

$\displaystyle{\bm{\nu}}_{{\rm u}}$
$\displaystyle={\rm{vec}}\left(\sum_{\ell=1}^{L_{1}}\frac{1}{\sigma_{\ell}^{2}}{\bf
D}_{\ell}^{\rm H}\tilde{{\bf
R}}_{\ell}+\sum_{\ell=L_{1}+1}^{N}\frac{1}{\bar{\sigma}_{\ell}^{2}}{\bf
D}_{\ell}^{\rm H}\bar{\bf r}_{\ell}\bar{\bf x}^{\rm H}\right).$ (32)

####  IV-B2 Optimize ${\bf H}_{{\rm g}}$

Similarly, the sub-problem of optimizing ${\bf H}_{{\rm g}}$ is given by

$\displaystyle{\mathcal{P}}{(\text{A}_{{\rm g}})}\;\min_{{\bf H}_{{\rm
g}}}\;\sum_{m=1}^{M}f_{{\rm g},m}({\bf h}_{{\rm g},m}),$

where

$\displaystyle f_{{\rm g},m}({\bf h}_{{\rm
g},m})=\sum_{k=1}^{K}\sum_{\ell=1}^{L_{1}}\frac{1}{\sigma_{\ell}^{2}}\left\|\tilde{{r}}_{\ell,m,k}-{\bf
h}_{{\rm g},m}^{\rm T}{\bf b}_{\ell,k}\right\|^{2}$
$\displaystyle+\sum_{\ell=L_{1}+1}^{N}\frac{1}{\bar{\sigma}_{\ell}^{2}}\left\|\bar{r}_{\ell,m}-\sum_{k=1}^{K}\bar{x}_{k}{\bf
h}_{{\rm g},m}^{\rm T}{\bf b}_{\ell,k}\right\|^{2}+{\bf h}_{{\rm g},m}^{\rm
H}{\bf C}_{{\rm g},m}{\bf h}_{{\rm g},m},$ (33) $\displaystyle\quad{\bf
b}_{\ell,k}={\rm diag}\left({\bm{\theta}}_{\ell}\right){\bf h}_{{\rm
u},k},\;{\text{for}}\;\ell=1,2,\cdots,N,$ (34) $\displaystyle\quad{\bf
C}_{{\rm g},m}=\sum_{k=1}^{K}{\rm{diag}}({\bf h}_{{\rm u},k})^{\rm H}{{\bf
C}_{m}^{(k)}}^{-1}{\rm{diag}}({\bf h}_{{\rm u},k}),\quad\forall m,$ (35)

$\tilde{{r}}_{\ell,m,k}$ denotes the entry in the $m$-th row and $k$-th column
of $\tilde{\bf R}_{\ell}$, and $\bar{r}_{\ell,m}$ denotes the $m$-th entry in
$\bar{\bf r}_{\ell}$.

######  Lemma 4 (Convexity of ${\mathcal{P}}{({\text{A}}_{{\rm g}})}$).

For any fixed ${\bf H}_{\rm u}$, the objective function of
${\mathcal{P}}{({\text{A}}_{{\rm g}})}$ is a convex quadratic function of
$\\{{\bf h}_{{\rm g},1},{\bf h}_{{\rm g},2},\cdots,{\bf h}_{{\rm g},M}\\}$.

###### Proof:

See Appendix E. ∎

Based on lemma 4, the optimal ${\bf h}_{{\rm g},m}$ is the root of the first
order derivative of $f_{{\rm g},m}({\bf h}_{{\rm g},m})$, which is given by

$\displaystyle{\bf h}_{{\rm g},m}^{\star}$ $\displaystyle={\bm{\Lambda}}_{{\rm
g},m}^{-1}{\bm{\nu}}_{{\rm g},m},$ (36)

for $m=1,2,\cdots,M$, where

$\displaystyle{\bm{\Lambda}}_{{\rm g},m}$ $\displaystyle={\bf C}_{{\rm
g},m}^{-1}+\sum_{k=1}^{K}\sum_{\ell=1}^{L_{1}}\frac{1}{\sigma_{\ell}^{2}}{\bf
b}_{\ell,k}^{\ast}{\bf b}_{\ell,k}^{\rm T}$ (37)
$\displaystyle\;+\sum_{\ell=L_{1}+1}^{N}\frac{1}{\bar{\sigma}_{\ell}^{2}}\left(\sum_{k=1}^{K}\bar{x}_{k}{\bf
b}_{\ell,k}^{\rm T}\right)^{\rm H}\left(\sum_{k=1}^{K}\bar{x}_{k}{\bf
b}_{\ell,k}^{\rm T}\right),$

and

$\displaystyle{\bm{\nu}}_{{\rm
g},m}=\sum_{k=1}^{K}\sum_{{\ell}=1}^{L_{1}}\frac{\tilde{{r}}_{\ell,m,k}}{\sigma_{\ell}^{2}}{\bf
b}_{\ell,k}^{\ast}+\sum_{\ell=L_{1}+1}^{N}\frac{\bar{r}_{\ell,m}}{\bar{\sigma}_{\ell}^{2}}\left(\sum_{k=1}^{K}\bar{x}_{k}^{\ast}{\bf
b}_{\ell,k}^{\ast}\right).$ (38)

####  IV-B3 Initial Estimation on ${\bf H}_{{\rm g}}$

The quality of the solution obtained by the alternative optimization depends
heavily on the initial point. Here, we propose an efficient estimator for
${\bf H}_{{\rm g}}$ to initialize the proposed alternative optimization
algorithm. In particular, we construct a special $\left\\{{\bf H}_{\rm g},{\bf
H}_{\rm u}\right\\}$ pair whose elements are given by

$\displaystyle{\bf h}_{{\rm g},m}$ $\displaystyle={\rm diag}({\bf H}_{\rm
r}\bar{\bf x}){\bf g}_{m},$ (39) $\displaystyle{\bf h}_{{\rm u},k}$
$\displaystyle={\rm diag}({\bf H}_{\rm r}\bar{\bf x})^{-1}{\bf h}_{{\rm
r},k}.$ (40)

Substituting the above ${\bf H}_{\rm g}$ into (18), we have

$\displaystyle\bar{\bf r}_{\ell}$ $\displaystyle={\bf H}_{\rm g}^{\rm
T}{\bm{\theta}}_{\ell}+\bar{\bf z}_{\ell},$ (41)

for $\ell=1,2,\cdots,N$. Therefore, ${\bf H}_{\rm g}$ can be initialized by
the least squares (LS) estimator as follows:

$\displaystyle\hat{\bf H}_{\rm g}=\left(\left[\bar{\bf r}_{1},\cdots,\bar{\bf
r}_{N}\right]\left[{\bm{\theta}}_{1},\cdots,{\bm{\theta}}_{N}\right]^{-1}\right)^{\rm
T}.$ (42)

Since $\hat{\bf H}_{\rm g}$ in (42) is unbiased and it has exploited most of
the available observations in all the $N+1$ training timeslots, it will give a
good initial point.

###  IV-C The Overall Proposed Algorithm

The overall proposed cascaded IRS channel estimation algorithm is summarized
in Algorithm 1. The convergence of the proposed alternating optimization
algorithm is analyzed in Lemma 5.

Algorithm 1 Proposed alternating optimization algorithm for the cascaded IRS
channel estimation

1: Initialize ${\bf H}_{{\rm g}}$ by (42).

2: repeat

3: Update ${\bf H}_{{\rm u}}$ by (29);

4: Update ${\bf h}_{{\rm g},m}$ by (36) for all $m$;

5: until $f_{{\rm A}}({\bf H}_{\rm g},{\bf H}_{\rm u})$ in (25) converges;

6: Output the cascaded channel ${\hat{\bf h}}_{{\rm I},k,m}={\rm{diag}}({\bf
h}_{{\rm u},k}){\bf h}_{{\rm g},m}$ for all $k$ and $m$.

######  Lemma 5 (Convergence of the Proposed Alternating Optimization
Algorithm).

The objective function $f_{{\rm A}}({\bf H}_{\rm g},{\bf H}_{\rm u})$ is non-
increasing in every step when ${\bf H}_{{\rm u}}$ or ${\bf H}_{{\rm g}}$ are
updated, and the optimization iterations in (29) and (36) converge to a local
optimum of ${\mathcal{P}}{(\text{A})}$.

###### Proof:

As shown in Section IV-B, the original problem is decomposed into two
unconstrained minimization problems whose objectives are convex quadratic
functions, and each subproblem has a unique optimal solution, which is derived
in (29) and (36). Therefore, the whole alternating optimization algorithm will
converge to a local optimum of the original problem
${\mathcal{P}}{(\text{A})}$ [37]. ∎

_Remark_ : The complexity for updating ${\bf H}_{{\rm u}}$ by (29) is ${\cal
O}(K^{3}N^{3}+KMN^{2})$. The complexity for updating ${\bf h}_{{\rm g},m}$ by
(36) is ${\cal O}(KN^{3})$, and thus the complexity to update ${\bf H}_{{\rm
g}}$ is ${\cal O}(KMN^{3})$. Therefore, the overall complexity of the solution
is ${\cal O}(IK^{3}N^{3}+IKMN^{3})$, where $I$ denotes the number of
iterations of the alternating optimization algorithm.101010We will show in
simulation that the algorithm will converge quickly in about two or three
iterations. In addition, the ${\cal O}(N^{3})$ complexity is costed by the
matrix inversion operation. However, since the matrices required inversion
operation are all Hermitian positive semi-definite matrices, they can be
implemented very efficiently by advanced algorithms such as the Cholesky-
decomposition-based algorithm [38].

##  V Training Phase Shifting Configuration

###  V-A Motivation of the Phase Shifting Configuration

The IRS steers the incident signal to different directions by configuring
different phase shifting vectors ${\bm{\theta}}_{\ell}$, as illustrated in
Figure 5. According to the protocol, the cascaded channel is scanned by $N$
spatial directions in $N$ training timeslots, and a proper design on
$\\{{\bm{\theta}}_{1},{\bm{\theta}}_{2},\cdots,{\bm{\theta}}_{N}\\}$
guarantees that the whole channel information in all directions is contained
by the received signals such that good channel estimation performance can be
achieved.

Figure 5: Illustration of the impact of different phase shiftings.

In the SU-MISO scenario, the overall received measurements after removing the
pilots and the BS-user channels is given by

$\displaystyle{\bf Y}_{\rm{SU}}$ $\displaystyle={\bf H}_{{\rm
I},1}{\bm{\Phi}}+{\bf Z},$ (43)

where
${\bm{\Phi}}=[{\bm{\theta}}_{1},{\bm{\theta}}_{2},\cdots,{\bm{\theta}}_{N}]$.
The LS estimator may adopted as follows [20]:

$\displaystyle\hat{\bf H}_{{\rm I},1}$ $\displaystyle={\bf
Y}_{\rm{SU}}{\bm{\Phi}}^{-1},$ (44)

where $\hat{\bf H}_{{\rm I},1}$ denotes the estimated cascaded channel. Then
${\bm{\Phi}}$ is optimized by minimizing the MSE:

$\displaystyle\min_{{\bm{\Phi}}}\;$
$\displaystyle{\rm{tr}}\left(\left({\bm{\Phi}}{\bm{\Phi}}^{\rm
H}\right)^{-1}\right)$ $\displaystyle{\bf s.t.}\;$
$\displaystyle|{\bm{\Phi}}_{i,j}|=1,\quad\forall i,j=1,2,\cdots,N.$

It is proved in [20] that the optimal value of the MSE is $1$, which can be
achieved by the DFT matrix such that ${\bm{\Phi}}={\bf F}$, where

$\displaystyle{\bf F}=\begin{bmatrix}1&1&1&\cdots&1\\\ 1&e^{-\jmath
2\pi\frac{1}{N}}&e^{-\jmath 2\pi\frac{2}{N}}&\cdots&e^{-\jmath
2\pi\frac{N-1}{N}}\\\ 1&e^{-\jmath 2\pi\frac{2}{N}}&e^{-\jmath
2\pi\frac{4}{N}}&\cdots&e^{-\jmath 2\pi\frac{2(N-1)}{N}}\\\
\vdots&\vdots&\vdots&\ddots&\vdots\\\ 1&e^{-\jmath
2\pi\frac{N-1}{N}}&e^{-\jmath 2\pi\frac{2(N-1)}{N}}&\cdots&e^{-\jmath
2\pi\frac{(N-1)^{2}}{N}}\end{bmatrix}.$ (45)

For the protocol extended from the SU case [20, 21, 22] shown in Figure 2, the
transmit signals from users are the same in different timeslots. Therefore,
the columns of ${\bf F}$ may be permuted to any orders, and the MSE will
remain the same. However, in our proposed protocol in Figure 4, the transmit
signals are different in stage I and stage II. In particular, the received
signals in stage I contribute to the estimation on both the common-link
variable and the user-specific variables in the cascaded channels, while the
signals in stage II contribute to the common-link variable only. Hence, the
phase shifting vectors
${\bm{\theta}}_{1},{\bm{\theta}}_{2},\cdots,{\bm{\theta}}_{L_{1}}$ in stage I
require additional design.

###  V-B Optimization Formulation of the Phase Shifting Configuration for MU-
MISO IRS Systems

As shown in (42), the initial estimation on the common-link variable $\hat{\bf
H}_{\rm g}$ in (42) is almost the same as the estimator for the SU case shown
in (44). Therefore, we still adopt the DFT-based phase shifting configuration
for all the $N$ timeslots. Additionally, an additional steering direction
${\bm{\vartheta}}\in{\mathbb{C}}^{N\times 1}$ is introduced for a more
flexible design:

$\displaystyle{\bm{\Phi}}$ $\displaystyle={\rm{diag}}({\bm{\vartheta}}){\bf
F}.$ (46)

Denote by $\vartheta_{n}$ the $n$-th element in ${\bm{\vartheta}}$. We have
$|\vartheta_{n}|=1$ for all $n=1,2,\cdots,N$. One can see that the value of
${\rm{tr}}\left(\left({\bm{\Phi}}{\bm{\Phi}}^{\rm H}\right)^{-1}\right)$ is
kept at $1$ for any ${\bm{\vartheta}}$. Denote by ${\bf f}_{\ell}$ the
$\ell$-th column of $\bf F$. The training phase shifting vector in the
$\ell$-th timeslot is given by

$\displaystyle{\bm{\theta}}_{\ell}$
$\displaystyle={\rm{diag}}({\bm{\vartheta}}){\bf f}_{\ell},$ (47)

where $\ell=1,2,\cdots,N$. The remaining task is to design ${\bm{\vartheta}}$.

However, it is difficult to design a straightforward objective function to
optimize ${\bm{\vartheta}}$ since the MSE of the estimated ${\bf H}_{\rm u}$
by the proposed algorithm is complicated. Considering that we have the
knowledge of the covariances ${\bf C}_{m}^{(k)}$ of the cascaded channels for
all $m$ and $k$, the average received power of the effective IRS channel from
user $k$ to the $m$-th BS antenna in timeslot $\ell$ can be denoted by a
function of ${\bm{\vartheta}}$, ${\bf f}_{\ell}$ and ${\bf C}_{m}^{(k)}$:

$\displaystyle Q_{\ell,k,m}$ $\displaystyle={\mathbb{E}}\left[|{\bf
g}_{m}^{\rm T}{\rm{diag}}({\bm{\theta}}_{\ell}){\bf h}_{{\rm
r},k}|^{2}\right]$ (48) $\displaystyle={\mathbb{E}}\left[|{\bf h}_{{\rm
I},k,m}^{\rm T}{\bm{\theta}}_{\ell}|^{2}\right]$
$\displaystyle={\mathbb{E}}\left[{\bm{\theta}}_{\ell}^{\rm H}\left({\bf
h}_{{\rm I},k,m}{\bf h}_{{\rm I},k,m}^{\rm
H}\right)^{\ast}{\bm{\theta}}_{\ell}\right]$
$\displaystyle={\bm{\theta}}_{\ell}^{\rm H}\left({\bf
C}_{m}^{(k)}\right)^{\ast}{\bm{\theta}}_{\ell}$
$\displaystyle={\bm{\vartheta}}^{\rm H}{\rm{diag}}({\bf f}_{\ell})^{\rm
H}\left({\bf C}_{m}^{(k)}\right)^{\ast}{\rm{diag}}({\bf
f}_{\ell}){\bm{\vartheta}}.$

Since ${\bf f}_{\ell}$ and ${\bf C}_{m}^{(k)}$ are known variables, the
summation of $Q_{\ell,k,m}$ over antennas $m=1,2,\cdots,M$, users
$k=1,2,\cdots,K$ and timeslots $\ell=1,2,\cdots,L_{1}$ is a function of
${\bm{\vartheta}}$, which is given by

$\displaystyle f_{\rm B}({\bm{\vartheta}})$
$\displaystyle=\sum_{\ell=1}^{L_{1}}\sum_{k=1}^{K}\sum_{m=1}^{M}Q_{\ell,k,m}$
(49)
$\displaystyle=\sum_{\ell=1}^{L_{1}}\sum_{k=1}^{K}\sum_{m=1}^{M}{\bm{\vartheta}}^{\rm
H}{\rm{diag}}({\bf f}_{\ell})^{\rm H}\left({\bf
C}_{m}^{(k)}\right)^{\ast}{\rm{diag}}({\bf f}_{\ell}){\bm{\vartheta}}$
$\displaystyle={\bm{\vartheta}}^{\rm
H}\left(\sum_{\ell=1}^{L_{1}}\sum_{k=1}^{K}\sum_{m=1}^{M}{\rm{diag}}({\bf
f}_{\ell})^{\rm H}\left({\bf C}_{m}^{(k)}\right)^{\ast}{\rm{diag}}({\bf
f}_{\ell})\right){\bm{\vartheta}}.$

We define

$\displaystyle{\bf
E}=\left(\sum_{\ell=1}^{L_{1}}\sum_{k=1}^{K}\sum_{m=1}^{M}{\rm{diag}}({\bf
f}_{\ell})^{\rm H}\left({\bf C}_{m}^{(k)}\right)^{\ast}{\rm{diag}}({\bf
f}_{\ell})\right).$ (50)

The optimization problem on ${\bm{\vartheta}}$ is formulated to maximize
$f_{\rm B}({\bm{\vartheta}})$:

$\displaystyle{\mathcal{P}}{(\text{B})}\quad\max_{{\bm{\vartheta}}}$
$\displaystyle f_{\rm B}({\bm{\vartheta}})={\bm{\vartheta}}^{\rm H}{\bf
E}{\bm{\vartheta}}$ $\displaystyle{\bf s.t.}$
$\displaystyle|\vartheta_{n}|=1,\quad\forall n=1,2,\cdots,N.$

###  V-C Solution for ${\mathcal{P}}{(\text{B})}$

${\mathcal{P}}{(\text{B})}$ is a non-convex problem due to the maximizing of a
convex objective function and the unit-modulus constraints. We solve
${\mathcal{P}}{(\text{B})}$ by the successive convex approximation (SCA)
algorithm. In particular, a surrogate problem, shown as follows, is
iteratively solved:

$\displaystyle{\mathcal{P}}{({\text{B}}_{i})}\quad\max_{{\bm{\vartheta}}}$
$\displaystyle{f}_{\rm B}^{(i)}({\bm{\vartheta}},\bar{\bm{\vartheta}})$
$\displaystyle{\bf s.t.}$ $\displaystyle|\vartheta_{n}|=1,\quad\forall
n=1,2,\cdots,N,$

where $i$ is the iteration index, $\bar{\bm{\vartheta}}$ is the solution of
the surrogate problem in the $(i-1)$-th iteration, and ${f}_{\rm
B}^{(i)}({\bm{\vartheta}},\bar{\bm{\vartheta}})$ is the first-order
approximation of ${f}_{\rm B}({\bm{\vartheta}})$ at $\bar{\bm{\vartheta}}$:

$\displaystyle{f}_{\rm B}^{(i)}({\bm{\vartheta}},\bar{\bm{\vartheta}})$
$\displaystyle=2{\rm Re}\left\\{{\bar{\bm{\vartheta}}}^{\rm H}{\bf
E}{\bm{\vartheta}}\right\\}-{\bar{\bm{\vartheta}}}^{\rm H}{\bf
E}{\bar{\bm{\vartheta}}}.$ (51)

One can see that ${f}_{\rm B}^{(i)}({\bm{\vartheta}},\bar{\bm{\vartheta}})$ is
a linear function of ${\bm{\vartheta}}$, and thus the optimal solution of
${\mathcal{P}}{({\text{B}}_{i})}$ is given by

$\displaystyle{\bm{\vartheta}}=e^{\jmath\angle({\bf
E}{\bar{\bm{\vartheta}}})}.$ (52)

The proof on the convergence of the SCA algorithm can be referred to in [39].

##  VI Numerical Examples

###  VI-A Simulation Setups

This section evaluates the performance of the proposed cascaded channel
estimation algorithm. In particular, we consider the indoor femtocell network
illustrated in Figure 6 in which $K=8$ users are randomly distributed in a 5 m
$\times$ 5 m square area and are served by one BS and one IRS. We generate the
channel coefficients according to the 3GPP ray-tracing model [40, Section 7.5]
using the model parameters for the Indoor-Office scenario [40, Table 7.5-6].
The system parameters for the simulations are summarized in Table I, in which
the path-loss is set according to the Indoor-Office pathloss model in [40,
Table 7.4.1-1].

Figure 6: The simulated IRS-aided $K$-user MISO communication scenario comprising of one $M$-antenna BS and one $N$-element IRS. TABLE I: Simulation Parameters Parameters | Values
---|---
Carrier frequency | 2.4 GHz
Transmission bandwidth |  $200$ kHz
Noise power spectral density |  $-170$ dBm/Hz
Path-loss for BS-IRS and IRS-user links (dB) | $40+17.3\lg d$
Path-loss for BS-user link (dB) | $30+31.9\lg d+\zeta$
Penetration loss $\zeta$ due to obstacle | 20 dB
Reflection efficiency of IRS | 0.8
Height of users | 1.5 m
Location of BS | (0, 0, 3m)
Location of IRS | (0, 10m, 3m)

In the simulation, we consider two baseline schemes to benchmark the proposed
scheme.

  * •

Baseline 1 [LMMSE using the protocol in Figure 2(b) [21, 22]]: This curve
illustrates the performance of the LMMSE estimator proposed in [21] and [22].
For simplicity, we assume that the BS-user channels have already been
perfectly estimated by this scheme. The protocol illustrated in Figure 2(b) is
adopted. In additon, for fair comparison, the number of training timeslots is
set as $\lceil\frac{N-1}{K}\rceil+\lceil\frac{N}{M}\rceil$ such that the total
pilot overhead is just slightly higher than that of the proposed scheme.

  * •

Baseline 2 [Bilinear alternating least squares (BALS) algorithm [17]]: In
[17], an iterative algorithm is proposed to estimate ${\bf G}$ and ${\bf
H}_{\rm r}$ by utilizing the PARAFAC decomposition, which adopts the same
channel estimation protocol as Baseline 1. Note that due to the ambiguity
issue (see Lemma 2 or [17, Section IV]), the BALS also cannot exactly
reconstruct ${\bf G}$ and ${\bf H}_{\rm r}$, and the actually estimated
variable is still the cascaded channel.

  * •

Baseline 3 [MAP modification for the BALS in [17]]: In this baseline, we make
a simple modification based on the MAP optimization in this paper to further
enhance the performance of the BALS algorithm in [17] by exploiting the prior
knowledge of the cascaded channels.

  * •

Baseline 4 [Selected On-off protocol [23]]: This curve illustrates the
performance of the estimation algorithm in [23] based on the selected on-off
channel estimation protocol shown in Figure 3. We assume that the covariances
of ${\bf h}_{{\rm u},k}={\rm diag}({\bf h}_{{\rm r},1})^{-1}{\bf h}_{{\rm
r},k}$ for all $k=1,2,\cdots,8$ are available, and the LMMSE estimator in [23,
Section V] is adopted. In addition, we always select user $1$ as the reference
user.

Note that the proposed scheme and Baselines 3 and 4 have the same pilot
overhead, i.e., $K+N+\lceil\frac{N}{M}\rceil(K-1)$. We focus on the evaluation
of the performance of cascaded channel estimation, and use the normalized MSE
(NMSE) as the evaluation metric, which is given by

${\rm{NMSE}}=\frac{\sum_{k=1}^{K}\sum_{m=1}^{M}{\mathbb{E}}\left[\left|{\bf
h}_{{\rm I},k,m}-{\hat{\bf h}}_{{\rm
I},k,m}\right|_{2}^{2}\right]}{\sum_{k=1}^{K}\sum_{m=1}^{M}{\mathbb{E}}\left[\left|{\bf
h}_{{\rm I},k,m}\right|_{2}^{2}\right]}.$ (53)

In addition, based on the proposed protocol, the BS-user channel estimation
can be independent of the cascaded channel estimation by applying the signal
pre-processing as shown in Section III-B. This signal pre-processing operation
will provide a theoretical $3$ dB gain for the BS-user channel estimation
compared to the conventional solution, which shuts down the IRS to estimate
the BS-user channel (see Appendix B), and thus we do not compare the
estimation performance for the BS-user channel in the simulations.

(a) $P_{\rm T}$ vs NMSE

(b) Convergence behavior when $P_{\rm T}=15$ dBm

Figure 7: The NMSE versus transmit power when $M=8$ and $N=32$.

###  VI-B Simulation Results

Figure 7(a) illustrates the NMSE of different schemes with respect to the
transmit power of users, in which the BS adopts a $4\times 2$ uniform planar
array (UPA), and the IRS adopts an $8\times 4$ UPA. Thus, we have $M=8$ and
$N=32$. The BALS algorithm in [17] achieves the worst performance since it
does not exploit the channel prior knowledge. The performance of the LMMSE
using the traditional protocol in [21] and [22] does not vary with the
increase of $P_{\rm T}$ since the main bottleneck is that the number of
training timeslots is smaller than $N$. Moreover, the performance of BALS-MAP
is better than that of LMMSE at a low SNR, but worsens as $P_{\rm T}$
increases since it will reduce to BALS when $P_{\rm T}$ is infinite. Based on
the above observations, we can draw a conclusion that the traditional channel
estimation protocol shown in Figure 2(b) is not effective for exploiting the
common-link structure, and thus we do not consider Baselines 2 and 3 in the
remaining simulations. On the other hand, it is seen that the proposed
protocol with the optimization-based channel estimation algorithm achieves
significant gain compared to all the baselines. In addition, the phase
shifting configuration by solving ${\mathcal{P}}{(\text{B})}$ achieves a more
than 3 dB gain by steering the reflected signals in Stage I to the direction
with a higher SNR compared to the random configuration baseline. Next, in
Figure 7(b), we fix the transmit power $P_{\rm T}$ at 15 dBm and show the
convergence behaviors of the proposed Algorithm 1 for
${\mathcal{P}}{(\text{A})}$. One can see that the proposed algorithm converges
quickly. Note that although the solution without phase-shift optimization
achieves a higher objective value, this does not imply it will have better
performance since the objective functions of the two curves are different due
to them adopting different ${\bm{\vartheta}}$ for the training phase shifting
configuration.

Figure 8: NMSE versus $M$, when $P_{\rm T}=20$ dBm.

In Figure 8, we simulate the performance of different BS antenna numbers $M$
when the BS adopts the uniform linear array and the IRS is still $8\times 4$
UPA. It is seen that the NMSE of all curves increases as $M$ increases since
the ratio of the channel unknowns to the training observations decreases as
$M$ increases. Moreover, the performance gain achieved by the phase shifting
configuration increases as $M$ increases. This is because when $M$ increases,
the number of training timeslots in Stage I of the proposed protocol
decreases, and the probability that the random configuration scheme steers to
the highest SNR direction becomes lower.

Figure 9: NMSE versus $N$, when $P_{\rm T}=20$ dBm. Figure 10: The CCDFs for
random user locations.

In Figure 9, we simulate the NMSE of different schemes for different IRS sizes
$N$. The BS is $4\times 2$ UPA, and the IRS is $N_{1}\times 8$ UPA in which
$N_{1}$ increases from $4$ to $10$. Note that as $N$ increases, the pilot
overhead increases according to the proposed protocol but the ratio of the
channel unknowns to the training observations is almost fixed. It is seen that
the NMSE of Baselines 1 and 4 varies only a little, while the NMSE of the
proposed scheme decreases as $N$ increases. This is because the channel
becomes more correlated as $N$ becomes large, and the proposed scheme has a
better capability of exploiting the channel prior knowledge.

Finally, we investigate the impact of user locations on the estimation
performance. In particular, we fix $N=8\times 4$ and $M=4\times 2$, and
generate $100$ snapshots for random user locations. For each snapshot, we
further generate $1000$ channel realizations with independent small-scale
fading to reduce the impact of other system parameters. Figure 10 plots the
complementary cumulative distribution functions (CCDFs) of the NMSE for
different snapshots. One can see that the performance gains of the proposed
scheme are irrespective of user locations. In addition, we further increase
$P_{\rm T}=40$ dBm for Baseline 4 (i.e., the selected-on-off-protocol-based
scheme [23]) such that it achieves a similar average NMSE to the proposed
scheme with $P_{\rm T}=20$ dBm. However, Baseline 4 achieves a much worse
outage performance. This is because the performance of the selected-on-off-
protocol-based scheme [23] highly depends on the channel quality of the
reference user, while the proposed scheme is much more robust since it does
not require selecting one reference user.

##  VII Conclusion

In this paper, we proposed a novel always-ON channel estimation protocol for
uplink cascaded channel estimation in IRS-assisted MU-MISO systems. In
contrast to the existing schemes, the pilot overhead required by the proposed
protocol is greatly reduced by exploiting the common-link structure. Based on
the protocol, we formulated an optimization-based joint channel estimation
problem that utilizes the combined statistical information of the cascaded
channels, and then we proposed an alternating optimization algorithm to solve
the problem with the local optimum solution. In addition, we optimized the
phase shifting configuration in the proposed protocol, which may further
enhance the channel estimation performance. The simulation results
demonstrated that the proposed protocol using the optimization based joint
channel estimation algorithm achieves a more than $15$ dB gain compared to the
benchmark. In addition, the proposed optimized phase shifting configuration
achieves a more than $3$ dB gain compared to the random configuration scheme.

##  Appendix A Proof of Lemma 1

Define the virtual reference channel by ${\bf H}_{\rm v}={\bf G}^{\rm T}{\rm
diag}({\bf H}_{\rm r}\bar{\bf x})$. Based on (18), we have

$\displaystyle\left[\bar{\bf r}_{1},\cdots,\bar{\bf r}_{N}\right]$
$\displaystyle={\bf H}_{\rm v}{\bm{\Phi}}.$ (54)

Thus ${\bf H}_{\rm v}$ can be perfectly estimated by:

$\displaystyle{\bf H}_{\rm v}=\left[\bar{\bf r}_{1},\cdots,\bar{\bf
r}_{N}\right]{\bm{\Phi}}^{-1}.$ (55)

We further define the $K$ relative channels by ${\bf h}_{{\rm A},k}={\rm
diag}({\bf H}_{\rm r}\bar{\bf x})^{-1}{\bf h}_{{\rm r},k}$, and ${\bf H}_{\rm
A}=[{\bf h}_{{\rm A},1},{\bf h}_{{\rm A},2},\cdots,{\bf h}_{{\rm A},K}]$. Then
the cascaded channels become ${\bf H}_{{\rm I},k}={\bf H}_{\rm v}{\rm
diag}({\bf h}_{{\rm A},k})$. Therefore, the remaining task is to perfectly
estimate ${\bf H}_{\rm A}$.

Based on the assumption on $\bf G$ and ${\bf H}_{\rm r}$, we have ${\bf
H}_{\rm v}={\bf F}_{\rm B}{\ddot{\bf G}}^{\rm T}{\bf F}_{\rm R}^{\rm T}{\rm
diag}({\bf F}_{\rm R}{\ddot{\bf H}}_{\rm r}\bar{\bf x})$ where ${\ddot{\bf
H}}_{\rm r}=[{\ddot{\bf h}}_{{\rm r},1},\cdots,{\ddot{\bf h}}_{{\rm r},K}]$.
Define ${\bf V}=[{\bf v}_{1},{\bf v}_{2},\cdots,{\bf v}_{M}]$ which is given
by

$\displaystyle{\bf v}_{m}={\rm diag}({\bf F}_{\rm R}{\ddot{\bf H}}_{\rm
r}\bar{\bf x}){\bf F}_{\rm R}{\ddot{\bf g}}_{m},$ (56)

and thus ${\bf H}_{\rm v}={\bf F}_{\rm B}{\bf V}^{\rm T}$. Using the
independence of ${\ddot{\bf g}}_{m}$ and ${\ddot{\bf h}}_{{\rm r},k}$, we have
${\mathbb{E}}[{\bf v}_{i}{\bf v}_{j}^{\rm H}]={\bm{0}}$. Since all ${\bf
v}_{m}$ ($m=1,2,\cdots,M$) follow joint multivariate normal distribution,
${\bf v}_{m}$ for all $m$ are pairwise independent to each other. In addition,
since ${\bf C}_{m}^{(k)}={\mathbb{E}}\left[{\bf h}_{{\rm I},k,m}{\bf h}_{{\rm
I},k,m}^{\rm H}\right]$ is full-rank for all $k$ and $m$, ${\mathbb{E}}[{\bf
v}_{m}{\bf v}_{m}^{\rm H}]$ is also full-rank for all $m$ with properly-
designed $\bar{\bf x}$.

Next, using ${\bf R}_{\ell}$ in (16) and (17), we have

$\displaystyle\tilde{\bf R}_{\ell}$ $\displaystyle={\bf F}_{\rm B}^{-1}{\bf
R}_{\ell}{\bf X}^{-1}$ (57) $\displaystyle={\bf F}_{\rm B}^{-1}{\bf G}^{\rm
T}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm r}$ $\displaystyle={\bf F}_{\rm
B}^{-1}{\bf H}_{\rm v}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm A}$
$\displaystyle={\bf V}^{\rm T}{\rm diag}({\bm{\theta}}_{\ell}){\bf H}_{\rm
A}.$

Stacking all $\tilde{\bf R}_{\ell}$, we have

$\displaystyle\begin{bmatrix}\tilde{\bf R}_{1}\\\ \vdots\\\ \tilde{\bf
R}_{L_{1}}\end{bmatrix}=\begin{bmatrix}{\bf V}^{\rm T}{\rm
diag}({\bm{\theta}}_{1})\\\ \vdots\\\ {\bf V}^{\rm T}{\rm
diag}({\bm{\theta}}_{L_{1}})\end{bmatrix}{\bf H}_{\rm A}.$ (58)

Define ${\bm{\Psi}}=[{\rm diag}({\bm{\theta}}_{1}){\bf V},\cdots,{\rm
diag}({\bm{\theta}}_{L_{1}}){\bf V}]^{\rm T}$. Now we need to prove
${\text{rank}}({\bm{\Psi}})={\rm{min}}\\{N,ML_{1}\\}$ with probability one,
and the whole proof is completed.

By permuting the columns of ${\bm{\Psi}}$, we have a new matrix
$\bar{\bm{\Psi}}=[{\rm diag}({\bf v}_{1})\bar{\bm{\Phi}},\cdots,{\rm
diag}({\bf v}_{M})\bar{\bm{\Phi}}]^{\rm T}$ where
$\bar{\bm{\Phi}}=[{\bm{\theta}}_{1},\cdots,{\bm{\theta}}_{L_{1}}]$. Then, it
is equivalent to prove that
${\text{rank}}(\bar{\bm{\Psi}})={\rm{min}}\\{N,ML_{1}\\}$ with probability
one. We prove it by induction. Define $\bar{\bm{\Psi}}_{m}=[{\rm diag}({\bf
v}_{1})\bar{\bm{\Phi}},\cdots,{\rm diag}({\bf v}_{m})\bar{\bm{\Phi}}]^{\rm
T}$. Since $\bar{\bm{\Phi}}$ is semi-orthogonal,
${\text{rank}}(\bar{\bm{\Psi}}_{1})=L_{1}$ with probability one. Let
${\text{rank}}(\bar{\bm{\Psi}}_{m-1})={\rm{min}}\\{N,(m-1)L_{1}\\}$ with
probability one, and the rest task is to prove
${\text{rank}}(\bar{\bm{\Psi}}_{m})={\rm{min}}\\{N,mL_{1}\\}$ with probability
one. We prove it by contradiction. Consider the case when
${\text{rank}}(\bar{\bm{\Psi}}_{m-1})=(m-1)L_{1}$ which is smaller than $N$
but ${\text{rank}}(\bar{\bm{\Psi}}_{m})<{\rm{min}}\\{N,mL_{1}\\}$. Since ${\rm
diag}({\bf v}_{m}){\bm{\theta}}_{i}$ and ${\rm diag}({\bf
v}_{m}){\bm{\theta}}_{j}$ are orthogonal for $i\neq j$, there shall exists
$\bf x$ which satisfies:

$\displaystyle[{\rm diag}({\bf v}_{1})\bar{\bm{\Phi}},\cdots,{\rm diag}({\bf
v}_{m-1})\bar{\bm{\Phi}}]{\bf x}={\rm diag}({\bf v}_{m}){\bm{\theta}}_{\ell},$
(59)

for some $1\leq\ell\leq L_{1}$ such that
${\text{rank}}(\bar{\bm{\Psi}}_{m})<{\rm{min}}\\{N,mL_{1}\\}$ is true.
However, since ${\bf v}_{m}$ is independent to ${\bf v}_{1},{\bf
v}_{2},\cdots,{\bf v}_{m-1}$ with full-rank covariance matrix, equation (59)
is inconsistent (which has no solution) with probability one, and the whole
proof is finished.

##  Appendix B Estimation on the BS-User Channel

Based on ${\bm{\theta}}_{0}=-{\bm{\theta}}_{1}$, we have

$\displaystyle{\bf R}_{0}$ $\displaystyle=\frac{1}{2}\left({{\bf Y}}_{0}+{{\bf
Y}}_{1}\right){\bf X}^{-1}$ (60) $\displaystyle={\bf H}_{\rm d}+\tilde{\bf
Z}_{0},$

where $\tilde{\bf Z}_{0}$ is the noise matrix consisting of $MK$ i.i.d.
complex Gaussian variables following
${\cal{CN}}({0},\frac{1}{2K}\sigma_{0}^{2})$. Let ${\bf r}_{{0},k}$ be the
$k$-th column of ${\bf R}_{0}$. The BS-user direct channel for the $k$-th user
can be estimated by the LMMSE estimator [41]:

$\displaystyle\hat{\bf h}_{{\rm d},k}={\bf C}_{{\rm d},k}\left({\bf C}_{{\rm
d},k}+\frac{\sigma_{0}^{2}}{2K}{\bf I}_{M}\right)^{-1}{\bf r}_{0,k},$ (61)

where ${\bf C}_{{\rm d},k}={\mathbb{E}}\left[{\bf h}_{{\rm d},k}{\bf h}_{{\rm
d},k}^{\rm H}\right]$ is the covariance matrix of the direct channel from the
BS to the $k$-th user. Note that, as shown in (60) and (61), the proposed
protocol may achieve a $3$ dB performance gain compared to the existing works,
[21, 22, 23], on the estimation of the BS-user channels since it exploits
doubled observation samples.

##  Appendix C Proof of Lemma 2

The objective function of ${\mathcal{P}}{(\text{A})}$ can be denoted by the
function of the cascaded channel coefficients $\\{{\bf h}_{{\rm I},k,m}\\}$,
as follows:

$\displaystyle f_{{\rm A}}({\bf H}_{\rm g},{\bf H}_{\rm u})=f_{{\rm
A}}(\\{{\bf h}_{{\rm I},k,m}\\})$ (62)
$\displaystyle=-\sum_{\ell=1}^{L_{1}}\sum_{k=1}^{K}{\frac{1}{\sigma_{\ell}^{2}}\left\|\tilde{{\bf
r}}_{\ell,k}-\sum_{m=1}^{M}{\bf h}_{{\rm I},k,m}^{\rm
T}{\bm{\theta}}_{\ell}\right\|^{2}}$
$\displaystyle\quad-\sum_{\ell=L_{1}+1}^{N}{\frac{1}{\bar{\sigma}_{\ell}^{2}}\left\|\bar{\bf
r}_{\ell}-\sum_{k=1}^{K}\sum_{m=1}^{M}\bar{x}_{k}{\bf h}_{{\rm I},k,m}^{\rm
T}{\bm{\theta}}_{\ell}\right\|^{2}}$
$\displaystyle\quad-\sum_{m=1}^{M}\sum_{k=1}^{K}{\bf h}_{{\rm I},k,m}^{\rm
H}{{\bf C}_{m}^{(k)}}^{-1}{\bf h}_{{\rm I},k,m}.$

Therefore, if one optimal solution $\left\\{{\bf H}_{\rm g}^{\star},{\bf
H}_{\rm u}^{\star}\right\\}\in{\cal A}$, any $\left\\{{\bf H}_{\rm g},{\bf
H}_{\rm u}\right\\}$ pair in set ${\cal A}$ is an optimal solution of
${\mathcal{P}}{(\text{A})}$, and the lemma is proved.

##  Appendix D Proof of Lemma 3

Denote $\ddot{\bf h}_{\rm u}={\rm{vec}}({\bf H}_{\rm u})$. The objective
function $f_{{\rm u}}$ in (26) is given by

$\displaystyle f_{{\rm u}}(\ddot{\bf h}_{\rm
u})=\sum_{\ell=1}^{L_{1}}{\frac{1}{\sigma_{\ell}^{2}}\left\|{\rm{vec}}(\tilde{{\bf
R}}_{\ell})-\left({\bf I}_{K}\otimes{\bf D}_{\ell}\right)\ddot{\bf h}_{\rm
u}\right\|^{2}_{2}}$ (63)
$\displaystyle\quad+\sum_{\ell=L_{1}+1}^{N}{\frac{1}{\bar{\sigma}_{\ell}^{2}}\left\|\bar{\bf
r}_{\ell}-\left({\bar{\bf x}}^{\rm T}\otimes{\bf D}_{\ell}\right)\ddot{\bf
h}_{\rm u}\right\|^{2}_{2}}+\ddot{\bf h}_{\rm u}^{\rm H}{\bf C}_{{\rm
u}}\ddot{\bf h}_{\rm u},$

where ${\bf C}_{{\rm u}}={\rm{blkdiag}}({\bf C}_{{\rm u},1},{\bf C}_{{\rm
u},2},\cdots,{\bf C}_{{\rm u},K})$. The second order derivative of $f_{{\rm
u}}(\ddot{\bf h}_{\rm u})$ is given by

$\displaystyle\frac{\partial^{2}f_{{\rm u}}(\ddot{\bf h}_{\rm
u})}{\partial\ddot{\bf h}_{\rm u}\partial\ddot{\bf h}_{\rm u}^{\rm H}}$
$\displaystyle=2{\bf C}_{{\rm u}}+2{\bf
I}_{K}\otimes\left(\sum_{\ell=1}^{L_{1}}\frac{1}{\sigma_{\ell}^{2}}{\bf
D}_{\ell}^{\rm H}{\bf D}_{\ell}\right)$ $\displaystyle\quad+2\left(\bar{\bf
x}^{\ast}\bar{\bf x}^{\rm
T}\right)\otimes\left(\sum_{\ell=L_{1}+1}^{N}\frac{1}{\bar{\sigma}_{\ell}^{2}}|\bar{x}_{k}|^{2}{\bf
D}_{\ell}^{\rm H}{\bf D}_{\ell}\right),$

which is a Hermitian positive semi-definite matrix. Thus, the lemma is proved.

##  Appendix E Proof of Lemma 4

The second order derivative of $f_{{\rm g},m}({\bf h}_{{\rm g},m})$ in (33) is
given by

$\displaystyle\frac{\partial^{2}f_{{\rm g},m}({\bf h}_{{\rm
g},m})}{\partial{\bf h}_{{\rm g},m}\partial{\bf h}_{{\rm g},m}^{\rm H}}=2{\bf
C}_{{\rm
g},m}^{-1}+2\sum_{k=1}^{K}\sum_{\ell=1}^{L_{1}}\frac{1}{\sigma_{\ell}^{2}}{\bf
b}_{\ell,k}^{\ast}{\bf b}_{\ell,k}^{\rm T}$ (64)
$\displaystyle\quad+2\sum_{\ell=L_{1}+1}^{N}\frac{1}{\bar{\sigma}_{\ell}^{2}}\left(\sum_{k=1}^{K}\bar{x}_{k}{\bf
b}_{\ell,k}^{\rm T}\right)^{\rm H}\left(\sum_{k=1}^{K}\bar{x}_{k}{\bf
b}_{\ell,k}^{\rm T}\right),$

which is a Hermitian positive semi-definite matrix. Thus, $f_{{\rm g},m}({\bf
h}_{{\rm g},m})$ is a convex quadratic function of ${\bf h}_{{\rm g},m}$.
Then, the objective function $\sum_{m=1}^{M}f_{{\rm g},m}({\bf h}_{{\rm
g},m})$ is a convex quadratic function of $\\{{\bf h}_{{\rm g},1},{\bf
h}_{{\rm g},2},\cdots,{\bf h}_{{\rm g},M}\\}$.

## References

  * [1] X. Tan, Z. Sun, D. Koutsonikolas, and J. M. Jornet, “Enabling indoor mobile millimeter-wave networks based on smart reflect-arrays,” in _Proc. IEEE INFOCOM_ , Apr. 2018, pp. 270–278.
  * [2] L. Li, T. J. Cui, W. Ji, S. Liu, J. Ding, X. Wan, Y. B. Li, M. Jiang, C.-W. Qiu, and S. Zhang, “Electromagnetic reprogrammable coding-metasurface holograms,” _Nature Commun._ , vol. 8, no. 1, p. 197, 2017.
  * [3] E. Björnson, Ö. Ö, and E. G. Larsson, “Reconfigurable intelligent surfaces: Three myths and two critical questions,” _IEEE Commun. Mag._ , vol. 58, no. 12, pp. 90–96, 2020.
  * [4] J. C. B. Garcia, A. Sibille, and M. Kamoun, “Reconfigurable intelligent surfaces: Bridging the gap between scattering and reflection,” _IEEE J. Sel. Areas Commun._ , vol. 38, no. 11, pp. 2538–2547, 2020.
  * [5] M. D. Renzo, A. Zappone, M. Debbah, M. S. Alouini, C. Yuen, J. d. Rosny, and S. Tretyakov, “Smart radio environments empowered by reconfigurable intelligent surfaces: How it works, state of research, and the road ahead,” _IEEE J. Sel. Areas Commun._ , vol. 38, no. 11, pp. 2450–2525, 2020.
  * [6] Q. Wu and R. Zhang, “Towards smart and reconfigurable environment: Intelligent reflecting surface aided wireless network,” _IEEE Commun. Mag._ , vol. 58, no. 1, pp. 106–112, 2019.
  * [7] C. Liaskos, S. Nie, A. Tsioliaridou, A. Pitsillides, S. Ioannidis, and I. Akyildiz, “A new wireless communication paradigm through software-controlled metasurfaces,” _IEEE Commun. Mag._ , vol. 56, no. 9, pp. 162–169, 2018.
  * [8] M. Di Renzo, M. Debbah, D.-T. Phan-Huy, A. Zappone, M.-S. Alouini, C. Yuen, V. Sciancalepore, G. C. Alexandropoulos, J. Hoydis, H. Gacanin _et al._ , “Smart radio environments empowered by reconfigurable AI meta-surfaces: an idea whose time has come,” _EURASIP J. Wireless Commun. and Netw._ , vol. 2019, no. 1, p. 129, 2019.
  * [9] Q. Wu and R. Zhang, “Intelligent reflecting surface enhanced wireless network via joint active and passive beamforming,” _IEEE Trans. Wireless Commun._ , vol. 18, no. 11, pp. 5394–5409, 2019.
  * [10] H. Guo, Y. Liang, J. Chen, and E. G. Larsson, “Weighted sum-rate maximization for reconfigurable intelligent surface aided wireless networks,” _IEEE Trans. Wireless Commun._ , vol. 19, no. 5, pp. 3064–3076, 2020.
  * [11] X. Mu, Y. Liu, L. Guo, J. Lin, and N. Al-Dhahir, “Exploiting intelligent reflecting surfaces in NOMA networks: Joint beamforming optimization,” _IEEE Trans. Wireless Commun._ , vol. 19, no. 10, pp. 6884–6898, 2020.
  * [12] G. Zhou, C. Pan, H. Ren, K. Wang, and A. Nallanathan, “A framework of robust transmission design for IRS-aided MISO communications with imperfect cascaded channels,” _IEEE Trans. Signal Process._ , vol. 68, pp. 5092–5106, 2020.
  * [13] S. Lin, B. Zheng, G. C. Alexandropoulos, M. Wen, M. D. Renzo, and F. Chen, “Reconfigurable intelligent surfaces with reflection pattern modulation: Beamforming design and performance analysis,” _IEEE Trans. Wireless Commun._ , vol. 20, no. 2, pp. 741–754, 2021.
  * [14] T. L. Jensen and E. De Carvalho, “An optimal channel estimation scheme for intelligent reflecting surfaces based on a minimum variance unbiased estimator,” in _Proc. IEEE ICASSP_ , 2020, pp. 5000–5004.
  * [15] C. You, B. Zheng, and R. Zhang, “Channel estimation and passive beamforming for intelligent reflecting surface: Discrete phase shift and progressive refinement,” _IEEE J. Sel. Areas Commun._ , vol. 38, no. 11, pp. 2604–2620, 2020.
  * [16] C. Zhang, Y. Jing, Y. Huang, and X. You, “Interleaved training for intelligent surface-assisted wireless communications,” _IEEE Signal Process. Lett._ , vol. 27, pp. 1774–1778, 2020.
  * [17] G. T. De Araujo, A. Dealmeida, and R. Boyer, “Channel estimation for intelligent reflecting surface assisted MIMO systems: A tensor modeling approach,” _IEEE J. Sel. Topics Signal Process._ , vol. 15, no. 3, pp. 789–802, 2021.
  * [18] D. Mishra and H. Johansson, “Channel estimation and low-complexity beamforming design for passive intelligent surface assisted MISO wireless energy transfer,” in _Proc. IEEE Int. Conf. Acoust., Speech, Signal Process._ , 2019, pp. 4659–4663.
  * [19] A. M. Elbir, A. Papazafeiropoulos, P. Kourtessis, and S. Chatzinotas, “Deep channel learning for large intelligent surfaces aided mm-wave massive MIMO systems,” _IEEE Wireless Commun. Lett._ , vol. 9, no. 9, pp. 1447–1451, 2020.
  * [20] Z. Zhou, N. Ge, Z. Wang, and L. Hanzo, “Joint transmit precoding and reconfigurable intelligent surface phase adjustment: A decomposition-aided channel estimation approach,” _IEEE Trans. Commun._ , vol. 69, no. 2, pp. 1228–1243, 2021.
  * [21] N. K. Kundu and M. R. McKay, “Channel estimation for reconfigurable intelligent surface aided MISO communications: From LMMSE to deep learning solutions,” _IEEE Open J. Commun. Soc._ , vol. 2, pp. 471–487, 2021\.
  * [22] H. Alwazani, A. Kammoun, A. Chaaban, M. Debbah, and M.-S. Alouini, “Intelligent reflecting surface-assisted multi-user MISO communication: Channel estimation and beamforming design,” _IEEE Open J. Commun. Soc._ , vol. 1, pp. 661–680, 2020.
  * [23] Z. Wang, L. Liu, and S. Cui, “Channel estimation for intelligent reflecting surface assisted multiuser communications: Framework, algorithms, and analysis,” _IEEE Trans. Wireless Commun._ , vol. 19, no. 10, pp. 6607–6620, 2020.
  * [24] C. Hu, L. Dai, S. Han, and X. Wang, “Two-timescale channel estimation for reconfigurable intelligent surface aided wireless communications,” _IEEE Trans. Commun._ , pp. 1–1, 2021.
  * [25] F. Liu, O. Tsilipakos, A. Pitilakis, A. C. Tasolamprou, M. S. Mirmoosa, N. V. Kantartzis, D.-H. Kwon, M. Kafesaki, C. M. Soukoulis, and S. A. Tretyakov, “Intelligent metasurfaces with continuously tunable local surface impedance for multiple reconfigurable functions,” _Physical Rev. Appl._ , vol. 11, no. 4, p. 044024, 2019.
  * [26] B. W. Frazier, T. M. Antonsen Jr, S. M. Anlage, and E. Ott, “Wavefront shaping with a tunable metasurface: Creating cold spots and coherent perfect absorption at arbitrary frequencies,” _Physical Rev. Research_ , vol. 2, no. 4, p. 043422, 2020.
  * [27] M. F. Imani, D. R. Smith, and P. del Hougne, “Perfect absorption in a metasurface-programmable complex scattering enclosure,” _arXiv preprint arXiv:2003.01766_ , 2020.
  * [28] B. Zheng, C. You, and R. Zhang, “Efficient channel estimation for double-IRS aided multi-user MIMO system,” _IEEE Trans. Commun._ , vol. 69, no. 6, pp. 3818–3832, 2021.
  * [29] Z. He and X. Yuan, “Cascaded channel estimation for large intelligent metasurface assisted massive MIMO,” _IEEE Wireless Commun. Lett._ , vol. 9, no. 2, pp. 210–214, 2020.
  * [30] X. Wei, D. Shen, and L. Dai, “Channel estimation for RIS assisted wireless communications: Part II-an improved solution based on double-structured sparsity,” _IEEE Commun. Lett._ , 2021.
  * [31] P. Wang, J. Fang, H. Duan, and H. Li, “Compressed channel estimation for intelligent reflecting surface-assisted millimeter wave systems,” _IEEE Signal Process. Lett._ , vol. 27, pp. 905–909, 2020.
  * [32] K. Ardah, S. Gherekhloo, and M. Haardt, “TRICE: A channel estimation framework for RIS-aided millimeter-wave MIMO systems,” _IEEE Signal Process. Lett._ , vol. 28, pp. 513–517, 2021.
  * [33] H. Liu, X. Yuan, and Y. J. A. Zhang, “Matrix-calibration-based cascaded channel estimation for reconfigurable intelligent surface assisted multiuser MIMO,” _IEEE J. Sel. Areas Commun._ , vol. 38, no. 11, pp. 2621–2636, 2020\.
  * [34] C. Liu, X. Liu, D. W. K. Ng, and J. Yuan, “Deep residual learning for channel estimation in intelligent reflecting surface-assisted multi-user communications,” _IEEE Trans. Wireless Commun._ , 2021.
  * [35] P. Cai, J. Zong, X. Luo, Y. Zhou, S. Chen, and H. Qian, “Downlink channel tracking for intelligent reflecting surface-aided FDD MIMO systems,” _IEEE Trans. Veh. Technol._ , vol. 70, no. 4, pp. 3341–3353, 2021.
  * [36] Z. Mao, M. Peng, and X. Liu, “Channel estimation for reconfigurable intelligent surface assisted wireless communication systems in mobility scenarios,” _China Commun._ , vol. 18, no. 3, pp. 29–38, 2021.
  * [37] Y. Xu and W. Yin, “A block coordinate descent method for regularized multiconvex optimization with applications to nonnegative tensor factorization and completion,” _SIAM J. Imag. Sci._ , vol. 6, no. 3, pp. 1758–1789, 2013.
  * [38] A. Krishnamoorthy and D. Menon, “Matrix inversion using Cholesky decomposition,” in _Proc. IEEE SPA_ , 2013, pp. 70–72.
  * [39] M. Razaviyayn, “Successive convex approximation: Analysis and applications,” _Ph.D. dissertation, Univ. of Minnesota_ , 2014.
  * [40] _Study on channel model for frequencies from 0.5 to 100 GHz_. 3GPP TR 38.901 version 14.0.0 Release 14, Nov. 2020\.
  * [41] S. M. Kay, _Fundamentals of statistical signal processing_. Prentice Hall PTR, 1993.