# Joint DL and UL Channel Estimation for Millimeter Wave MIMO Systems Using Tensor Modeling

Authors: Paulo R. B. Gomes; André L. F. de Almeida; João Paulo C. L. da Costa; Rafael T. de Sousa
Journal: Wireless Communications and Mobile Computing
 (2019)

---

## Abstract

In this paper, we address the problem of joint downlink (DL) and uplink (UL) channel estimation for millimeter wave (mmWave) multiple-input multiple-output (MIMO) systems. Assuming a closed-loop and multifrequency-based channel training framework in which pilot signals received by multiple antenna mobile stations (MSs) are coded and spread in the frequency domain via multiple adjacent subcarriers, we propose two tensor-based semiblind receivers by capitalizing on the multilinear structure and sparse feature of the received signal at the BS equipped with a hybrid analog-digital beamforming (HB) architecture. As a first processing stage, the joint estimation of the compressed DL and UL channel matrices can be obtained in an iterative way by means of an alternating least squares (ALS) algorithm that capitalizes on a parallel factors model for the received signals. Alternatively, for more restricted scenarios, a closed-form solution is also proposed. From the estimated effective channel matrices, the users’ channel parameters such as angles of departure (AoD), angles of arrival (AoA), and path gains are then estimated in a second processing stage by solving independent compressed sensing (CS) problems (one for each MS). In contrast to the classical approach in the literature, in which the DL and UL channel estimation problems are usually considered as two separate problems, our idea is to jointly estimate both the DL and UL channels as a single problem by concentrating most of the processing burden for channel estimation at the BS side. Simulation results demonstrate that the proposed receivers achieve a performance close to the classical approach that is applied on DL and UL communication links separately, with the advantage of avoiding complex computations for channel estimation at the MS side as well as dedicated feedback channels for each MS, which are attractive features for massive MIMO systems.


---

## Body

## 1. Introduction

In recent years, millimeter wave (mmWave) massive multiple-input multiple-output (MIMO) technology has been a subject of increasing interest in both academia and industry for future wireless standards due to its great potential to provide substantial gains in data rates and energy efficiency. However, due to the severe path loss over the mmWave frequency bands, large antenna arrays should be deployed at the base station (BS) and mobile stations (MSs) to provide sufficient beamforming gain in mmWave MIMO scenarios [1]. In this context, the implementation of fully digital beamforming architectures becomes prohibitive due its expensive cost, hardware constraints, and power consumption of high-resolution analog-to-digital converters (ADC) and digital-to-analog converters (DAC) per antenna port [2, 3]. To overcome these practical limitations, hybrid analog-digital beamforming (HB) architectures that split the signal processing between analog and digital domains using a reduced number of radio frequency (RF) chains (assumed to be smaller than the number of antennas) have been investigated [4, 5]. In the HB architectures, the digital part performs baseband processing using microprocessors, while the analog part can be implemented at the RF domain using different analog approaches such as phase-shifter networks [6], switches [3, 7], or lenses [8].To fully benefit from the beamforming gains in mmWave MIMO systems, an accurate channel estimation is crucial to realize the hybrid precoding designs in which the analog part is used to improve the signal power, while the digital part is designed to suppress interuser interferences [9–11]. For this purpose, several channel estimation techniques such as [12–16] have been proposed. The authors of [12] proposed an iterative method based on the least squares estimation (LSE) concept and sparse message passing (SMP) algorithm. In this method, the location of nonzero entries of the channel vector is detected through the SMP, while the LSE is used for estimating the channel coefficients at each iteration. Zhu et al. [13] proposed an auxiliary beam pair design for mmWave channel estimation in which the best auxiliary beam pair is fed back to the transmitter via a feedback channel. The method proposed by Ghauch et al. [14] consists of a subspace-based approach that exploits the channel reciprocity in time-division duplexing (TDD) MIMO systems for hybrid precoding design. It iteratively estimates the dominant singular modes of the channel instead of the entire channel. In contrast, the works [15, 16] explore the angular sparsity of mmWave channels and use compressed sensing (CS) theory to estimate only the channel parameters from which the mmWave channel can be reconstructed. In [15], the angular spreads over the angle of arrival (AoA) and angle of departure (AoD) are considered in the channel modeling, while the low-rank structure of the channel is exploited to reduce the number of samples needed to recover the mmWave channel. Similarly, but disregarding the angular spreads in the spatial domains, the method in [16] uses the 2D unitary ESPRIT algorithm for spatial parameters estimation, while the path gains are estimated by means of the LS criterion.The researchers [17–21] have proposed CS-based and tensor-based channel estimators for mmWave MIMO systems, respectively. They assume the conventional channel training framework, where the DL and UL channel estimation problems are treated separately (as two decoupled procedures at the MS and BS, respectively). In particular, for frequency division duplexing (FDD) systems, where channel estimation is usually carried at the power-limited MS side, computational complexity plays a significant role due to the large number of channel coefficients to be estimated. An interesting approach to deal with this problem exploits the poor scattering nature of the mmWave channels via CS techniques [22, 23]. For instance, in [17, 18], the intrinsic sparse feature of the mmWave channel is exploited and CS-based channel estimation algorithms are formulated. However, the adaptive algorithm proposed in [17] can be applied to estimate the DL or UL channel separately. In [19], a layered pilot transmission scheme is proposed to UL channel estimation, while [20, 21] exploit the DL communication of wideband mmWave channels. The main idea in [21] is to divide the overall channel estimation problem into three smaller CS subproblems via tensor-based modeling to estimate the channel parameters (AoDs, AoAs, and delays) with less computational complexity. The system model is formulated as a parallel factors (PARAFAC) decomposition [24], and tensor-based algorithms combined with CS tools are proposed to solve the channel estimation problem. However, Zhou et al. [19] only considered the UL channel estimation, while Zhou et al. and Araújo and de Almeida [20, 21] focused on the DL case. On the other hand, Shen et al. [25] proposed a closed-loop based training framework, where the MSs directly feed the received pilots back to the BS without channel estimation. Therein, a simplified approach is adopted, in which the UL channel is modeled as an additive white Gaussian noise (AWGN) term. Moreover, in [25], only the DL channel could be estimated. Different from [25], we are interested in a joint estimation of the UL and DL channels, and a more realistic multipath channel model that fits to mmWave MIMO scenario and HB architecture is considered.In this paper, we study the problem of joint DL and UL channel estimation in the context of mmWave MIMO systems that employ HB architectures. Initially, we propose a novel closed-loop and multifrequency-based channel training framework in which the pilot signals received by multiple MSs are coded and spread in the frequency domain and then fed back to the BS over the same UL resources. Making use of the proposed framework for channel estimation, the received closed-loop signal at the BS can be modeled as a three-way array (i.e., a third-order tensor) that follows a PARAFAC model. By capitalizing both on the multidimensional and sparse structures of the received signal, we propose two tensor-based semiblind receivers for joint DL and UL channel estimation. The first receiver is an iterative solution based on the alternating least squares (ALS) algorithm [26]. The second is a closed-form solution based on the least squares Khatri-Rao factorization (LS-KRF) algorithm [27]. In the proposed receivers, we first obtain joint estimates of the compressed DL and UL channel matrices. Then, we exploit the sparse representation of the DL and UL channels to individually recover the channel parameters (AoDs, AoAs, and path gains) of each user via CS-based techniques. The proposed framework allows concentrating most of the processing burden for channel estimation at the BS side, i.e., avoiding unnecessary computational overhead for channel estimation at the power-limited MS side. Our simulation results reveal that the proposed receivers achieve a performance close to the classical framework that treats the estimation of DL and UL channels as separate problems.In summary, the main contributions of this paper can be listed as follows:(i)
We propose a novel closed-loop and multifrequency-based channel training framework for channel estimation that focuses jointly on the DL and UL communication links. The proposed framework concentrates the processing burden for joint channel estimation at the BS, avoiding processing with high computational cost at the MS side.(ii)
We show that, by making use of the proposed framework for channel estimation, the received closed-loop signal can be modeled as a third-order tensor that follows a PARAFAC model. Then, we formulate two tensor-based semiblind receivers (iterative and closed-form ones) for joint DL and UL channel estimation by capitalizing on a tensor structure of the received closed-loop signal.(iii)
We study the identifiability issues under which the DL and UL channel matrices can be jointly and uniquely estimated using the proposed receivers. Useful lower bounds on the number of subcarriers required to accomplish the joint channel estimation are derived.The rest of this paper is structured as follows: In Section2, we provide as a presentation complement some important tensor definitions, tensor algebra operations, and a brief overview on the PARAFAC decomposition. In Section 3, we present the proposed channel training framework, and the system and channel models. Section 4 formulates the two proposed tensor-based semiblind receivers for joint DL and UL channel estimation. In Section 5, we analyze the identifiability conditions of the proposed receivers. Simulation results are provided in Section 6. Conclusions and perspectives for future work are drawn in Section 7.


### 1.1. Notation and Properties

Scalars, column vectors, matrices, and tensors are denoted by nonbold lowercase lettersa, bold lowercase letters a, bold uppercase letters A, and calligraphic uppercase letters A, respectively. The superscripts ⋅T, ⋅∗, ⋅H, and ⋅† denote the transpose, complex conjugate, conjugate transpose, and pseudoinverse operations. ⋅F represents the Frobenius norm of a matrix or tensor. The i,r-th entry of A is denoted by Ai,r. The operator diaga converts a into a diagonal matrix, while DiA consists in a diagonal matrix formed by the i-th row of A. vecA converts A to a vector a by stacking its columns on top of each other, while unvecI×Ra converts a∈ℂIR to a matrix A∈ℂI×R. vecdA converts the diagonal elements of A into a vector. ∘ denotes the outer product operator. The Kronecker and Khatri-Rao products are denoted by ⊗ and ⋄, respectively. The Khatri-Rao product between the matrices A=a1,…,aR∈ℂI×R and B=b1,…,bR∈ℂJ×R corresponds to a column-wise Kronecker product, i.e.,(1)A⋄B=a1⊗b1,a2⊗b2,…,aR⊗bR∈ℂIJ×R.We shall make use of the following two properties of the Kronecker and Khatri-Rao products:(2)AC⋄BD=A⊗BC⋄D,(3)vecABCT=C⋄AvecdB,(4)a⊗b=vecb∘a∈ℂIJ,where B is assumed to be a diagonal matrix in (3). In both cases, the matrices have compatible dimensions.


## 1.1. Notation and Properties

Scalars, column vectors, matrices, and tensors are denoted by nonbold lowercase lettersa, bold lowercase letters a, bold uppercase letters A, and calligraphic uppercase letters A, respectively. The superscripts ⋅T, ⋅∗, ⋅H, and ⋅† denote the transpose, complex conjugate, conjugate transpose, and pseudoinverse operations. ⋅F represents the Frobenius norm of a matrix or tensor. The i,r-th entry of A is denoted by Ai,r. The operator diaga converts a into a diagonal matrix, while DiA consists in a diagonal matrix formed by the i-th row of A. vecA converts A to a vector a by stacking its columns on top of each other, while unvecI×Ra converts a∈ℂIR to a matrix A∈ℂI×R. vecdA converts the diagonal elements of A into a vector. ∘ denotes the outer product operator. The Kronecker and Khatri-Rao products are denoted by ⊗ and ⋄, respectively. The Khatri-Rao product between the matrices A=a1,…,aR∈ℂI×R and B=b1,…,bR∈ℂJ×R corresponds to a column-wise Kronecker product, i.e.,(1)A⋄B=a1⊗b1,a2⊗b2,…,aR⊗bR∈ℂIJ×R.We shall make use of the following two properties of the Kronecker and Khatri-Rao products:(2)AC⋄BD=A⊗BC⋄D,(3)vecABCT=C⋄AvecdB,(4)a⊗b=vecb∘a∈ℂIJ,where B is assumed to be a diagonal matrix in (3). In both cases, the matrices have compatible dimensions.


## 2. Tensor Preliminaries

In order to facilitate the presentation of the proposed receivers, we provide below a brief overview on some important tensor definitions and tensor algebra operations. We also introduce the PARAFAC decomposition and its different representation forms.


### 2.1. Tensor Definitions and Basic Operations

Throughout this paper, the definitions and operations involving tensors are in accordance with [28, 29]. A tensor is defined here as a multidimensional array. The order of a tensor corresponds to the number of dimensions. It can be seen as a generalization of a matrix to higher-order dimensions. For instance, a scalar is a tensor of order 0, a vector is a tensor of order 1, and a matrix is a tensor of order 2. An n-mode fiber of a tensor is a vector obtained by varying the n-th index and keeping all the other indexes fixed. Slices are two-dimensional sections of a tensor, obtained by fixing all but two indices. The operator A⊔nB denotes the concatenation of two matrices along the n-th mode of a tensor. The 1 mode, 2 mode, and 3 mode unfolding matrices of the third-order tensor X∈ℂI1×I2×I3, denoted by X1∈ℂI1×I2I3, X2∈ℂI2×I1I3, and X3∈ℂI3×I1I2, are obtained by collecting all the 1 mode, 2 mode, and 3 mode fibers to be columns of the resulting matrices, respectively. The n-mode product between the tensor X and a matrix A∈ℂRn×In is denoted by Y=X×An, which is equivalent in a matrix fashion to Yn=AXn n=1,2,3.


### 2.2. PARAFAC Decomposition

By definition, the PARAlell FACtor (PARAFAC) analysis decomposition of a third-order tensorX∈ℂI1×I2×I3, introduced by [24], is the factorization of X in a sum of R third-order rank-one tensors each one being formed by the outer product of three vectors. Mathematically, the PARAFAC decomposition of X∈ℂI1×I2×I3 is given by(5)X=∑r=1Rar1∘ar2∘ar3,where R is the rank of the PARAFAC decomposition and is defined as the minimum number of rank-one tensors for which X holds exactly. The vector arn∈ℂIn denotes the r-th column of the factor matrix An=a1n,…,aRn∈ℂIn×R along the n-th mode n=1,2,3.The PARAFAC decomposition can also be represented in terms of the frontal slices ofX as follows:(6)Xi3=A1Di3A3A2T∈ℂI1×I2,for i3=1,…,I3.By usingn-mode product notation, equations (5) and (6) can be written as(7)X=ℐ3,R×1A1×2A2×3A3,where ℐ3,R denotes a third-order identity tensor of size R×R×R. Its elements are equal to 1 when all indices are equal and 0 elsewhere.The 1 mode, 2 mode, and 3 mode unfolding matrices ofX admit the following factorizations with respect to the factor matrices An n=1,2,3:(8)X1=A1A3⋄A2T,(9)X2=A2A3⋄A1T,(10)X3=A3A2⋄A1T.


## 2.1. Tensor Definitions and Basic Operations

Throughout this paper, the definitions and operations involving tensors are in accordance with [28, 29]. A tensor is defined here as a multidimensional array. The order of a tensor corresponds to the number of dimensions. It can be seen as a generalization of a matrix to higher-order dimensions. For instance, a scalar is a tensor of order 0, a vector is a tensor of order 1, and a matrix is a tensor of order 2. An n-mode fiber of a tensor is a vector obtained by varying the n-th index and keeping all the other indexes fixed. Slices are two-dimensional sections of a tensor, obtained by fixing all but two indices. The operator A⊔nB denotes the concatenation of two matrices along the n-th mode of a tensor. The 1 mode, 2 mode, and 3 mode unfolding matrices of the third-order tensor X∈ℂI1×I2×I3, denoted by X1∈ℂI1×I2I3, X2∈ℂI2×I1I3, and X3∈ℂI3×I1I2, are obtained by collecting all the 1 mode, 2 mode, and 3 mode fibers to be columns of the resulting matrices, respectively. The n-mode product between the tensor X and a matrix A∈ℂRn×In is denoted by Y=X×An, which is equivalent in a matrix fashion to Yn=AXn n=1,2,3.


## 2.2. PARAFAC Decomposition

By definition, the PARAlell FACtor (PARAFAC) analysis decomposition of a third-order tensorX∈ℂI1×I2×I3, introduced by [24], is the factorization of X in a sum of R third-order rank-one tensors each one being formed by the outer product of three vectors. Mathematically, the PARAFAC decomposition of X∈ℂI1×I2×I3 is given by(5)X=∑r=1Rar1∘ar2∘ar3,where R is the rank of the PARAFAC decomposition and is defined as the minimum number of rank-one tensors for which X holds exactly. The vector arn∈ℂIn denotes the r-th column of the factor matrix An=a1n,…,aRn∈ℂIn×R along the n-th mode n=1,2,3.The PARAFAC decomposition can also be represented in terms of the frontal slices ofX as follows:(6)Xi3=A1Di3A3A2T∈ℂI1×I2,for i3=1,…,I3.By usingn-mode product notation, equations (5) and (6) can be written as(7)X=ℐ3,R×1A1×2A2×3A3,where ℐ3,R denotes a third-order identity tensor of size R×R×R. Its elements are equal to 1 when all indices are equal and 0 elsewhere.The 1 mode, 2 mode, and 3 mode unfolding matrices ofX admit the following factorizations with respect to the factor matrices An n=1,2,3:(8)X1=A1A3⋄A2T,(9)X2=A2A3⋄A1T,(10)X3=A3A2⋄A1T.


## 3. System and Channel Models

In this section, we introduce the proposed closed-loop and multifrequency channel training framework. Then, we formulate our DL and UL signal models. In addition, the considered mmWave massivo MIMO channel model is also presented.


### 3.1. Downlink Signal Model

Consider a wireless communication system operating in the FDD mode, where a BS equipped withNBS antennas serves simultaneously U MSs equipped with NMS antennas. We assume that the BS employs a hybrid beamforming architecture using MRF RF chains. Due to the different instants of time dedicated to DL and UL communications, the beamforming matrix associated with DL transmission is denoted by W=WRFWBB∈ℂNBS×MRF, while the beamforming matrix associated with UL reception is denoted by F=FRFFBB∈ℂNBS×MRF. Note that equal beamforming matrices can also be considered in the transmission and reception phases without loss of generality. The BS transmits a length-T pilot sequence sp∈ℂT over the p-th spatial direction using the beamforming vector wp∈ℂNBS (p=1,…,P and P≤MRF). The received signal at the u-th MS over P different directions is given by(11)Yu=HuWS+VuDL∈ℂNMS×T,where Hu∈ℂNMS×NBS denotes the DL channel matrix associated with the u-th MS, W=w1,w2,…,wP∈ℂNBS×P denotes the transmission beamforming matrix, S=s1,s2,…,sPT∈ℂP×T concatenates the pilot sequences to be sent by each transmission beam. The matrix VuDL∈ℂNMS×T is the additive white Gaussian noise (AWGN) term at the u-th MS.During the training phase, we assume identity matrices for the digital beamforming matrices, while the analog beamforming matrices have constant unit modulus entries with random phases. Thus, the entries ofW and F are chosen uniformly from a unit circle scaled by a constant 1/NBS, i.e.,(12)Wi,j=1NBSejϑi,j,Fi,j=1NBSejφi,j,where ϑi,j and φi,j∈−π,π follow a uniform distribution. Since this work deals with the channel estimation problem, the optimum design of the beamforming matrices is not addressed here.


### 3.2. Uplink Signal Model

The pilot signal (11) received at the u-th MS is fed back to the BS after a multifrequency coding operation (i.e., no channel estimation is done at the MS side). More specifically, we assume that Yu (u=1,…,U) is coded and spread in the frequency domain across K adjacent subcarriers over which the fading channel is assumed to be constant. The coded signal of the u-th MS transmitted at the k-th subcarrier can be expressed as(13)Yk,u=diagck,uYu∈ℂNMS×T,where ck,u∈ℂNMS denotes a known code vector associated with the k-th subcarrier and used by u-th MS. It is worth noting that the coding vectors used by the different MSs do not need to be orthogonal. As will be discussed later, the linear independence assumption is enough. In practice, this means that these codes can be locally generated at each MS as pseudorandom sequences, i.e., no prior signaling and coordination between MSs is necessary [30].In the UL communication, the BS employsQ beamforming vectors fq∈ℂNBS (q=1,…,Q and Q≤MRF) to receive the coded uplink pilot signals over a set of Q different spatial directions. The received closed-loop signal at the BS associated with the k-th subcarrier is then given by(14)Xk=FH∑u=1UGuYk,u+VUL=FHGediagckYeT+FHVUL∈ℂQ×T,where Ge=G1,G2,…,GU∈ℂNBS×UNMS denotes an extended version of the UL channel matrix that concatenates the U UL channel matrices Gu∈ℂNBS×NMS (u=1,…,U) of all MSs, Ye=Y1T,Y2T,…,YUT∈ℂT×UNMS denotes an extended matrix that concatenates the feedback signals sent by all MSs, ck=ck,1T,ck,2T,…,ck,UTT∈ℂUNMS is an extended code vector that contains the coding vectors of all MSs with respect to the k-th subcarrier, and FHVUL represents the filtered noise term at the output of the RF chains.


### 3.3. Conventional× Proposed Channel Training Framework

The conventional framework for channel estimation, summarized in Figure1, assumes channel reciprocity in TDD or treat the DL and UL channel estimation as two separated problems in FDD, i.e., solved independently at the MS and BS, respectively. For the DL channel estimation, the BS first sends pilot signals to all MSs. At the MS side, the DL channel estimation can be performed by means of the state-of-the-art least squares (LS), minimum mean square error (MMSE), or CS-based estimators. Then, the estimated DL channel is reported back to the BS via dedicated UL resources [30]. To solve the UL channel estimation problem, a pilot signal is sent to the BS by each MS. Finally, the UL channels of all MSs are estimated by the BS. In practice, the DL channel estimated from UL pilots under the reciprocity assumption may not be accurate due to radio frequency distortions or a different carrier frequency such as in FDD. In addition, the conventional framework for channel estimation may imply a high computational complexity at the MS side, especially for power-limited devices.Figure 1
Conventional training framework. The DL and UL channel estimation problems are solved independently. The BS first transmits pilot sequences. Then, the DL channel is estimated at the MS side. The estimated DL channel is fed back to the BS via dedicated uplink resources. The UL channel is estimated at the BS side. The text in blue refer to DL communication, while the text in red refer to UL communication.In the proposed framework, summarized in Figure2, no processing for channel estimation is performed at the MSs side. In contrast to the conventional approach, the received pilot signals at each MS are fed back to the BS after a multifrequency coding operation across a set of adjacent subcarriers. After this closed-loop procedure, the joint DL and UL channel estimation can be performed at the BS from the received signal given in (14). We can note that the proposed framework alleviates computational overhead due to channel estimation at the power-limited devices, by shifting this processing burden to the BS side. Furthermore, it also relaxes channel reciprocity assumptions since the DL and UL channels can be estimated jointly from (14).Figure 2
Proposed closed-loop and multifrequency-based training framework. The DL and UL channels are jointly estimated. The BS first transmits pilot sequences. The MSs encode the received pilots and then feed them back to the BS. The BS jointly estimates the DL and UL channels. The text in blue refer to DL communication, while the text in red refer to UL communication.After estimating the UL and DL channels, the BS may report the estimated parameters (AoDs, AoAs, and path gains) to the MS. Then, each MS can rebuild an estimation of the DL channel (according to procedure presented in Section4.5) before decoding the information data. This approach reduces the overhead since a number of DL channel parameters reported to the MS is much smaller with respect to the size of the DL channel matrix in an mmWave MIMO system.


### 3.4. Channel Model

In (11) and (14), we consider a general formulation in which the DL and UL channels are completely independent. In other words, the channels do not share any reciprocity in the angular or path gain domains. We also assume that the UL channels are constant across the K adjacent subcarriers used in the multifrequency coding operation. Due to the severe path loss, mmWave channels can be modeled by a narrow-band clustered channel model with few Lu dominant paths between the u-th MS and the BS. The DL channel matrix Hu∈ℂNMS×NBS associated with the u-th MS can be written as [19](15)Hu=∑l=1Luαu,laMSθu,laBSϕu,lT,where αu,l denotes the complex path gain of the u-th MS related to the l-th path in the DL communication. The path gains are modeled as circular symmetric Gaussian random variables with zero mean and unit variance. aMSθu,l∈ℂNMS and aBSϕu,l∈ℂNBS are the antenna array response vectors evaluated at the angle of arrival θu,l and angle of departure ϕu,l uniformly distributed in the interval 0,2π. Throughout this paper, we assume uniform linear arrays (ULAs) at the BS and MSs. However, the proposed method can be applied to arbitrary array geometries without loss of generality. For ULA configurations with interantennas spacing equals to d=λ/2, where λ denotes the wavelength of the signal, the array response vectors at the MS and BS can be formulated as(16)aMSθu,l=1NMS1,ejπcosθu,l,…,ejπNMS−1cosθu,lT,aBSϕu,l=1NBS1,ejπcosϕu,l,…,ejπNBS−1cosϕu,lT.In matrix form,Hu can be rewritten as(17)Hu=AMSdiagαABST,where α=NMSNBS/Luαu,1,αu,2,…,αu,LuT∈ℂLu denotes the vector that contains the Lu path gains in the DL. The array response matrices AMS∈ℂNMS×Lu and ABS∈ℂNBS×Lu at the MS and BS are expressed as(18)AMS=aMSθu,1,aMSθu,2,…,aMSθu,Lu,ABS=aBSϕu,1,aBSϕu,2,…,aBSϕu,Lu.The UL channel matrixGu∈ℂNBS×NMS from the u-th MS to the BS can be represented in a similar way. We define Gu as follows(19)Gu=ABSdiagβAMST,where ABS∈ℂNBS×Mu and AMS∈ℂNMS×Mu are now functions of the spatial parameters in the UL, while β=NMSNBS/Muβu,1,βu,2,…,βu,MuT∈ℂMu denotes the vector that contains the Mu path gains of the UL channel.


## 3.1. Downlink Signal Model

Consider a wireless communication system operating in the FDD mode, where a BS equipped withNBS antennas serves simultaneously U MSs equipped with NMS antennas. We assume that the BS employs a hybrid beamforming architecture using MRF RF chains. Due to the different instants of time dedicated to DL and UL communications, the beamforming matrix associated with DL transmission is denoted by W=WRFWBB∈ℂNBS×MRF, while the beamforming matrix associated with UL reception is denoted by F=FRFFBB∈ℂNBS×MRF. Note that equal beamforming matrices can also be considered in the transmission and reception phases without loss of generality. The BS transmits a length-T pilot sequence sp∈ℂT over the p-th spatial direction using the beamforming vector wp∈ℂNBS (p=1,…,P and P≤MRF). The received signal at the u-th MS over P different directions is given by(11)Yu=HuWS+VuDL∈ℂNMS×T,where Hu∈ℂNMS×NBS denotes the DL channel matrix associated with the u-th MS, W=w1,w2,…,wP∈ℂNBS×P denotes the transmission beamforming matrix, S=s1,s2,…,sPT∈ℂP×T concatenates the pilot sequences to be sent by each transmission beam. The matrix VuDL∈ℂNMS×T is the additive white Gaussian noise (AWGN) term at the u-th MS.During the training phase, we assume identity matrices for the digital beamforming matrices, while the analog beamforming matrices have constant unit modulus entries with random phases. Thus, the entries ofW and F are chosen uniformly from a unit circle scaled by a constant 1/NBS, i.e.,(12)Wi,j=1NBSejϑi,j,Fi,j=1NBSejφi,j,where ϑi,j and φi,j∈−π,π follow a uniform distribution. Since this work deals with the channel estimation problem, the optimum design of the beamforming matrices is not addressed here.


## 3.2. Uplink Signal Model

The pilot signal (11) received at the u-th MS is fed back to the BS after a multifrequency coding operation (i.e., no channel estimation is done at the MS side). More specifically, we assume that Yu (u=1,…,U) is coded and spread in the frequency domain across K adjacent subcarriers over which the fading channel is assumed to be constant. The coded signal of the u-th MS transmitted at the k-th subcarrier can be expressed as(13)Yk,u=diagck,uYu∈ℂNMS×T,where ck,u∈ℂNMS denotes a known code vector associated with the k-th subcarrier and used by u-th MS. It is worth noting that the coding vectors used by the different MSs do not need to be orthogonal. As will be discussed later, the linear independence assumption is enough. In practice, this means that these codes can be locally generated at each MS as pseudorandom sequences, i.e., no prior signaling and coordination between MSs is necessary [30].In the UL communication, the BS employsQ beamforming vectors fq∈ℂNBS (q=1,…,Q and Q≤MRF) to receive the coded uplink pilot signals over a set of Q different spatial directions. The received closed-loop signal at the BS associated with the k-th subcarrier is then given by(14)Xk=FH∑u=1UGuYk,u+VUL=FHGediagckYeT+FHVUL∈ℂQ×T,where Ge=G1,G2,…,GU∈ℂNBS×UNMS denotes an extended version of the UL channel matrix that concatenates the U UL channel matrices Gu∈ℂNBS×NMS (u=1,…,U) of all MSs, Ye=Y1T,Y2T,…,YUT∈ℂT×UNMS denotes an extended matrix that concatenates the feedback signals sent by all MSs, ck=ck,1T,ck,2T,…,ck,UTT∈ℂUNMS is an extended code vector that contains the coding vectors of all MSs with respect to the k-th subcarrier, and FHVUL represents the filtered noise term at the output of the RF chains.


## 3.3. Conventional× Proposed Channel Training Framework

The conventional framework for channel estimation, summarized in Figure1, assumes channel reciprocity in TDD or treat the DL and UL channel estimation as two separated problems in FDD, i.e., solved independently at the MS and BS, respectively. For the DL channel estimation, the BS first sends pilot signals to all MSs. At the MS side, the DL channel estimation can be performed by means of the state-of-the-art least squares (LS), minimum mean square error (MMSE), or CS-based estimators. Then, the estimated DL channel is reported back to the BS via dedicated UL resources [30]. To solve the UL channel estimation problem, a pilot signal is sent to the BS by each MS. Finally, the UL channels of all MSs are estimated by the BS. In practice, the DL channel estimated from UL pilots under the reciprocity assumption may not be accurate due to radio frequency distortions or a different carrier frequency such as in FDD. In addition, the conventional framework for channel estimation may imply a high computational complexity at the MS side, especially for power-limited devices.Figure 1
Conventional training framework. The DL and UL channel estimation problems are solved independently. The BS first transmits pilot sequences. Then, the DL channel is estimated at the MS side. The estimated DL channel is fed back to the BS via dedicated uplink resources. The UL channel is estimated at the BS side. The text in blue refer to DL communication, while the text in red refer to UL communication.In the proposed framework, summarized in Figure2, no processing for channel estimation is performed at the MSs side. In contrast to the conventional approach, the received pilot signals at each MS are fed back to the BS after a multifrequency coding operation across a set of adjacent subcarriers. After this closed-loop procedure, the joint DL and UL channel estimation can be performed at the BS from the received signal given in (14). We can note that the proposed framework alleviates computational overhead due to channel estimation at the power-limited devices, by shifting this processing burden to the BS side. Furthermore, it also relaxes channel reciprocity assumptions since the DL and UL channels can be estimated jointly from (14).Figure 2
Proposed closed-loop and multifrequency-based training framework. The DL and UL channels are jointly estimated. The BS first transmits pilot sequences. The MSs encode the received pilots and then feed them back to the BS. The BS jointly estimates the DL and UL channels. The text in blue refer to DL communication, while the text in red refer to UL communication.After estimating the UL and DL channels, the BS may report the estimated parameters (AoDs, AoAs, and path gains) to the MS. Then, each MS can rebuild an estimation of the DL channel (according to procedure presented in Section4.5) before decoding the information data. This approach reduces the overhead since a number of DL channel parameters reported to the MS is much smaller with respect to the size of the DL channel matrix in an mmWave MIMO system.


## 3.4. Channel Model

In (11) and (14), we consider a general formulation in which the DL and UL channels are completely independent. In other words, the channels do not share any reciprocity in the angular or path gain domains. We also assume that the UL channels are constant across the K adjacent subcarriers used in the multifrequency coding operation. Due to the severe path loss, mmWave channels can be modeled by a narrow-band clustered channel model with few Lu dominant paths between the u-th MS and the BS. The DL channel matrix Hu∈ℂNMS×NBS associated with the u-th MS can be written as [19](15)Hu=∑l=1Luαu,laMSθu,laBSϕu,lT,where αu,l denotes the complex path gain of the u-th MS related to the l-th path in the DL communication. The path gains are modeled as circular symmetric Gaussian random variables with zero mean and unit variance. aMSθu,l∈ℂNMS and aBSϕu,l∈ℂNBS are the antenna array response vectors evaluated at the angle of arrival θu,l and angle of departure ϕu,l uniformly distributed in the interval 0,2π. Throughout this paper, we assume uniform linear arrays (ULAs) at the BS and MSs. However, the proposed method can be applied to arbitrary array geometries without loss of generality. For ULA configurations with interantennas spacing equals to d=λ/2, where λ denotes the wavelength of the signal, the array response vectors at the MS and BS can be formulated as(16)aMSθu,l=1NMS1,ejπcosθu,l,…,ejπNMS−1cosθu,lT,aBSϕu,l=1NBS1,ejπcosϕu,l,…,ejπNBS−1cosϕu,lT.In matrix form,Hu can be rewritten as(17)Hu=AMSdiagαABST,where α=NMSNBS/Luαu,1,αu,2,…,αu,LuT∈ℂLu denotes the vector that contains the Lu path gains in the DL. The array response matrices AMS∈ℂNMS×Lu and ABS∈ℂNBS×Lu at the MS and BS are expressed as(18)AMS=aMSθu,1,aMSθu,2,…,aMSθu,Lu,ABS=aBSϕu,1,aBSϕu,2,…,aBSϕu,Lu.The UL channel matrixGu∈ℂNBS×NMS from the u-th MS to the BS can be represented in a similar way. We define Gu as follows(19)Gu=ABSdiagβAMST,where ABS∈ℂNBS×Mu and AMS∈ℂNMS×Mu are now functions of the spatial parameters in the UL, while β=NMSNBS/Muβu,1,βu,2,…,βu,MuT∈ℂMu denotes the vector that contains the Mu path gains of the UL channel.


## 4. Proposed Tensor-Based Semiblind Receivers for Joint DL and UL Channel Estimation

Our aim is to jointly estimate the DL and UL channel matricesHu and Gu of each MS (u=1,…,U) by solving a multiuser channel estimation problem at the BS. To this end, we first rewrite the received closed-loop signal (14) using a tensor formalism. Then, by capitalizing on its multidimensional structure, we formulate two tensor-based semiblind receivers to initially obtain the estimates for the compressed channel matrices. As a final step, we exploit the sparse representation of the DL and UL channels to estimate their parameters by decoupling the multiuser channel estimation problem into multiple single-user ones solved in a parallel way via separate CS problems.


### 4.1. PARAFAC Modeling

According to (6) and (7), the noiseless term in the received closed-loop signal (14) can be interpreted as the k-th frontal slice of the following third-order PARAFAC decomposition:(20)X=ℐ3,U⋅NMS×1FHGe×2Ye×3C∈ℂQ×T×K,obtained by concatenating the K signal matrices Xkk=1K associated with the different adjacent subcarriers along the third mode of X, i.e.,(21)X=X1⊔3X2⊔3…⊔3XK.By analogy with (7), the following correspondence holds:(22)A1,A2,A3⟷FHGe,Ye,C,I1,I2,I3,R⟷Q,T,K,UNMS.The three dimensions, or modes, ofX stands for the number of receive beams, pilot sequence length, and number of subcarriers. The matrix F can be seen as a compression matrix associated with the first mode of X which reduces the size of the first mode from NBS to Q RF chains (i.e., number of beams). The k-th row of the multifrequency coding matrix C∈ℂK×UNMS contains the code values used by the U MSs at the k-th subcarrier, i.e.,(23)C=c1,c2,…,cKT∈ℂK×UNMS.Estimates of the channel parameters (AoDs, AoAs, and path gains) that build up the channel matricesHuu=1U and Guu=1U can be obtained by fitting the noisy version of X to a PARAFAC decomposition. In the following, we formulate the first stage of the proposed receivers that consists of estimating the factor matrices FHGe∈ℂQ×UNMS and Ye∈ℂT×UNMS in an iterative or closed-form way. Once the factor matrices are estimated, the second stage of the proposed receivers is to solve U independent CS problems that yield to channel parameter estimation of each MS, as will be shown later.


### 4.2. First Stage: Bilinear Alternating Least Squares (B-ALS Receiver)

According to (8)–(10), we can obtain the following representations for the unfolding matrices X1∈ℂQ×TK, X2∈ℂT×QK, and X3∈ℂK×QT of X in terms of its factor matrices:(24)X1=ΦC⋄YeT,(25)X2=YeC⋄ΦT,(26)X3=CYe⋄ΦT,where we denote Φ=FHGe for the simplicity of representation.Since the multifrequency coding matrixC is assumed to be known at the BS while the DL and UL channel-state information are not available, the proposed B-ALS receiver consists of estimating Φ and Ye in an alternating way from X1 and X2 by optimizing, respectively, the following two nonlinear LS problems:(27)Φ^=argminΦX1−ΦC⋄YeTF2,Y^e=argminYeX2−YeC⋄ΦTF2.The solutions of which are given byΦ^=X1C⋄Y^eT† and Y^e=X2C⋄Φ^T†, respectively.Each iteration of the bilinear ALS-PARAFAC algorithm contains only two LS updating steps. At each step, one factor matrix is updated, while the other is assumed fixed to its value obtained in the previous step [31]. This procedure is repeated until the convergence of the algorithm, denoted by(28)εi=X1−Φ^C⋄Y^eTF2,the residual error between the received signal tensor and the reconstructed signal tensor at the i-th iteration. We declare that the first stage has converged at the i-th iteration when(29)εi−εi−1≤σ,where σ is a threshold. In our computational simulations, we set σ=10−6. Convergence to the global minimum is always achieved within a few iterations due to the knowledge of the frequency spread matrix at the BS. The proposed B-ALS receiver is summarized in Algorithm 1.Algorithm 1:Pseudocode of the B-ALS receiver.
(1)
Seti=0;
Initialize randomly the factor matrixY^ei=0;
(2)
i
⟵
i
+
1;
(3)
FromX1, obtain a LS estimate of Φ^i:
Φ^i=X1C⋄Y^ei−1T†;
(4)
FromX2, obtain a LS estimate of Y^ei:
Y^ei=X2C⋄Φ^iT†;
(5)
Repeat steps 2–4 until convergence. The convergence is achieved whenεi−εi−1≤10−6, where εi is the residual error computed in the i-th iteration.


### 4.3. Alternative Closed-Form Solution to the First Stage

In contrast to the B-ALS receiver, the second proposed receiver named LS-KRF is an alternative closed-form solution that can be employed in particular cases in whichK≥UNMS. The idea is to filter the received signal tensor X by exploiting the knowledge of the multifrequency coding matrix C and then solve a set of rank 1 approximation problems.Initially, by multiplying both sides ofX3T in (26) by the pseudoinverse of CT from the right-hand side, we obtain(30)Ye⋄Φ=X3TCT†=y1⊗φ1,…,yUNMS⊗φUNMS∈ℂQT×UNMS.According to property in (4), the i-th column of (30) i=1,…,UNMS can be rewritten as(31)yi⊗φi=vecφi∘yi,which denotes the vectorization operation of the rank 1 matrix Ψi=φi∘yi∈ℂQ×T. By defining UiΣiViH as the singular value decomposition (SVD) of Ψi, estimates for φi∈ℂQ and yi∈ℂT i=1,…,UNMS can be obtained by truncating the SVD of Ψi to a rank 1 approximation, i.e., [27](32)φ^i=σ1u1,y^i=σ1v1∗,where u1∈ℂQ and v1∈ℂT are the corresponding first left and right singular vectors of Ui and Vi, respectively. σ1 denotes the largest singular value of the matrix Σi. Final estimates for Φ^ and Y^e are obtained by repeating this SVD computation i=1,…,UNMS times in parallel, one for each column of (30). The pseudocode of the LS-KRF receiver is summarized in Algorithm 2.Algorithm 2:Pseudocode of the LS-KRF receiver.
(1)
Apply theunvecQ×T operator in the i-th column of (30) to obtain the rank 1 matrix Ψi∈ℂQ×T;
(2)
Compute the SVDΨi=UiΣiViH. Then, obtain the estimates for the i-th columns of Φ^ and Y^e as:
φ^i=σ1u1andy^i=σ1v1∗,
whereu1∈ℂQ and v1∈ℂT are the first left and right singular vectors of Ui and Vi, respectively, and σ1 is the largest singular value;
(3)
Repeat steps 1-2 for all columns of (30).


### 4.4. Second Stage: Sparse Formulation to DL and UL Channel Parameters Estimation

Once the matricesΦ^ and Y^e are estimated, the second stage of the proposed receivers consists in estimating the channel parameters (AoDs, AoAs, and path gains) to reconstruct the channel matrices H^u and G^u related to each MS. Thanks to the knowledge of the multifrequency coding matrix C, the estimated factor matrices are not affected by permutation of columns ambiguity. Therefore, the first stage provides automatic separation of the compressed users’ channels. From Φ^ and Y^e, the multiuser channel estimation problem can be decoupled into U single-user ones as formulated below.Let us rewrite the block representation ofΦ^, defined as(33)Φ^=Φ^1,Φ^2,…,Φ^U∈ℂQ×UNMS,where(34)Φ^u=FHG^u∈ℂQ×NMS,u=1,…,U.where G^u∈ℂNBS×NMS denotes the estimated UL channel matrix related to the u-th MS. By replacing G^u for (19), and then vectorizing (34) according to the property in (3), we have(35)φ^u=vecΦ^u=A^MS⋄FHA^BSβ^.Using the property in (2), we straightforwardly obtain(36)φ^u=INMS⊗FHA^MS⋄A^BSβ^∈ℂQNMS,where INMS denotes an identity matrix of size NMS×NMS.The same procedure can directly be applied in theu-th block Y^uT of the estimated factor matrix:(37)Y^e=Y^1T,Y^2T,…,Y^UT∈ℂT×UNMS,where(38)Y^uT=STWTH^uT+VuDLT∈ℂT×NMS.H
^
u
∈
ℂ
N
MS
×
N
BS denotes the estimated DL channel matrix defined above in (17). Similar to (36), we obtain the following vector formulation:(39)y^u=INMS⊗STWTA^MS⋄A^BSα^+vuDL∈ℂTNMS,where y^u=vecY^uT and vuDL=vecVuDLT.From (36) and (39), two independent CS problems can be formulated to jointly estimate the parameters of the DL and UL channels of the u-th MS. We assume that grid quantization errors are neglected, i.e., the AoDs and AoAs are drawn from a uniform angle grid of N points contained in the set 0,2π/N,…,2πN−1/N, with N≫Lu and N≫Mu. Based on this assumption, we can obtain the following sparse formulations for the vectors y^u and φ^u, respectively:(40)y^u=INMS⊗STWTΣDα¯,(41)φ^u=INMS⊗FHΣDβ¯,where ΣD∈ℂNMSNBS×N2 denotes the known dictionary matrix used to solve the sparse signal recovery problem, defined as(42)ΣD=A¯MS⊗A¯BS∈ℂNMSNBS×N2.The matricesA¯MS∈ℂNMS×N and A¯BS∈ℂNBS×N that make the dictionary are given by(43)A¯MS=aMS0,aMS2πN,…,aMS2πN−1N,A¯BS=aBS0,aBS2πN,…,aBS2πN−1N,and contain all points of the uniform angle grid. The left-hand side matrices INMS⊗STWT and INMS⊗FH are called measurement matrices of the sparse problems. α¯∈ℂN2 and β¯∈ℂN2 are sparse vectors obtained by augmenting the vector gains α and β with zero elements, respectively.Estimates for the parameters of the channel matricesHu and Gu can be obtained by applying CS algorithms in the estimated sparse vectors (40) and (41). Many efficient algorithms such as orthogonal matching pursuit (OMP) [22], structured compressive sampling matching pursuit (S-CoSaMP) [25], and fast iterative shrinkage-thresholding (FISTA) [23], to name a few, can be used to solve these two sparse signal recovery problems. In a simplified view, the estimates for the path gains α^ and β^ correspond to nonzero entries of the estimated sparse vectors α¯ and β¯, while estimation for the spatial parameters (AoDs and AoAs) are obtained by selecting the columns of the dictionary matrix ΣD related to the positions of the estimated path gains in the sparse vector. In our numerical results, we adopt the OMP algorithm to estimate the UL and DL channel parameters for simplicity reasons, although any state-of-the-art CS-based algorithm is equally applicable to solve problems (40) and (41).Remark:.
Compared to (36), the sparse signal recovery problem formulated from (39) naturally incorporates the DL noise contribution in its structure. For this reason, the proposed closed-loop framework for channel estimation can lead to some performance degradation in the DL channel estimation compared to UL channel estimation. Therefore, we can observe a trade-off between DL channel estimation accuracy performed by the BS and reduction of the processing cost for channel estimation at the MS side. This discussion is reinforced by means of numerical simulations in Section 6.


### 4.5. Joint DL and UL Channel Estimation

Finally, from the estimated channel parameters (AoDs, AoAs, and path gains), the BS can construct the estimated DL and UL channel matricesH^u and G^u of the u-th MS according to relations (17) and (19) as follows:(44)H^u=A^MSdiagα^A^BST,(45)G^u=A^BSdiagβ^A^MST.As previously presented in Section4.4, the proposed receivers decouple the multiuser channel estimation problem into 2U single-user ones (U problems dedicated to each communication link) that can be solved independently for each MS. Since 2U digital processing units are available at the BS, the second stage of the proposed receivers can be computed in parallel. Therefore, its processing delay can be kept constant (i.e., it does not increase with the number of MSs), when the BS is equipped with multiple (at least 2U) digital processing units. The parallelized processing for the second stage of the proposed receivers is illustrated in Figures 3 and 4. In addition, the overall pseudocode of the proposed two-stage tensor-based receivers for joint DL and UL channel estimation is summarized in Algorithm 3.Figure 3
Illustration of the parallelized processing for the estimation of U uplink channels. Theu-th block A^u of the estimated factor matrix A^ is forwarded and processed for a dedicated processor.Figure 4
Illustration of the parallelized processing for the estimation of the U downlink channels. Theu-th block Y^uT of the estimated factor matrix B^ is forwarded and processed for a dedicated processor.Algorithm 3:Overall pseudocode of the proposed receivers for joint DL and UL channel estimation.
1
First Stage. estimation of the compressed channel matrices
(1.1)
From the received signal tensorX in (20), obtain the estimated factor matrices Φ^ and Y^e via B-ALS or LS-KRF described in Algorithms 1 and 2, respectively,
(2)
Second Stage. parameters estimation and channels reconstruction
(2.1)
FromΦ^=Φ^1,Φ^2,…,Φ^UandY^e=Y^1T,Y^2T,…,Y^UT, obtain the estimates of the channel parameters (AoDs, AoAs, and path gains) for each MS by applying a CS recovery algorithm (e.g., OMP) to problems (40) and (41) independently
(2.2)
The BS constructs the estimated DL and UL channel matricesH^u and G^u according to (44) and (45), respectively.


## 4.1. PARAFAC Modeling

According to (6) and (7), the noiseless term in the received closed-loop signal (14) can be interpreted as the k-th frontal slice of the following third-order PARAFAC decomposition:(20)X=ℐ3,U⋅NMS×1FHGe×2Ye×3C∈ℂQ×T×K,obtained by concatenating the K signal matrices Xkk=1K associated with the different adjacent subcarriers along the third mode of X, i.e.,(21)X=X1⊔3X2⊔3…⊔3XK.By analogy with (7), the following correspondence holds:(22)A1,A2,A3⟷FHGe,Ye,C,I1,I2,I3,R⟷Q,T,K,UNMS.The three dimensions, or modes, ofX stands for the number of receive beams, pilot sequence length, and number of subcarriers. The matrix F can be seen as a compression matrix associated with the first mode of X which reduces the size of the first mode from NBS to Q RF chains (i.e., number of beams). The k-th row of the multifrequency coding matrix C∈ℂK×UNMS contains the code values used by the U MSs at the k-th subcarrier, i.e.,(23)C=c1,c2,…,cKT∈ℂK×UNMS.Estimates of the channel parameters (AoDs, AoAs, and path gains) that build up the channel matricesHuu=1U and Guu=1U can be obtained by fitting the noisy version of X to a PARAFAC decomposition. In the following, we formulate the first stage of the proposed receivers that consists of estimating the factor matrices FHGe∈ℂQ×UNMS and Ye∈ℂT×UNMS in an iterative or closed-form way. Once the factor matrices are estimated, the second stage of the proposed receivers is to solve U independent CS problems that yield to channel parameter estimation of each MS, as will be shown later.


## 4.2. First Stage: Bilinear Alternating Least Squares (B-ALS Receiver)

According to (8)–(10), we can obtain the following representations for the unfolding matrices X1∈ℂQ×TK, X2∈ℂT×QK, and X3∈ℂK×QT of X in terms of its factor matrices:(24)X1=ΦC⋄YeT,(25)X2=YeC⋄ΦT,(26)X3=CYe⋄ΦT,where we denote Φ=FHGe for the simplicity of representation.Since the multifrequency coding matrixC is assumed to be known at the BS while the DL and UL channel-state information are not available, the proposed B-ALS receiver consists of estimating Φ and Ye in an alternating way from X1 and X2 by optimizing, respectively, the following two nonlinear LS problems:(27)Φ^=argminΦX1−ΦC⋄YeTF2,Y^e=argminYeX2−YeC⋄ΦTF2.The solutions of which are given byΦ^=X1C⋄Y^eT† and Y^e=X2C⋄Φ^T†, respectively.Each iteration of the bilinear ALS-PARAFAC algorithm contains only two LS updating steps. At each step, one factor matrix is updated, while the other is assumed fixed to its value obtained in the previous step [31]. This procedure is repeated until the convergence of the algorithm, denoted by(28)εi=X1−Φ^C⋄Y^eTF2,the residual error between the received signal tensor and the reconstructed signal tensor at the i-th iteration. We declare that the first stage has converged at the i-th iteration when(29)εi−εi−1≤σ,where σ is a threshold. In our computational simulations, we set σ=10−6. Convergence to the global minimum is always achieved within a few iterations due to the knowledge of the frequency spread matrix at the BS. The proposed B-ALS receiver is summarized in Algorithm 1.Algorithm 1:Pseudocode of the B-ALS receiver.
(1)
Seti=0;
Initialize randomly the factor matrixY^ei=0;
(2)
i
⟵
i
+
1;
(3)
FromX1, obtain a LS estimate of Φ^i:
Φ^i=X1C⋄Y^ei−1T†;
(4)
FromX2, obtain a LS estimate of Y^ei:
Y^ei=X2C⋄Φ^iT†;
(5)
Repeat steps 2–4 until convergence. The convergence is achieved whenεi−εi−1≤10−6, where εi is the residual error computed in the i-th iteration.


## 4.3. Alternative Closed-Form Solution to the First Stage

In contrast to the B-ALS receiver, the second proposed receiver named LS-KRF is an alternative closed-form solution that can be employed in particular cases in whichK≥UNMS. The idea is to filter the received signal tensor X by exploiting the knowledge of the multifrequency coding matrix C and then solve a set of rank 1 approximation problems.Initially, by multiplying both sides ofX3T in (26) by the pseudoinverse of CT from the right-hand side, we obtain(30)Ye⋄Φ=X3TCT†=y1⊗φ1,…,yUNMS⊗φUNMS∈ℂQT×UNMS.According to property in (4), the i-th column of (30) i=1,…,UNMS can be rewritten as(31)yi⊗φi=vecφi∘yi,which denotes the vectorization operation of the rank 1 matrix Ψi=φi∘yi∈ℂQ×T. By defining UiΣiViH as the singular value decomposition (SVD) of Ψi, estimates for φi∈ℂQ and yi∈ℂT i=1,…,UNMS can be obtained by truncating the SVD of Ψi to a rank 1 approximation, i.e., [27](32)φ^i=σ1u1,y^i=σ1v1∗,where u1∈ℂQ and v1∈ℂT are the corresponding first left and right singular vectors of Ui and Vi, respectively. σ1 denotes the largest singular value of the matrix Σi. Final estimates for Φ^ and Y^e are obtained by repeating this SVD computation i=1,…,UNMS times in parallel, one for each column of (30). The pseudocode of the LS-KRF receiver is summarized in Algorithm 2.Algorithm 2:Pseudocode of the LS-KRF receiver.
(1)
Apply theunvecQ×T operator in the i-th column of (30) to obtain the rank 1 matrix Ψi∈ℂQ×T;
(2)
Compute the SVDΨi=UiΣiViH. Then, obtain the estimates for the i-th columns of Φ^ and Y^e as:
φ^i=σ1u1andy^i=σ1v1∗,
whereu1∈ℂQ and v1∈ℂT are the first left and right singular vectors of Ui and Vi, respectively, and σ1 is the largest singular value;
(3)
Repeat steps 1-2 for all columns of (30).


## 4.4. Second Stage: Sparse Formulation to DL and UL Channel Parameters Estimation

Once the matricesΦ^ and Y^e are estimated, the second stage of the proposed receivers consists in estimating the channel parameters (AoDs, AoAs, and path gains) to reconstruct the channel matrices H^u and G^u related to each MS. Thanks to the knowledge of the multifrequency coding matrix C, the estimated factor matrices are not affected by permutation of columns ambiguity. Therefore, the first stage provides automatic separation of the compressed users’ channels. From Φ^ and Y^e, the multiuser channel estimation problem can be decoupled into U single-user ones as formulated below.Let us rewrite the block representation ofΦ^, defined as(33)Φ^=Φ^1,Φ^2,…,Φ^U∈ℂQ×UNMS,where(34)Φ^u=FHG^u∈ℂQ×NMS,u=1,…,U.where G^u∈ℂNBS×NMS denotes the estimated UL channel matrix related to the u-th MS. By replacing G^u for (19), and then vectorizing (34) according to the property in (3), we have(35)φ^u=vecΦ^u=A^MS⋄FHA^BSβ^.Using the property in (2), we straightforwardly obtain(36)φ^u=INMS⊗FHA^MS⋄A^BSβ^∈ℂQNMS,where INMS denotes an identity matrix of size NMS×NMS.The same procedure can directly be applied in theu-th block Y^uT of the estimated factor matrix:(37)Y^e=Y^1T,Y^2T,…,Y^UT∈ℂT×UNMS,where(38)Y^uT=STWTH^uT+VuDLT∈ℂT×NMS.H
^
u
∈
ℂ
N
MS
×
N
BS denotes the estimated DL channel matrix defined above in (17). Similar to (36), we obtain the following vector formulation:(39)y^u=INMS⊗STWTA^MS⋄A^BSα^+vuDL∈ℂTNMS,where y^u=vecY^uT and vuDL=vecVuDLT.From (36) and (39), two independent CS problems can be formulated to jointly estimate the parameters of the DL and UL channels of the u-th MS. We assume that grid quantization errors are neglected, i.e., the AoDs and AoAs are drawn from a uniform angle grid of N points contained in the set 0,2π/N,…,2πN−1/N, with N≫Lu and N≫Mu. Based on this assumption, we can obtain the following sparse formulations for the vectors y^u and φ^u, respectively:(40)y^u=INMS⊗STWTΣDα¯,(41)φ^u=INMS⊗FHΣDβ¯,where ΣD∈ℂNMSNBS×N2 denotes the known dictionary matrix used to solve the sparse signal recovery problem, defined as(42)ΣD=A¯MS⊗A¯BS∈ℂNMSNBS×N2.The matricesA¯MS∈ℂNMS×N and A¯BS∈ℂNBS×N that make the dictionary are given by(43)A¯MS=aMS0,aMS2πN,…,aMS2πN−1N,A¯BS=aBS0,aBS2πN,…,aBS2πN−1N,and contain all points of the uniform angle grid. The left-hand side matrices INMS⊗STWT and INMS⊗FH are called measurement matrices of the sparse problems. α¯∈ℂN2 and β¯∈ℂN2 are sparse vectors obtained by augmenting the vector gains α and β with zero elements, respectively.Estimates for the parameters of the channel matricesHu and Gu can be obtained by applying CS algorithms in the estimated sparse vectors (40) and (41). Many efficient algorithms such as orthogonal matching pursuit (OMP) [22], structured compressive sampling matching pursuit (S-CoSaMP) [25], and fast iterative shrinkage-thresholding (FISTA) [23], to name a few, can be used to solve these two sparse signal recovery problems. In a simplified view, the estimates for the path gains α^ and β^ correspond to nonzero entries of the estimated sparse vectors α¯ and β¯, while estimation for the spatial parameters (AoDs and AoAs) are obtained by selecting the columns of the dictionary matrix ΣD related to the positions of the estimated path gains in the sparse vector. In our numerical results, we adopt the OMP algorithm to estimate the UL and DL channel parameters for simplicity reasons, although any state-of-the-art CS-based algorithm is equally applicable to solve problems (40) and (41).Remark:.
Compared to (36), the sparse signal recovery problem formulated from (39) naturally incorporates the DL noise contribution in its structure. For this reason, the proposed closed-loop framework for channel estimation can lead to some performance degradation in the DL channel estimation compared to UL channel estimation. Therefore, we can observe a trade-off between DL channel estimation accuracy performed by the BS and reduction of the processing cost for channel estimation at the MS side. This discussion is reinforced by means of numerical simulations in Section 6.


## 4.5. Joint DL and UL Channel Estimation

Finally, from the estimated channel parameters (AoDs, AoAs, and path gains), the BS can construct the estimated DL and UL channel matricesH^u and G^u of the u-th MS according to relations (17) and (19) as follows:(44)H^u=A^MSdiagα^A^BST,(45)G^u=A^BSdiagβ^A^MST.As previously presented in Section4.4, the proposed receivers decouple the multiuser channel estimation problem into 2U single-user ones (U problems dedicated to each communication link) that can be solved independently for each MS. Since 2U digital processing units are available at the BS, the second stage of the proposed receivers can be computed in parallel. Therefore, its processing delay can be kept constant (i.e., it does not increase with the number of MSs), when the BS is equipped with multiple (at least 2U) digital processing units. The parallelized processing for the second stage of the proposed receivers is illustrated in Figures 3 and 4. In addition, the overall pseudocode of the proposed two-stage tensor-based receivers for joint DL and UL channel estimation is summarized in Algorithm 3.Figure 3
Illustration of the parallelized processing for the estimation of U uplink channels. Theu-th block A^u of the estimated factor matrix A^ is forwarded and processed for a dedicated processor.Figure 4
Illustration of the parallelized processing for the estimation of the U downlink channels. Theu-th block Y^uT of the estimated factor matrix B^ is forwarded and processed for a dedicated processor.Algorithm 3:Overall pseudocode of the proposed receivers for joint DL and UL channel estimation.
1
First Stage. estimation of the compressed channel matrices
(1.1)
From the received signal tensorX in (20), obtain the estimated factor matrices Φ^ and Y^e via B-ALS or LS-KRF described in Algorithms 1 and 2, respectively,
(2)
Second Stage. parameters estimation and channels reconstruction
(2.1)
FromΦ^=Φ^1,Φ^2,…,Φ^UandY^e=Y^1T,Y^2T,…,Y^UT, obtain the estimates of the channel parameters (AoDs, AoAs, and path gains) for each MS by applying a CS recovery algorithm (e.g., OMP) to problems (40) and (41) independently
(2.2)
The BS constructs the estimated DL and UL channel matricesH^u and G^u according to (44) and (45), respectively.


## 5. Identifiability Issues

In this section, we examine the identifiability issues under which the compressed DL and UL channel matricesΦ and Ye can be jointly and uniquely recovered using the proposed receivers.


### 5.1. B-ALS Receiver

Unique LS solutions for the compressed DL and UL channel matricesΦ and Ye obtained from (24) and (25) require that C⋄YeT∈ℂUNMS×KT and C⋄ΦT∈ℂUNMS×KQ have full row rank to be right invertible. Hence, the following two conditions must be satisfied:(46)KT≥UNMS,KQ≥UNMS.Combining these conditions yields the following lower bound on the number of subcarriers required for the multifrequency coding at the MS:(47)K≥maxUNMST,UNMSQ,where x denotes the smallest integer number that is greater or equal to x.


### 5.2. LS-KRF Receiver

The LS-KRF receiver requires that the following necessary and sufficient uniqueness condition be satisfied:(48)K≥UNMS.Note that this condition indicates that the application of the LS-KRF receiver requires a more restricted scenario compared to the proposed B-ALS receiver since the number of frequency resources (subcarriers) increases with the number of antennas at the MSs and active MSs. On the contrary, the LS-KRF receiver is a closed-form solution in contrast to the iterative B-ALS receiver.


## 5.1. B-ALS Receiver

Unique LS solutions for the compressed DL and UL channel matricesΦ and Ye obtained from (24) and (25) require that C⋄YeT∈ℂUNMS×KT and C⋄ΦT∈ℂUNMS×KQ have full row rank to be right invertible. Hence, the following two conditions must be satisfied:(46)KT≥UNMS,KQ≥UNMS.Combining these conditions yields the following lower bound on the number of subcarriers required for the multifrequency coding at the MS:(47)K≥maxUNMST,UNMSQ,where x denotes the smallest integer number that is greater or equal to x.


## 5.2. LS-KRF Receiver

The LS-KRF receiver requires that the following necessary and sufficient uniqueness condition be satisfied:(48)K≥UNMS.Note that this condition indicates that the application of the LS-KRF receiver requires a more restricted scenario compared to the proposed B-ALS receiver since the number of frequency resources (subcarriers) increases with the number of antennas at the MSs and active MSs. On the contrary, the LS-KRF receiver is a closed-form solution in contrast to the iterative B-ALS receiver.


## 6. Simulation Results

In this section, we present a set of simulation results to evaluate the performance of the proposed joint DL and UL channel estimator. We compare the proposed channel training framework with the conventional training framework illustrated in Figure1, where the CS-based OMP algorithm [22] is applied at both MSs and BS to estimate the channel parameters in a decoupled way. The OMP algorithm is also considered as the second stage of our algorithm, according Section 4.4. The MSs and the BS employ uniform linear arrays with half-wavelength-spaced antennas. We set NBS=32, NMS=16, U=2, N=64, and equal signal-to-noise ratio (SNR) for the DL and UL communications in all experiments. The obtained results are averaged over 1000 independent Monte Carlo runs. At each run, the DL and UL channel matrices with Lu=3 and Mu=3 paths per user and HB matrices are generated in accordance with equations (17), (19), and (12), respectively. The pilot signal S is a binary phase shift keying (BPSK) modulated matrix, and the multifrequency coding matrix C has random coefficients following a uniform distribution.The receiver’s performance is evaluated in terms of the normalized mean square error (NMSE) measures between the estimated and true DL and UL channel matrices:(49)NMSEH^=∑u=1UHu−H^uF2∑u=1UHuF2,NMSEG^=∑u=1UGu−G^uF2∑u=1UGuF2.In our experiments, we evaluate the accuracy of channel estimation in terms of the NMSE metric for different values of signal-to-noise ratio (SNR), number of transmission (P) and reception (Q) beams, number of training subcarriers (K), and length of the pilot sequences (T).Figures5 and 6 show the NMSE as a function of the number of transmission (P) and reception (Q) beams for different values of SNR, and fixed values U=2, T=16, and K=25. According Figure 6, the proposed method outperforms the classical framework to the UL channel estimation, while the DL performance is worse in all the simulated SNR ranges as shown in Figure 5. From this experiment, we can observe the trade-off between DL channel estimation accuracy and computational complexity. In other words, the proposed framework concentrates most of the processing burden for channel estimation at the BS side, while a better performance of the DL channel estimation comes at the expense of a high computational cost to complex channel estimation processing at the MS side when the conventional framework is utilized. On the contrary, the performance loss at DL is compensated with more accurate estimations at UL. In addition, the NMSE performance is not influenced by the number of RF chains when P and Q are greater than 12. This result reveals that the proposed framework provides a good channel estimation accuracy even when the BS is equipped with a few number of RF chains, which is the case in HB architectures.Figure 5
NMSE vs. number of transmission beams P for the DL channel estimation:U=2, T=16, and K=25.Figure 6
NMSE vs. number of reception beamsQ for the UL channel estimation: U=2, T=16, and K=25.In Figures7 and 8, the NMSE performance is evaluated as a function of the number of subcarriers (K). An increase of K leads to an improved performance only until K=32 subcarriers in the DL channel estimation, while for the UL channel estimation this value is approximately equal to K=64 subcarriers. This result shows that the proposed closed-loop channel training framework can operate with few frequency resources to jointly estimate the MS channels with high accuracy.Figure 7
NMSE vs. number of subcarriers K for the DL channel estimation:U=2 and T=16.Figure 8
NMSE vs. number of subcarriers K for the UL channel estimation:U=2 and T=16.Figures9 and 10 depict the NMSE performance in terms of the length of the pilot sequence (T). In this experiment, we also set U=2, P=8, and K=25. Here, we conclude that short pilot sequences are necessary to estimate the DL and UL channels from the proposed method. Low variability in the NMSE values is observed when T>15. For a massive MIMO scenario, this result implies in a substantial reduction in the pilot overhead to joint DL and UL channel estimation.Figure 9
NMSE vs. length of the training sequenceT for the DL channel estimation: U=2, P=8, and K=25.Figure 10
NMSE vs. length of the training sequenceT for the UL channel estimation: U=2, P=8, and K=25.


## 7. Conclusion and Perspectives

In this paper, we have addressed the joint DL and UL channel estimation problem for multiuser FDD massive MIMO systems with HB architecture. As contributions of this work, we firstly proposed a novel closed-loop and multifrequency-based channel training framework that concentrates most of the processing burden for channel estimation at the BS side. We have shown that making use of the proposed framework, the received closed-loop signal follows a third-order PARAFAC model, which can be exploited by two tensor-based semiblind receivers followed by compressed sensing recovery of the channel parameters. Additionally, we have also provided an identifiability study. We have compared our proposed approach with the conventional channel training framework, where the DL and UL channel estimation problems are treated as two decoupled problems, i.e., solved by the MSs and BS, separately. Compared to the conventional framework, the proposed receivers have shown a superior performance in the estimation of the UL channel, while the performance of the DL channel estimation exhibits some degradation. It is worth noting that such a degradation is the price to pay for the complexity reduction at the MS by transferring the processing burden associated with the DL channel estimation to the BS. Perspectives include the extension of the proposed modeling to frequency- and time-selective channels.



---

*Source: 4858137-2019-09-15.xml*