Asymptotic Analysis of RZF over Double Scattering Channels with MMSE Estimation

TL;DR

This paper provides an asymptotic analysis of RZF precoding in multi-user MISO systems with double scattering channels, deriving deterministic approximations for SINR and ergodic rate, and revealing the maximum number of supported users.

Contribution

It introduces the first MMSE channel estimate for double scattering models and derives tight deterministic approximations for system performance metrics in large-scale regimes.

Findings

Maximum supported users equals the number of scatterers.

Asymptotic approximations closely match simulations for moderate sizes.

Supports large-scale MIMO gains with double scattering channels.

Abstract

This paper studies the ergodic rate performance of regularized zero-forcing (RZF) precoding in the downlink of a multi-user multiple-input single-output (MISO) system, where the channel between the base station (BS) and each user is modeled by the double scattering model. This non-Gaussian channel model is a function of both the antenna correlation and the structure of scattering in the propagation environment. This paper makes the preliminary contribution of deriving the minimum-mean-square-error (MMSE) channel estimate for this model. Then under the assumption that the users are divided into groups of common correlation matrices, this paper derives deterministic approximations of the signal-to-interference-plus-noise ratio (SINR) and the ergodic rate, which are almost surely tight in the limit that the number of BS antennas, the number of users and the number of scatterers in each group grow infinitely large. The derived results are expressed in a closed-form for the special case of multi-keyhole channels. Simulation results confirm the close match provided by the asymptotic analysis for moderate system dimensions. We show that the maximum number of users that can be supported simultaneously, while realizing large-scale MIMO gains, is equal to the number of scatterers.

Equations

$$y_{k,g}={\text{h}}_{k,g}^H\text{x}+n_{k,g},k=1,\dots ,K_g,g=1,\dots ,G,$$

$$\text{x}=\sum _{g=1}^G\sum _{k=1}^{K_g}\sqrt{p_{k,g}}{\text{g}}_{k,g}{s}_{k,g},$$

$$\mathbb{E}[||\text{x}|{|}^2]=\mathbb{E}[\text{tr }(\text{P}{\text{G}}^H\text{G})]\le \bar{P},$$

$${\text{h}}_{k,g}=\sqrt{S_g}\left(\frac{1}{\sqrt{S_g}}{\text{R}}_{B{S}_g}^{1/2}{\text{W}}_g{\bar{\text{S}}}_g^{1/2}\right){\tilde{\text{w}}}_{k,g},$$

$${\text{y}}_{k,g}^{tr}={\text{h}}_{k,g}+\frac{1}{\sqrt{{\rho }_{tr}}}{\text{n}}_{k,g}^{tr},$$

$${\hat{\text{h}}}_{k,g}=d_g{\text{R}}_{B{S}_g}{\text{Q}}_g{\text{y}}_{k,g}^{tr},$$

$$R_{k,g}=\log (1+{\gamma }_{k,g}),$$

$${\gamma }_{k,g}=\frac{p_{k,g}|\mathbb{E}[{\text{h}}_{k,g}^H{\text{g}}_{k,g}]{|}^2}{p_{k,g}\text{var}[{\text{h}}_{k,g}^H{\text{g}}_{k,g}]+{\sum }_{(k^{\prime },g^{\prime })\ne (k,g)}{p}_{k^{\prime },g^{\prime }}\mathbb{E}[|{\text{h}}_{k,g}^H{\text{g}}_{k^{\prime },g^{\prime }}{|}^2]+{\sigma }^2}.$$

$$R_{sum}=\sum _{g=1}^G\sum _{k=1}^{K_g}{R}_{k,g}.$$

$$\text{G}=\zeta ({\hat{\text{H}}}^H\hat{\text{H}}+K\alpha {\text{I}}_N{)}^{-1}{\hat{\text{H}}}^H,$$

$${\zeta }^2=\bar{P}/\mathbb{E}[\text{tr}(\text{P}\hat{\text{H}}({\hat{\text{H}}}^H\hat{\text{H}}+K\alpha {\text{I}}_N{)}^{-2}{\hat{\text{H}}}^H)]=\bar{P}/\Theta ,$$

$${\gamma }_{k,g\text{RZF}}=\frac{p_{k,g}|\mathbb{E}[{\text{h}}_{k,g}^H\hat{\text{V}}{\hat{\text{h}}}_{k,g}]{|}^2}{\mathbb{E}[{\text{h}}_{k,g}^H\hat{\text{V}}{\hat{\text{H}}}_{[k,g]}^H{\text{P}}_{[k,g]}{\hat{\text{H}}}_{[k,g]}\hat{\text{V}}{\text{h}}_{k,g}]+p_{k,g}\text{var }[{\text{h}}_{k,g}^H\hat{\text{V}}{\hat{\text{h}}}_{k,g}]+\frac{\Theta }{\rho }},$$

$$0<\text{lim inf }\frac{S_g}{N}\le \text{lim sup }\frac{S_g}{N}<\infty ,0<\text{lim inf }\frac{K_g}{N}\le \text{lim sup }\frac{K_g}{N}<\infty .$$

$$X_N-X_N^o\underset{N\to \infty }{\overset{a.s.}{\to }}0.$$

$${\gamma }_{k,g\text{RZF}}-{\gamma }_{k,g\text{RZF}}^o\underset{N\to \infty }{\overset{a.s.}{\to }}0.$$

$${\text{h}}_{k,g}=\sqrt{S_g}{\text{Z}}_g{\tilde{\text{w}}}_{k,g},$$

$${\hat{\text{h}}}_{k,g}={\mathbf{\Phi }}_g^{1/2}{\bar{\text{q}}}_{k,g},$$

$${\mathbf{\Phi }}_g=d_g^2{\text{R}}_{B{S}_g}{\text{Q}}_g\left({\text{Z}}_g{\text{Z}}_g^H+\frac{1}{{\rho }_{tr}}{\text{I}}_N\right){\text{Q}}_g^H{\text{R}}_{B{S}_g}^H.$$

$$m_g(\alpha )=\frac{1}{K}{d}_g^2\left(\sum _{j=1}^{S_g}\frac{{\bar{s}}_{g,j}{\delta }_g({\bar{\text{D}}}_g,\alpha )}{1+\frac{K_g}{K}{d}_g^2{\bar{s}}_{g,j}{\bar{m}}_g(\alpha ){\delta }_g({\bar{\text{D}}}_g,\alpha )}+\frac{S_g}{{\rho }_{tr}}{\delta }_g({\tilde{\text{D}}}_g,\alpha )\right),$$

$${\bar{m}}_g(\alpha )=\frac{1}{1+m_g(\alpha )},$$

$${\delta }_g({\text{R}}_g,\alpha )=\frac{1}{S_g}\text{tr }{\text{R}}_g\bar{\text{T}}(\alpha ),$$

$$\text{where, }\bar{\text{T}}(\alpha )={\left(\sum _{i=1}^G\frac{{\bar{\text{D}}}_i}{S_i}\left(\sum _{j=1}^{S_i}\frac{\frac{K_i}{K}{d}_i^2{\bar{s}}_{i,j}{\bar{m}}_i(\alpha )}{1+\frac{K_i}{K}{d}_i^2{\bar{s}}_{i,j}{\bar{m}}_i(\alpha ){\delta }_i({\bar{\text{D}}}_i,\alpha )}\right)+\frac{K_i}{K}{d}_i^2\frac{{\bar{m}}_i(\alpha )}{{\rho }_{tr}}{\tilde{\text{D}}}_i+\alpha {\text{I}}_N\right)}^{-1},$$

$$\frac{1}{K}\text{tr}(\text{U}{\hat{\text{C}}}^{-1}(\alpha ))-\frac{1}{K}\text{tr}(\text{U}\bar{\text{T}}(\alpha ))\underset{N\to \infty }{\overset{a.s.}{\to }}0.$$

$$\frac{1}{K}{\hat{\text{h}}}_{k,g}^H{\hat{\text{C}}}_{[k,g]}^{-1}{\hat{\text{h}}}_{k,g}-m_g(\alpha )\underset{N\to \infty }{\overset{a.s.}{\to }}0,$$

$$\frac{1}{K}{\text{h}}_{k,g}^H{\hat{\text{C}}}_{[k,g]}^{-1}{\hat{\text{h}}}_{k,g}-h_g(\alpha )\underset{N\to \infty }{\overset{a.s.}{\to }}0,$$

$$\mathbf{\chi }(\alpha )-({\text{I}}_N-\bar{\text{J}}(\alpha ){)}^{-1}\bar{\text{v}}(\alpha )\underset{N\to \infty }{\overset{a.s.}{\to }}0,$$

$$[\bar{\text{J}}(\alpha ){]}_{g,i}=d_g^2\frac{K_i}{K}\frac{1}{(1+m_i(\alpha ){)}^2}\left({\bar{\beta }}_{g,i}({\text{R}}_{B{S}_g}{\text{Q}}_g,{\text{Q}}_g^H{\text{R}}_{BS}^H,\alpha )+\frac{1}{{\rho }_{tr}}{\tilde{\beta }}_i({\text{R}}_{B{S}_g}{\text{Q}}_g,{\text{Q}}_g^H{\text{R}}_{BS}^H,\alpha )\right)$$

$${\bar{\text{v}}}_g(\alpha )=d_g^2\frac{1}{K}\sum _{i=1}^G\sum _{l=1}^{K_i}\frac{p_{l,i}}{(1+m_i(\alpha ){)}^2}\left({\bar{\beta }}_{g,i}({\text{R}}_{B{S}_g}{\text{Q}}_g,{\text{Q}}_g^H{\text{R}}_{BS}^H,\alpha )+\frac{1}{{\rho }_{tr}}{\tilde{\beta }}_i({\text{R}}_{B{S}_g}{\text{Q}}_g,{\text{Q}}_g^H{\text{R}}_{BS}^H,\alpha )\right),$$

$$\text{where }{\bar{\beta }}_{g,i}(\text{A},\text{B},\alpha )={\bar{\beta }}_{g,i}^1(\text{A},\text{B},\alpha )+{\bar{\beta }}_{g,i}^2(\text{A},\text{B},\alpha ),$$

$$\displaystyle\bar{\beta}^{1}_{g,i}(\textbf{A},\textbf{B},\alpha)=\begin{cases}\frac{d_{i}^{2}}{S_{i}K}\sum_{j=1}^{S_{i}}\frac{\bar{s}_{i,j}}{(1+\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}(\alpha)\bar{s}_{i,j}\delta_{i}(\bar{\textbf{D}}_{i},\alpha))^{2}}\sum_{n=1}^{S_{g}}\frac{\bar{s}_{g,n}u^{\prime}_{g}(\textbf{A}\textbf{R}_{BS_{g}}\textbf{B},\bar{\textbf{D}}_{i},\alpha)}{\left(1+\frac{K_{g}}{K}d_{g}^{2}\bar{m}_{g}(\alpha)\bar{s}_{g,n}\delta_{g}(\bar{\textbf{D}}_{g},\alpha)\right)^{2}},&\text{ if }g\neq i,\\ \frac{d_{i}^{2}}{K}\sum_{j=1}^{S_{i}}\frac{\bar{s}_{i,j}}{\left(1+\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}(\alpha)\bar{s}_{i,j}\delta_{i}(\bar{\textbf{D}}_{i},\alpha)\right)^{2}}\Big{(}\underset{n\neq j}{\sum_{n=1}^{S_{i}}}\frac{\frac{1}{S_{i}}\bar{s}_{i,n}u^{\prime}_{i}(\textbf{A}\textbf{R}_{BS_{i}}\textbf{B},\bar{\textbf{D}}_{i},\alpha)}{\left(1+\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}(\alpha)\bar{s}_{i,n}\delta_{i}(\bar{\textbf{D}}_{i},\alpha)\right)^{2}}\\ +\bar{s}_{i,j}\delta_{i}(\textbf{R}_{BS_{i}}\textbf{Q}_{i}\textbf{R}_{BS_{i}}\textbf{B},\alpha)\delta_{i}(\textbf{A}\textbf{R}_{BS_{i}}\textbf{Q}^{H}_{i}\textbf{R}^{H}_{BS_{i}},\alpha)\Big{)},&\text{ if }g=i,\end{cases}$$

$${\bar{\beta }}_{g,i}^2(\text{A},\text{B},\alpha )=\frac{d_i^2}{K{\rho }_{tr}}\sum _{j=1}^{S_g}\frac{{\bar{s}}_{g,j}{u}_g^{\prime }(\text{A}{\text{R}}_{B{S}_g}\text{B},{\tilde{\text{D}}}_i,\alpha )}{(1+\frac{K_g}{K}{d}_g^2{\bar{m}}_g(\alpha ){\bar{s}}_{g,j}{\delta }_g({\bar{\text{D}}}_g,\alpha ){)}^2},$$

$${\tilde{\beta }}_i(\text{A},\text{B},\alpha )=\frac{d_i^2}{K}\left(\sum _{j=1}^{S_i}\frac{{\bar{s}}_{i,j}{u}_i^{\prime }({\bar{\text{D}}}_i,\text{A}\text{B},\alpha )}{(1+\frac{K_i}{K}{d}_i^2{\bar{m}}_i(\alpha ){\bar{s}}_{i,j}{\delta }_i({\bar{\text{D}}}_i,\alpha ){)}^2}+\frac{S_i}{{\rho }_{tr}}{u}_i^{\prime }({\tilde{\text{D}}}_i,\text{A}\text{B},\alpha )\right),$$

$$u_g^{\prime }({\text{R}}_g,\text{L},\alpha )=\frac{1}{S_g}\text{tr }{\text{R}}_g\left[\bar{\text{T}}\left(\sum _{z=1}^G\frac{{\bar{\text{D}}}_z(d_z^2{K}_z{\bar{m}}_z(\alpha ){)}^2}{S_z{K}^2}\text{tr}({\bar{\text{S}}}_z{\text{W}}_z(\alpha {)}^2{\bar{\text{S}}}_z)u_z^{\prime }({\bar{\text{D}}}_z,\text{L},\alpha )+\text{L}\right)\bar{\text{T}}\right],$$

$${\text{W}}_i(\alpha )={\left({\text{I}}_{S_i}+\frac{K_i}{K}{d}_i^2{\bar{m}}_i(\alpha ){\delta }_i({\bar{\text{D}}}_i,\alpha ){\bar{\text{S}}}_i\right)}^{-1},$$

$${\text{u}}^{\prime }(\bar{\text{D}},\text{L},\alpha )=({\text{I}}_N-\text{J}(\bar{\text{D}},\alpha ){)}^{-1}\text{v}(\bar{\text{D}},\text{L},\alpha ),$$

Full text

I Introduction

Large-scale multiple-input multiple-output (MIMO) is widely considered as a promising technology for next generation wireless communication systems [1, 2, 3, 4, 5, 6, 7, 8, 9]. However, most works on this subject share the underlying assumption of rich scattering conditions and, thus, work with full rank Rayleigh or Rician fading channel matrices. Although the use of full rank channel matrices facilitates the derivation of closed-form capacity bounds and approximations, these models do not capture the characteristics of realistic propagation environments, where the presence of spatial correlation and poor scattering conditions significantly affects the system performance [10, 11, 12].

Many works have already considered correlated Rayleigh fading channel models in their analysis of large-scale and Massive MIMO systems [2, 7, 3, 4, 5, 8]. The authors in [2] show that in the limit of an infinite number of antennas, the effects of noise, channel estimation errors, and interference vanish, while pilot contamination remains the only performance limitation. However, very recently the authors in [3] show in the same limit that the capacity increases without bound even under pilot contamination, provided that the spatial correlation matrices of the pilot contaminating users are asymptotically linearly independent and multi-cell minimum-mean-square-error (MMSE) precoding/combining techniques are employed.

Despite the important role played by spatial correlation in determining large-scale MIMO performance, low rank channels have been observed in MIMO systems that have low antenna correlation at both ends of the transmission link [13]. In fact, Gesbert et al. show in [10] that MIMO capacity is governed by both the spatial correlation at the communication ends and the structure of scattering in the propagation environment. Motivated by this, the authors devise a “double scattering channel model”, which utilizes the geometry of the propagation environment to model spatial correlation, rank deficiency and limited scattering. Unlike the commonly utilized MIMO channel models, the double scattering model is non-Gaussian, rendering the theoretical analysis of techniques designed using this model extremely hard.

I-A Related Literature

The main literature related to the double scattering channel model is scarce and is mainly represented by [14, 15, 16, 17, 18, 19, 13, 20, 21]. The authors in [14] study its diversity order and show that a MIMO system with $t$ transmit (Tx) and $r$ receive (Rx) antennas and $s$ scatterers achieves a diversity of order $trs/\max(t,r,s)$. The authors in [16] analyze the ergodic MIMO capacity taking into account the presence of double scattering and keyhole effects and show that the use of multiple antennas in keyhole channels only offers diversity gains, but no spatial multiplexing gains. Closed-form upper-bounds on the sum-capacity are obtained in [19] for the MIMO multiple access channel with double scattering fading.

A few papers have analyzed the double scattering model in large-scale MIMO settings [13, 21, 22]. The authors in [13] study this model without Tx and Rx correlation using tools from free probability theory and derive implicit expressions for the asymptotic mutual information (MI) and signal-to-interference-plus-noise ratio (SINR) of the MMSE detector. The authors in [21] study a MIMO multiple access system with double-scattering channels and use random matrix theory (RMT) tools to derive an almost surely tight deterministic approximation of the MI. However, none of these papers consider practical large-scale multi-user MIMO systems with linear signal processing schemes and channel estimation. The authors are only aware of [23], that deals with these aspects using numerical examples instead of theoretical analysis.

I-B Main Contributions

The focus of this work is on the downlink (DL) of a single-cell multi-user multiple-input single-output (MISO) system. We consider double scattering fading between the BS and the users in its most general form. The users are further divided into $G$ groups, such that the users in the same group are characterized by common correlation matrices. This assumption is needed to provide some level of analytical tractability for an, otherwise, extremely difficult to study non-Gaussian channel model. We stress here that having variation across the correlation matrices is important in realizing the fundamental performance of large-scale MIMO systems as shown in [3, 4, 9]. However, these works base their analysis on the correlated Rayleigh channel model, which encourages us to study the benefits of having a large number of antennas under the more realistic double scattering model, where we focus on studying the effects of the scattering conditions (instead of the structure of the correlation matrices) on the system performance.

We consider a realistic large-scale MIMO system where the BS obtains channel state information (CSI) from uplink pilot transmissions and applies MMSE estimation technique. Under the assumption that the number of BS antennas, scatterers and users grow large, we derive asymptotically tight deterministic approximations of the SINR and the ergodic rate with regularized zero-forcing (RZF) precoding, which are accurate for moderate system dimensions as well, as shown through simulations. The deterministic approximations are expressed in a closed-form for the multi-keyhole channel and some important insights into the impact of different system parameters on the sum rate performance are drawn. We show that the number of users scheduled simultaneously should be less than the number of scatterers to see ‘large-scale MIMO gains’.

I-C Outline and Notation

The rest of the paper is organized as follows. Section II presents the transmission model and introduces the double scattering channel model along with its MMSE estimate. In Section III, the asymptotically tight deterministic equivalents of the SINR and user rates with RZF precoding are derived. Simulation results are provided in Section IV and Section V concludes the paper.

The following notation is used. Boldface lower-case and upper-case characters denote vectors and matrices respectively. The operator $\text{tr }(\textbf{X})$ denotes the trace of a matrix X. The spectral norm of a matrix X is denoted by $||\textbf{X}||$ and the $N\times N$ diagonal matrix of entries $\{x_{n}\}$ is denoted by $\textbf{X}=\text{diag}(x_{1},x_{2},\dots,x_{N})$. A random vector $\textbf{x}\sim\mathcal{CN}(\textbf{m},\boldsymbol{\Phi})$ is complex Gaussian distributed with mean m and covariance matrix $\boldsymbol{\Phi}$. The notation $\xrightarrow[]{a.s.}$ denotes almost sure convergence.

II System Model

Consider a single-cell multi-user MISO system where a BS equipped with $N$ antennas serves $K$ single-antenna users, who are divided into $G$ groups of $K_{g}$, $g=1,\dots,G$, users such that the users in the same group experience similar propagation conditions. In this section, we outline the transmission model, and introduce the double scattering channel along with its MMSE estimate.

II-A Transmission Model

The received complex baseband signal $y_{k,g}$ at user $k$ in group $g$ is given as,

[TABLE]

where $\textbf{h}_{k,g}^{H}\in\mathbb{C}^{1\times N}$ is the channel vector from the BS to user $k$ in group $g$, $\textbf{x}\in\mathbb{C}^{N\times 1}$ is the Tx signal vector and $n_{k,g}\sim\mathcal{CN}(0,\sigma^{2})$ is the receiver noise. The Tx signal vector x is given as,

[TABLE]

where $\textbf{g}_{k,g}\in\mathbb{C}^{N\times 1}$, $p_{k,g}\geq 0$ and $s_{k,g}\sim\mathcal{CN}(0,1)$ are the precoding vector, signal power and data symbol for user $k$ in group $g$ respectively. The precoding vectors satisfy,

[TABLE]

where $\bar{P}>0$ is the average total Tx power, $\textbf{P}=\text{diag}(p_{1,1},\dots,p_{K_{1},1},p_{1,2},\dots,p_{K_{G},G})$ $\in\mathbb{R}^{K\times K}$ and $\textbf{G}=[\textbf{G}_{1},\dots,\textbf{G}_{G}]\in\mathbb{C}^{N\times K}$ is the precoding matrix, where $\textbf{G}_{g}=[\textbf{g}_{1,g},\dots,\textbf{g}_{K_{g},g}]$ $\in\mathbb{C}^{N\times K_{g}}$.

II-B Double Scattering Channel Model

A main contribution of this paper is to apply the double scattering channel model proposed in [10] to a multi-user MISO system. This non-Gaussian model has rank that is determined by both the spatial correlation between the antennas at the BS and the structure of scattering in the propagation environment. The double scattering channel vector $\textbf{h}_{k,g}$ is given as [10],

[TABLE]

where $S_{g}$ is the number of scatterers at the Tx and the Rx sides in group $g$, $\textbf{R}_{BS_{g}}\in\mathbb{C}^{N\times N}$ is the correlation matrix between the BS antennas and the $S_{g}$ Tx scatterers, $\bar{\textbf{S}}_{g}\in\mathbb{C}^{S_{g}\times S_{g}}$ is the correlation matrix between the $S_{g}$ Tx and Rx scatterers, $\textbf{W}_{g}\in\mathbb{C}^{N\times S_{g}}$ is a standard complex Gaussian matrix that describes the small-scale fading between the BS and the scattering cluster at the Tx side and $\tilde{\textbf{w}}_{k,g}\sim\mathcal{CN}(\boldsymbol{0},\frac{1}{S_{g}}\textbf{I}_{S_{g}})\in\mathbb{C}^{S_{g}\times 1}$ describes the small-scale fading between the user $k$ in group $g$ and the scattering cluster at the Rx side. Note that we can assume $\bar{\textbf{S}}_{g}$ to be diagonal without any loss of generality for the statistics of the received signal.

Schematic of the double scattering channel is shown in Fig. 1, where $\sigma_{t,g}$ and $\sigma_{s,g}$ represent the angular spread of the radiated signal from the BS array and the Tx scatterers respectively in group $g$, and $\mu_{t,g}$ and $\mu_{s,g}$ represent the mean angle of departure (AoD) of the radiated signal from the BS array and the Tx scatterers respectively in group $g$, where $\mu_{t,g}=\mu_{s,g}$.

II-C Channel Estimation

During a dedicated uplink training phase, the users transmit mutually orthogonal pilot sequences that allow the BS to compute the MMSE estimates $\hat{\textbf{h}}_{k,g}$ of the channel vectors $\textbf{h}_{k,g}$. After correlating the received training signal with the pilot sequence of user $k$ in group $g$, the BS estimates the channel vector $\textbf{h}_{k,g}$ based on the received observation, $\textbf{y}^{tr}_{k,g}\in\mathbb{C}^{N\times 1}$, given as,

[TABLE]

where $\textbf{n}^{tr}_{k,g}\sim\mathcal{CN}(\boldsymbol{0},\textbf{I}_{N})$ and $\rho_{tr}>0$ is the effective training SNR, assumed to be given here.

*Lemma 1: * The MMSE estimate $\hat{\textbf{h}}_{k,g}\in\mathbb{C}^{N\times 1}$ of the channel $\textbf{h}_{k,g}$ in (4) is given as,

[TABLE]

where $d_{g}=\frac{1}{S_{g}}(\text{tr }\bar{\textbf{S}}_{g})$ and $\textbf{Q}_{g}=\left(d_{g}\textbf{R}_{BS_{g}}+\frac{1}{\rho_{tr}}\textbf{I}_{N}\right)^{-1}$.

The proof is postponed to Appendix A.

We stress that Lemma 1 has been derived for a non-Gaussian channel. Under the orthogonality property of the MMSE estimate, $\textbf{h}_{k,g}=\hat{\textbf{h}}_{k,g}+\tilde{\textbf{h}}_{k,g}$, where $\tilde{\textbf{h}}_{k,g}$ is the uncorrelated channel estimation error. Since our focus is on a single-cell system so the derived channel estimate is not corrupted by pilot contamination from adjacent cells. Therefore the considered system is referred to as a large-scale MIMO system instead of a Massive MIMO system.

II-D Achievable Rates

Since the BS does not have CSI, so it utilizes the MMSE estimate in (6) to implement digital precoding. The analysis presented in this paper is, therefore, a function of $\rho_{tr}$, which determines the error that is incurred in channel estimation. However, to compute the SINR in the downlink, the user needs to have CSI too. A well-known feature of large-scale MIMO systems is channel hardening, which means that the effective useful channel $\textbf{h}_{k,g}^{H}\textbf{g}_{k,g}$ of a user converges to its average value when $N,S$ grow large. Hence, it is sufficient for each user to have only the statistical CSI and the resulting performance loss vanishes in the large system limit. Using this idea, an ergodic achievable user rate can be computed using a technique from [24], widely applied to large-scale MIMO systems in [1, 25, 2, 8]. The main idea is to decompose $y_{k,g}$ in (1) as $\sqrt{p_{k,g}}\mathbb{E}[\textbf{h}_{k,g}^{H}\textbf{g}_{k,g}]s_{k,g}+\sqrt{p_{k,g}}(\textbf{h}_{k,g}^{H}\textbf{g}_{k,g}-\mathbb{E}[\textbf{h}_{k,g}^{H}\textbf{g}_{k,g}])s_{k,g}+\sum_{(k^{\prime},g^{\prime})\neq(k,g)}\sqrt{p_{k^{\prime},g^{\prime}}}\textbf{h}_{k,g}^{H}\textbf{g}_{k^{\prime},g^{\prime}}s_{k^{\prime},g^{\prime}}+{n}_{k,g}$ and assume that the average effective channel $\mathbb{E}[\textbf{h}_{k,g}^{H}\textbf{g}_{k,g}]$ is perfectly known at the corresponding user. By treating the inter-user interference and channel uncertainty as worst-case Gaussian noise, the user $k$ in group $g$ can achieve the ergodic rate,

[TABLE]

without knowing the instantaneous values of $\textbf{h}^{H}_{k,g}\textbf{g}_{k,g}$. The parameter $\gamma_{k,g}$ can be interpreted as the effective average SINR of user $k$ in group $g$ and is defined as,

[TABLE]

The ergodic achievable sum rate is given as,

[TABLE]

This paper considers RZF precoding, which is a state-of-the-art heuristic precoding scheme with a simple closed-form expression given as [2, 7],

[TABLE]

where $\hat{\textbf{H}}=[\hat{\textbf{H}}_{1}^{H},\hat{\textbf{H}}_{2}^{H},\dots,\hat{\textbf{H}}_{G}^{H}]^{H}\in\mathbb{C}^{K\times N}$, where $\hat{\textbf{H}}_{g}=[\hat{\textbf{h}}_{1,g},\hat{\textbf{h}}_{2,g},$ $\dots,\hat{\textbf{h}}_{K_{g},g}]^{H}\in\mathbb{C}^{K_{g}\times N}$, $\alpha$ is the regularization parameter and $\zeta$ ensures that the power constraint in (3) is satisfied as,

[TABLE]

where $\Theta=\mathbb{E}[\text{tr}(\textbf{P}\hat{\textbf{H}}\hat{\textbf{V}}^{2}\hat{\textbf{H}}^{H})]$, where $\hat{\textbf{V}}=(\hat{\textbf{H}}^{H}\hat{\textbf{H}}+K\alpha\textbf{I}_{N})^{-1}$. The SINR in (8) is now defined as,

[TABLE]

where $\rho=\frac{\bar{P}}{\sigma^{2}}$, $\hat{\textbf{H}}_{[k,g]}=[\hat{\textbf{H}}_{1}^{H},\dots,\hat{\textbf{H}}_{g-1}^{H},\hat{\textbf{h}}_{1,g},\dots,\hat{\textbf{h}}_{k-1,g},\hat{\textbf{h}}_{k+1,g},\dots,\hat{\textbf{h}}_{K_{g},g},\dots,\hat{\textbf{H}}_{G}^{H}]^{H}\in\mathbb{C}^{K-1\times N}$ and $\textbf{P}_{[k,g]}=\text{diag}(p_{1,1},\dots,p_{K_{g-1},g-1},p_{1,g},\dots,p_{k-1,g},p_{k+1,g},\dots,p_{K_{g},g},\dots,p_{K_{G},G})\in\mathbb{C}^{K-1,K-1}$.

III Main Results

The performance of the users are characterized by their ergodic rates, $R_{k,g}$, which depend on their SINRs. These ergodic rates, with RZF precoding, are difficult to compute for finite system dimensions. However, they tend to deterministic quantities in the large ($N,K$) regime, as shown in [2, 8, 7]. These quantities depend only on the statistics of the channels and are referred to as deterministic equivalents. These deterministic equivalents are almost surely (a.s.) tight in the asymptotic limit. By a.s. tight in the asymptotic limit we mean that as the system dimensions grow to infinity, the approximations yielded by these deterministic equivalents tend to the actual values with probability one. Deterministic equivalents were first proposed by Hachem et al. in [26], who showed their ability to capture important system performance indicators.

In the sequel, the asymptotic limit, denoted as $N\rightarrow\infty$, represents the following assumption.

A-1. For all $g$, $N$, $S_{g}$, $K_{g}$ and $K$ tend to infinity such that,

[TABLE]

Also, in the sequel, the deterministic equivalent of a sequence of RVs $X_{N}$ is represented by the deterministic sequence $X^{o}_{N}$, which approximates $X_{N}$ such that,

[TABLE]

The objective of this section is to derive the deterministic equivalent $\gamma^{o}_{k,g\text{RZF}}$ of the SINR $\gamma_{k,g\text{RZF}}$ with RZF precoding, such that,

[TABLE]

Note that $\gamma_{k,g\text{RZF}}$ is a function of both $\tilde{\textbf{w}}_{k,g}$ and $\textbf{W}_{g}$. From now on, these random vectors and matrices must be understood as sequences of vectors and matrices of growing dimensions. For the sake of simplicity, their dependence on $N$, $S$ and $K$ is not explicitly shown.

We stress again that the result in (14) must be understood in the way that, for each given set of system parameters $N$, $S$, and $K$, we provide approximations of the SINR and the user rates that become increasingly tight as $N$, $S$ and $K$ grow large. We will show later by simulations that these approximations are very accurate for moderate system dimensions - an observation made in [2, 7, 3, 8, 26] as well. The asymptotic analysis requires the following two assumptions.

A-2. For all $g$, $\text{lim sup}_{N}||\textbf{R}_{BS_{g}}||<\infty$ and $\text{lim sup}_{S_{g}}||\bar{\textbf{S}}_{g}||<\infty$.

A-3. For all $g$, $\text{lim inf}_{N}\frac{1}{N}\text{tr }\textbf{R}_{BS_{g}}>0$ and $\text{lim inf}_{S_{g}}\frac{1}{N}\text{tr }\bar{\textbf{S}}_{g}>0$.

We first present the Fubini theorem, which is a key mathematical idea behind our derivations.

III-A Fubini Theorem

Most works derive deterministic equivalents for random matrix models created from sums of independent random matrices. In many practical cases, like for the double scattering model, it is necessary to consider products of random matrices. To handle such models, we rely on the concept of iterative deterministic equivalents which utilizes the Fubini theorem [27, 21].

A consequence of the Fubini theorem is that if we have a function $f_{N}(\textbf{H}^{\prime}_{N},\textbf{H}^{\prime\prime}_{N})$ of two independent random sequences of matrices, $(\textbf{H}^{\prime\prime}_{N})_{N\geq 1}$ and $(\textbf{H}^{\prime}_{N})_{N\geq 1}$, then we can condition on one sequence, let’s say $(\textbf{H}^{\prime}_{N})_{N\geq 1}$, and find the deterministic equivalent $\tilde{g}_{N}(\textbf{H}^{\prime}_{N})$ for $f_{N}$. If it can be proved that this deterministic equivalent holds true for every $(\textbf{H}^{\prime}_{N})_{N\geq 1}$ generated by a space $\Omega$, then it is also valid for the random sequence $(\textbf{H}^{\prime}_{N},\textbf{H}^{\prime\prime}_{N})_{N\geq 1}$ [27]. However $\tilde{g}_{N}(\textbf{H}^{\prime}_{N})$ is still random due to dependence on $(\textbf{H}^{\prime}_{N})_{N\geq 1}$, so we need to get a deterministic equivalent $g_{N}$ for it.

This is the main mathematical idea behind our derivation of the deterministic equivalent of $\gamma_{k,g\text{RZF}}$ which is a function of the random matrices $\textbf{W}_{g}$ and the random vectors $\tilde{\textbf{w}}_{k,g}$. First, we interpret the double-scattering channel model in (4) as,

[TABLE]

where $\textbf{Z}_{g}=\frac{1}{\sqrt{S_{g}}}\textbf{R}_{BS_{g}}^{1/2}\textbf{W}_{g}\bar{\textbf{S}}_{g}^{1/2}$, and we assume $\textbf{Z}_{g}$ to be deterministic. Under this setting, the estimate of the double scattering model in (6) can be interpreted as,

[TABLE]

where $\bar{\textbf{q}}_{k,g}\sim\mathcal{CN}(\boldsymbol{0},\textbf{I}_{N})$ and $\boldsymbol{\Phi}_{g}$ is the covariance matrix of the channel estimate given as,

[TABLE]

We obtain the deterministic equivalent of $\gamma_{k,g\text{RZF}}$ in terms of certain fixed point equations that depend on the deterministic matrices $\textbf{Z}_{g}$’s using RMT results from [28]. We then extend the analysis by allowing $\textbf{Z}_{g}$’s to be random based on the Fubini theorem. In this second step, we derive the deterministic equivalents of the fixed point equations under the actual random $\textbf{Z}_{g}$’s.

To summarize, the “randomness" related to $\textbf{W}_{g}$ is first removed and then introduced later after we have the deterministic equivalent of the SINR as a function of the random vector $\tilde{\textbf{w}}_{k,g}$ and the deterministic matrix $\textbf{W}_{g}$. From this construction, the resulting deterministic approximation is referred to as an iterative deterministic equivalent.

III-B New and Useful Results

The two theorems in this section represent the major contributions of this work as they are required to cope with the non-Gaussian channel model in (4) and its estimate in (6) and form the mathematical basis of the subsequent large system analysis of RZF precoding. They provide deterministic equivalents of normalized traces and quadratic forms involving single and double occurrences of the resolvent matrix $\hat{\textbf{C}}^{-1}(\alpha)$, which is defined as $\hat{\textbf{C}}^{-1}(\alpha)=\left(\frac{1}{K}\hat{\textbf{H}}^{H}\hat{\textbf{H}}+\alpha\textbf{I}_{N}\right)^{-1}$.

The resolvent is of constant use in this paper as it arises from the expression of $\hat{\textbf{V}}$ in the RZF precoder in (10) and appears in all the terms of the SINR in (12). The trace of the resolvent is characterized by the quantities $(m_{g}(\alpha),\bar{m}_{g}(\alpha),\delta_{g}(\textbf{R}_{g},\alpha))$, $1\leq g\leq G$ in the asymptotic limit. The deterministic approximation of the SINR explicitly depends on these quantities, which is why they are introduced first in the following theorem.

Theorem 1: Consider the resolvent matrix, $\hat{\textbf{C}}^{-1}(\alpha)$, where the columns of $\hat{\textbf{H}}^{H}$ are distributed as (6). Then the following system of $3G$ implicit equations in $m_{g}(\alpha),\bar{m}_{g}(\alpha)$ and $\delta_{g}(\textbf{R}_{g},\alpha)$,

[TABLE]

has a unique solution satisfying $(m_{g}(\alpha),\bar{m}_{g}(\alpha)$, $\delta_{g}(\textbf{R}_{g},\alpha))>0$ for all $g$ and $\alpha>0$, where $\textbf{R}_{g}$ is an arbitrary matrix with a uniformly bounded spectral norm, $\bar{\textbf{D}}_{g}=\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{R}_{BS_{g}}\textbf{Q}^{H}_{g}\textbf{R}^{H}_{BS_{g}}$ and $\tilde{\textbf{D}}_{g}=\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{Q}^{H}_{g}\textbf{R}^{H}_{BS_{g}}$. Let U be any deterministic matrix with a uniformly bounded spectral norm. Under assumptions A-1, A-2, A-3 and for $\alpha>0$, we have

[TABLE]

Proof: The proof of Theorem 1 is postponed to Appendix B.

Next we present the deterministic equivalents of $\frac{1}{K}\hat{\textbf{h}}^{H}_{k,g}\hat{\textbf{C}}_{[k,g]}^{-1}\hat{\textbf{h}}_{k,g}$ and $\frac{1}{K}{\textbf{h}}^{H}_{k,g}\hat{\textbf{C}}_{[k,g]}^{-1}\hat{\textbf{h}}_{k,g}$, where $\hat{\textbf{C}}_{[k,g]}=\frac{1}{K}\hat{\textbf{H}}_{[k,g]}^{H}\hat{\textbf{H}}_{[k,g]}+\alpha\textbf{I}_{N}$, where $\hat{\textbf{H}}_{[k,g]}$ is defined in Section II-D. These quantities arise in the expression of $\gamma_{k,g\text{RZF}}$ and appear repeatedly in our analysis.

Lemma 2: Under the setting of assumptions A-1, A-2 and A-3 and for $\alpha>0$,

[TABLE]

where $m_{g}$ is obtained as the unique solution of (18) and $h_{g}=\frac{1}{K}d_{g}\left(\sum_{j=1}^{S_{g}}\frac{\bar{s}_{g,j}\delta_{g}(\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{R}_{BS_{g}},\alpha)}{1+\frac{K_{g}}{K}d_{g}^{2}\bar{s}_{g,j}\bar{m}_{g}(\alpha)\delta_{g}(\bar{\textbf{D}}_{g},\alpha)}\right)$.

Proof: The proof of Lemma 2 is provided in Appendix C.

Theorem 1 and Lemma 2 approximate quantities with one occurrence of the resolvent matrix. However, with RZF, random terms involving two resolvents arise in the expression of the SINR, a case which is out of the scope of Theorem 1. We therefore develop the result for this case.

Theorem 2: Define $\boldsymbol{\chi}(\alpha)=[\chi_{1}(\alpha),\dots,\chi_{G}(\alpha)]^{T}$, where $\chi_{g}(\alpha)$ =$\frac{1}{K^{2}}\text{tr }\boldsymbol{\Phi}_{g}\hat{\textbf{C}}^{-1}(\alpha)\hat{\textbf{H}}^{H}\textbf{P}\hat{\textbf{H}}\hat{\textbf{C}}^{-1}(\alpha)$, where $\boldsymbol{\Phi}_{g}$ is given by (17). Under the setting of assumptions A-1, A-2 and A-3 and for $\alpha>0$,

[TABLE]

where,

[TABLE]

[TABLE]

where for arbitrary matrices $\textbf{R}_{g}$ and L having uniformly bounded spectral norm,

[TABLE]

and $\textbf{u}^{\prime}(\bar{\textbf{D}},\textbf{L},\alpha)=[u^{\prime}_{1}(\bar{\textbf{D}}_{1},\textbf{L},\alpha),u^{\prime}_{2}(\bar{\textbf{D}}_{2},\textbf{L},\alpha),\dots,u^{\prime}_{G}(\bar{\textbf{D}}_{G},\textbf{L},\alpha)]^{T}$, which can be expressed as a system of linear equations as follows,

[TABLE]

for $g,i=1,\dots,G$.

Proof: The proof of Theorem 2 is provided in Appendix D.

III-C Deterministic Approximation of the SINR

Based on the results in Theorem 1, Lemma 2 and Theorem 2, the deterministic equivalent of the SINR in (12) can be derived and is presented in the next theorem.

Theorem 3: Under the setting of assumptions A-1, A-2 and A-3 and for $\alpha>0$, the downlink SINR of user $k$ in group $g$ defined in (12) converges almost surely as,

[TABLE]

where,

[TABLE]

where,

[TABLE]

Proof: The proof of Theorem 3 is given in Appendix F.

Corollary 1: Assume that A-1, A-2 and A-3 hold true and $\alpha>0$. Then the individual downlink rates $R_{k,g}$ of users converge as,

[TABLE]

where,

[TABLE]

where $\gamma^{o}_{k,g\text{RZF}}$ is given by (33).

Proof: The proof follows from the application of the continuous mapping theorem [29] to the logarithm function and the almost sure convergence of $\gamma_{k,g\text{RZF}}$ in (32).

An approximation of the average system sum rate can be obtained by replacing $R_{k,g}$ in (9) with its asymptotic approximation as,

[TABLE]

These asymptotic expressions will be shown to provide good approximations even for moderate system sizes by the means of simulations in Section IV. This means they can be used for evaluating the performance of practical systems without relying on time-consuming Monte-Carlo simulations. These expressions are also useful in performing different optimization tasks, for example, determining the optimal $\alpha$ and power allocation as shown in [7].

The deterministic equivalents depend only on the ‘slowly varying’ covariance matrices $\textbf{R}_{BS_{g}}$ and $\bar{\textbf{S}}_{g}$, instead of the instantaneous channels that vary very fast. These correlation matrices, although huge in size in the large $N,S$ regime, can be computed at the BS using knowledge of only the large-scale channel statistics. In particular, each correlation matrix, can be computed using estimates of the angular spreads, $\sigma_{t,g},\sigma_{s,g}$, and the mean AoD, $\mu_{t,g}=\mu_{s,g}$. These parameters are illustrated in Fig. 1. The spreads can be locally estimated while $\mu_{t,g}$ depends on the line of sight (LoS) angles of the users within a group, that have to be fed to the BS. However, the estimation of these large-scale parameters must be performed only once per coherence period (rather than at the same pace as the small-scale fading) and the number of required estimated parameters depends on $G$ and $K$ (much less than $N$ and $S$ in practical systems). Discussion on the estimation of large-scale parameters can be found in [30, 31], where the authors reason that the estimation of channel covariance matrices imposes a low CSI feedback overhead.

We would also point out here that the considered RZF precoding, known for its high performance, also imposes a high computational complexity in large-scale MIMO systems due to the matrix inversion operation. It is therefore interesting to study the double scattering model with other precoding schemes that reduce this prohibitively high computational complexity. One way is to employ truncated polynomial expansion (TPE) precoding proposed in [8], which replaces the matrix inversion by a TPE. The TPE precoding can perform close to RZF, with a reduced computational burden. Theorem 1 and Theorem 2 developed in this work will play a direct role in the asymptotic analysis of the SINR of the double scattering model with TPE.

Despite being useful, the expressions in this work are quite involved and do not provide direct insights into the interplay between different system parameters. We therefore focus on Rayleigh product channel and obtain some simplifications under different operating conditions.

III-D Rayleigh Product Channel

A special case of the double-scattering channel is the Rayleigh product channel which does not exhibit any form of correlation [18]. For this model, popularly known as the multi-keyhole channel, we find that Theorem 3 can be given in a closed-form as shown in the next corollary.

Corollary 2: For $G=1$, let $S_{1}=S$, $K_{1}=K$, and assume $\bar{\textbf{S}}_{1}=\textbf{I}_{S}$ and $\textbf{R}_{BS_{1}}=\textbf{I}_{N}$. Then $\gamma^{o}_{k\text{RZF}}$ in Theorem 3 can be given in a closed-form as,

[TABLE]

and $\bar{m}\in\left(0,1\right)$ is given as the unique root to,

[TABLE]

Sketch of Proof: One can show by straightforward but tedious calculations that the fundamental equations in (18) can be reduced to a single polynomial equation (III-D) in $\bar{m}$, for a single group and $\bar{\textbf{S}}_{1}=\textbf{I}_{S_{1}}$ and $\textbf{R}_{BS_{1}}=\textbf{I}_{N}$. Simplifying (120) for the setting of Corollary 2, writing it as a polynomial in $F(\textbf{R},\textbf{L})$ and applying implicit function theorem with respect to $l$ yields the expression of $u^{\prime}$. After some tedious calculations, we show that Theorem 2 can be given in a closed form and the SINR expressions in Theorem 3 can be greatly simplified.

Corollary 3: Under the setting of Corollary 2, let $\frac{S}{N},\frac{S}{K}\rightarrow\infty$. Then $\gamma^{o}_{k\text{RZF}}$ defined in Theorem 3 approaches the limit $\gamma^{o}_{k\text{RZF},\frac{S}{N},\frac{S}{K}\rightarrow\infty}$ which is given by,

[TABLE]

where $h=\frac{a}{\alpha}\left(\frac{N}{K}-1+\bar{m}\right)$, $\bar{\beta}=\frac{a^{4}\left(\frac{N}{K}-1+\bar{m}\right)(1+\rho_{tr})}{\alpha(\alpha\rho_{tr}+a^{2}\bar{m}+a^{2}\bar{m}\rho_{tr})}$, $\bar{m}=\frac{1-\frac{N}{K}-c+\sqrt{(c+\frac{N}{K}-1)^{2}+4c}}{2}$ and $c=\frac{\rho_{tr}\alpha}{a^{2}(1+\rho_{tr})}$.

Note that as $\frac{S}{N},\frac{S}{K}\rightarrow\infty$, h behaves as a Rayleigh fading channel, whose SINR is given as (45). Setting, $\rho_{tr}\rightarrow\infty$, we obtain the following corollary for the perfect CSI case.

Corollary 4: Under the setting of Corollary 3, let $\rho_{tr}\rightarrow\infty$. Then $\gamma^{o}_{k\text{RZF},\frac{S}{N},\frac{S}{K}\rightarrow\infty}$ is given by,

[TABLE]

where $\bar{m}=\frac{1-\frac{N}{K}-\alpha+\sqrt{(\alpha+\frac{N}{K}-1)^{2}+4\alpha}}{2}$.

This result has also been obtained in Corollary 2 of [7], where the authors derive the deterministic equivalent of the SINR under RZF, with the channels modeled as correlated Rayleigh.

Corollary 5: Under the setting of Corollary 2, let $\frac{N}{S}\rightarrow\infty$ and $\frac{N}{K}\rightarrow\infty$ with $K>S$. Then $\gamma^{o}_{k\text{RZF}}$ defined in Theorem 3 approaches the limit $\gamma^{o}_{k\text{RZF},\frac{N}{S},\frac{N}{K}\rightarrow\infty}$ which is given by,

[TABLE]

Corollary 6: If $K\leq S$ in the setting of Corollary 5, then $\bar{m}\rightarrow 0$ and as a consequence, the interference term $\Upsilon_{k}^{o}\rightarrow 0$ and the power normalization term $\xi^{o}(\textbf{I}_{N},\textbf{I}_{N})\rightarrow 0$ in Theorem 3.

Remark 1: Corollary 6 implies for $\frac{N}{S},\frac{N}{K}\rightarrow\infty$ with $K\leq S$, $\gamma^{o}_{k\text{RZF},\frac{N}{S},\frac{N}{K}\rightarrow\infty}$ grows unboundedly.

Two very important observations are in order now. First, the performance of a large-scale MIMO system is limited by the number of scatterers in the propagation environment. Recent works that show the spectral efficiency of large-scale MIMO systems to be unlimited base their analysis on correlated Rayleigh channels and neglect the impact of limited scattering, that may exist in practice. Even though some of these works realize that poor scattering will deteriorate the performance but they do not account for it in the theoretical analysis. This is the first paper that derives a bound on the asymptotic ergodic rate that can be achieved under limited scattering.

The second observation is that the maximum number of users that can be served simultaneously while achieving large-scale MIMO gains is less than or equal to the number of scatterers. In fact, we find through Corollary 6 that if $K\leq S$, $\gamma^{o}_{k\text{RZF},\frac{N}{S},\frac{N}{K}\rightarrow\infty}$ grows unboundedly with $N$. This ‘large-scale MIMO gain’ is in accordance with the results in [2] and [1], which when simplified for a single cell case, reveal that $\gamma_{k\text{RZF}}\rightarrow\infty$ for $N\rightarrow\infty$, $K/N\rightarrow 0$. The inherent assumption in these works that rely on the correlated Rayleigh channel is that there is enough scattering to serve all the $K$ users. However, as $K$ exceeds $S$, we start to see a significant performance loss. For example, with uniform power allocation, serving $K=S+1$ users would result in the SINR of each user approaching $S$ as $N\rightarrow\infty$, according to Corollary 5. If we serve $K=S+2$ users, the SINR will approach $S/2$, resulting in a $3\rm{dB}$ loss. Thus, the value of $S$ is what determines the maximum number of users that can be served while realizing large-scale MIMO gains.

IV Simulation Results

Under the double scattering model, the correlation matrices $\textbf{R}_{BS_{g}}$ and $\bar{\textbf{S}}_{g}$ are given as $\textbf{R}_{BS_{g}}=\textbf{G}(\mu_{t,g},\sigma_{t,g},d_{t},S_{g})$ and $\bar{\textbf{S}}_{g}=\textbf{G}(\mu_{s,g},\sigma_{s,g},d_{s,g},S_{g})$, where $\textbf{G}(\mu,\sigma,d,n)$ is defined as [21],

[TABLE]

All parameters have already been defined in Fig. 1. The parameter values are set as $G=4$, $K=128$, $N=128$, $S_{g}$=$\{130,140,134,144\}$, $\mu_{t,g}=\mu_{s,g}$=$\{-\pi/3,-\pi/9,\pi/9,\pi/3\}$, $\sigma_{t,g}$=$\{\pi/5,\pi/6,\\ \pi/5,\pi/7\}$ and $\sigma_{s,g}$=$\{\pi/6,\pi/6,\pi/6,\pi/6\}$. Also, $d_{t}=0.5$ and $d_{s,g}=2$ for all $g$. We assume an equal number of users in each group, i.e. $K_{g}=K/G$, with uniform power allocation, $\textbf{P}=\textbf{I}_{K}$. Fig. 3 compares the downlink system sum rate $R_{sum}=\sum_{g=1}^{G}\sum_{k=1}^{K_{g}}\log(1+\gamma_{{k,g}_{\text{RZF}}})$ obtained using 2000 Monte-Carlo realizations of the SINR in (12) to the asymptotic approximation provided in (39), where $\gamma^{o}_{k,g\text{RZF}}$ is given by (33). It can be seen that the asymptotic result yields a very good approximation for moderate system dimensions.

To explain the decrease in sum-rate at high values of $\rho$ for $\rho_{tr}=2\rm{dB}$, note that the goal of the regularization parameter $\alpha$ in the RZF precoder in (10) is to improve the conditioning of the random matrix $\hat{\textbf{H}}^{H}\hat{\textbf{H}}+K\alpha\textbf{I}_{N}$ and hence it should depend on the channel parameters, including the training SNR $\rho_{tr}$ and the downlink SNR $\rho$. Finding the optimal $\alpha$ that maximizes $\gamma^{o}_{k,g\text{RZF}}$ is outside the scope of this work. Now for the chosen $\alpha=1/\rho$ (optimal under perfect CSI), the sum rate is decreasing at high SNR values for $\rho_{tr}=2\rm{dB}$. This is because $\alpha$ does not account for $\rho_{tr}$ and thus the matrix $\hat{\textbf{H}}^{H}\hat{\textbf{H}}+K\alpha\textbf{I}_{N}$ in the RZF precoder becomes ill-conditioned as the quality of the estimate deteriorates and $\alpha$ decreases (due to increase in $\rho$). This results in a loss in the performance. This effect has also been observed in ([7], Fig. 2). In fact, the authors there obtain the optimal $\alpha^{*}$ as a solution of a fixed point equation that depends on $\rho_{tr}$ and show that for $\alpha=\alpha^{*}$, the sum rate continues to increase with $\rho$ for any quality of channel estimate.

Furthermore, note that the mismatch between the theoretical and the Monte-Carlo results starts to increase for high SNR values due to the slower convergence of $\gamma_{k,g\text{RZF}}$ to its deterministic approximation as well documented in RMT literature [7, 32]. A better approximation at high SNR can be seen in systems with higher values of $N,K$ and $S$.

Next we plot in Fig. 3 the average achievable rate for $G=11$, $N=121$, $S_{g}=120$ $\forall g$, and $K=60$ users having mean AoDs generated between $-15^{o}$ to $15^{o}$ and spreads between $3^{o}$ and $6^{o}$. We employ the fixed quantization method from [31] to group the users, based on the criteria of the minimum chordal distance between the users’ correlation eigenspaces and the group subspaces. The Monte-Carlo and deterministic average rates when common correlation matrices are assumed are shown to be close. In fact, the average (and sum) rate now continues to increase for $\rho_{tr}=2\rm{dB}$ because for $K<N$, $\hat{\textbf{H}}^{H}\hat{\textbf{H}}+K\alpha\textbf{I}_{N}$ is well-conditioned. The Monte-Carlo average rate under the actual different correlation matrices of the users is plotted in black to highlight the importance of having variation in correlation matrices. However, for us to analyze this difficult non-Gaussian model, it was necessary to assume common correlation matrices within groups. The performance loss due to this assumption can be reduced by increasing the number of groups.

Fig. 5 studies the effect of the number of scatterers on the system sum rate for a single group with multi-keyhole channel, i.e. $G=1$, $\textbf{R}_{BS}=\textbf{I}_{N}$, $\bar{\textbf{S}}=\textbf{I}_{S}$, with $N=128$, $K=32$, $\alpha=1/\rho$ and $\rho_{tr}=10\rm{dB}$. The downlink sum rate in (39) plotted using the closed-form expression of $\gamma^{o}_{k\text{RZF}}$ in Corollary 2 is close to the Monte-Carlo result even for a very low number of scatterers. The spatial multiplexing gains are seen to increase linearly with $S$. However, for $S>N$, the gains start to decrease since the degrees of freedom are limited by the number of antennas at the BS. The limiting sum rate as $S/N,S/K\rightarrow\infty$ is also plotted using the SINR in (45). As the number of scatterers increases, the performance approaches to that of a Rayleigh fading channel.

Finally, we study the performance of the Rayleigh product channel in Fig. 5 as the number of BS antennas increases for a single group. It can be seen that the performance for the case when $S=30$, $K=40$ saturates to a limiting sum rate, given by substituting (47) in (39), which is also plotted in black on the figure. This result confirms that it is useless to deploy more antennas when the number of scatterers in the environment is limited. The other curves highlight that when $K\leq S$, we see the ‘large-scale MIMO effect’ and the performance continues to grow with $N$ because the channel has enough degrees of freedom to support all the users.

V Conclusion

In this paper, we studied a large-scale multi-user MISO system with double-scattering channels that are more realistic than the commonly used Gaussian channels. We first derived the MMSE estimate for this channel. Then under the assumption of per-group channel correlation matrices, we derived the deterministic approximations of the SINR and ergodic rates with RZF precoding, that are tight in the large system limit. Simulation results showed a close match between the asymptotic and the Monte-Carlo simulated sum rate and provided insights into the performance of multi-keyhole channels. We showed that large-scale MIMO gains can only be realized when the number of scheduled users is less than the number of scatterers in the environment.

Appendix A Proof of Lemma 1

The MMSE estimate $\hat{\textbf{h}}_{k,g}$ of $\textbf{h}_{k,g}$ is computed as $\hat{\textbf{h}}_{k,g}=\textbf{F}_{g}\textbf{y}^{tr}_{k,g}$ where $\textbf{F}_{g}$ is obtained as the solution of $\text{min }_{\textbf{F}_{g}}\mathbb{E}[||\textbf{h}_{k,g}-\textbf{F}_{g}\textbf{y}^{tr}_{k,g}||^{2}]$ resulting in,

[TABLE]

Using the expression of $\textbf{y}_{k,g}^{tr}$ in (5) and the independence of $\textbf{h}_{k,g}$ and $\textbf{n}_{k,g}$ we get,

[TABLE]

Conditioning the expectation on $\textbf{W}_{g}$ first, we obtain,

[TABLE]

Therefore $\hat{\textbf{h}}_{k,g}=\frac{1}{S_{g}}(\text{tr }\bar{\textbf{S}}_{g})\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{y}^{tr}_{k,g}$, where $\textbf{Q}_{g}=\left(\frac{1}{S_{g}}(\text{tr }\bar{\textbf{S}}_{g})\textbf{R}_{BS_{g}}+\frac{1}{\rho_{tr}}\textbf{I}_{N}\right)^{-1}$.

Appendix B Proof of Theorem 1

As a starting point we assume $\textbf{Z}_{g}$ to be deterministic such that $\text{lim sup}_{N}||\textbf{Z}_{g}\textbf{Z}_{g}^{H}||<\infty$. Under this setting, we use the expression of $\hat{\textbf{h}}_{k,g}$ given by (16) in the resolvent matrix, $\hat{\textbf{C}}^{-1}$, to have,

[TABLE]

where the entries of $\bar{\textbf{Q}}_{g}\in\mathbb{C}^{N\times K_{g}}\sim\mathcal{CN}(0,\frac{1}{K_{g}})$ and $\tilde{\boldsymbol{\Phi}}_{g}=\frac{K_{g}}{K}\boldsymbol{\Phi}_{g}$, where $\boldsymbol{\Phi}_{g}$ is given by (17). Note that $\tilde{\boldsymbol{\Phi}}^{1/2}_{g}\bar{\textbf{Q}}_{g}\textbf{I}_{K_{g}}\bar{\textbf{Q}}_{g}^{H}\tilde{\boldsymbol{\Phi}}^{1/2^{H}}_{g}$ is a double scattering model with deterministic $\tilde{\boldsymbol{\Phi}}^{1/2}_{g}$s. Under this setting, the deterministic equivalent of $\frac{1}{K}\text{tr }\textbf{U}\hat{\textbf{C}}^{-1}$ is obtained using Corollary 1 from [28] as,

[TABLE]

where $(e_{g},\bar{e}_{g})$ are given as a unique solution to the following set of implicit equations,

[TABLE]

[TABLE]

Now for the actual double-scattering channel model, $\textbf{Z}_{g}$s are random and modeled as

$\frac{1}{\sqrt{S_{g}}}\textbf{R}_{BS_{g}}^{1/2}\textbf{W}_{g}\bar{\textbf{S}}_{g}^{1/2}$. Using the Fubini Theorem [27] we can extend the result in (55) to random $\textbf{Z}_{g}$s by finding deterministic equivalents, denoted as $m_{g}$ and $\bar{m}_{g}$, for $e_{g}$ and $\bar{e}_{g}$ respectively.

To do this, we first define quantities $e_{g,i,j}$ and $\bar{e}_{g,i,j}$ for $i=1,\dots,G$, $j=1,\dots,S_{i}$ as $\bar{e}_{g,i,j}=\frac{1}{1+e_{g,i,j}}$ and $e_{g,i,j}=\frac{1}{K_{g}}\text{tr }\tilde{\boldsymbol{\Phi}}_{g,i,j}(\sum_{l=1}^{G}\bar{e}_{l,i,j}\tilde{\boldsymbol{\Phi}}_{l,i,j}+\alpha\textbf{I}_{N})^{-1}$, where,

[TABLE]

It can be shown following the techniques used in Appendix E of [27] that for all $g,i,j$,

[TABLE]

Now using the expression of $\boldsymbol{\Phi}_{g}$ from (17) in the expression of $e_{g}$ in (57) we have,

[TABLE]

[TABLE]

where $\tilde{\textbf{D}}_{g}=\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{Q}^{H}_{g}\textbf{R}^{H}_{BS_{g}}$. First we treat $T_{1}$ using (61) and matrix inversion lemma (Lemma 1 in [2]) to remove the dependence of $\left(\sum_{i=1}^{G}\bar{e}_{i}\tilde{\boldsymbol{\Phi}}_{i}+\alpha\textbf{I}_{N}\right)^{-1}$ on the vector $\textbf{z}_{g,j}$ as,

[TABLE]

Then from trace lemma and rank-one perturbation lemma (Lemmas 3 and 5 in [33]) we have,

[TABLE]

where $\bar{\textbf{D}}_{g}=\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{R}_{BS_{g}}\textbf{Q}^{H}_{g}\textbf{R}^{H}_{BS_{g}}$.

Notice that $\tilde{\boldsymbol{\Phi}}_{i}$ is still a function of $\textbf{Z}_{i}$ which is random. To remove the dependence of $T_{1}$ on $\textbf{Z}_{i}$, we need the deterministic equivalent of $\frac{1}{S_{g}}\text{tr }\bar{\textbf{D}}_{g}\left(\sum_{i=1}^{G}\bar{e}_{i}\tilde{\boldsymbol{\Phi}}_{i}+\alpha\textbf{I}_{N}\right)^{-1}$. Using the expression of $\tilde{\boldsymbol{\Phi}}_{i}$ in (17) and Corollary 1 from [28], we have the following convergence,

[TABLE]

where ${\textbf{R}}_{g}$ is any deterministic matrix that satisfies $\text{lim sup}_{N}||\textbf{R}_{g}||<\infty$ and,

[TABLE]

such that $(f_{g}({\textbf{R}}_{g}),\bar{f}_{g})\geq 0$. Substituting $\bar{f}_{g}$ in $f_{g}({\textbf{R}}_{g})$, we have,

[TABLE]

Consequently we get the following convergence for $T_{1}$ and $T_{2}$,

[TABLE]

Plugging $T_{1}$ and $T_{2}$ into the expression of $e_{g}$ in (63) yields,

[TABLE]

where $\epsilon_{g}\xrightarrow[N\rightarrow\infty]{a.s.}0$. Consider the deterministic counterpart of $(e_{g}(\alpha),\bar{e}_{g}(\alpha),f_{g}(\textbf{R}_{g},\alpha))$ as $(m_{g}(\alpha),\\ \bar{m}_{g}(\alpha),\delta_{g}(\textbf{R}_{g},\alpha))$ defined in (18) and define $\Upsilon_{1}=\text{max}_{g}|e_{g}(\alpha)-m_{g}(\alpha)|$, $\Upsilon_{2}=\text{max }_{g}|\bar{e}_{g}(\alpha)-\bar{m}_{g}(\alpha)|$, $\Upsilon_{3}=\text{max}_{g}|f_{g}({\textbf{R}}_{g},\alpha)-\delta_{g}({\textbf{R}}_{g},\alpha)|$ and $\epsilon=\text{max}_{g}|\epsilon_{g}|$. It can be shown using the techniques from Appendix E of [27] that,

[TABLE]

for $\alpha$ sufficiently large and $\epsilon\xrightarrow[N\rightarrow\infty]{a.s.}0$. This result can be extended to all $\alpha$ by Vitali convergence theorem [34].

Now combining (55) with the result in (70) for $\textbf{R}_{g}=\textbf{U}$ will yield,

[TABLE]

where $\bar{\textbf{T}}(\alpha)$ is given by (19). This completes the proof of Theorem 1. The uniqueness of the solution $(m_{g}(\alpha),\bar{m}_{g}(\alpha),\delta_{g}(\textbf{R}_{g},\alpha))$ can be proved by showing that the G-variate function in (18) is a standard interference function as done for Theorem 2 in [35].

Appendix C Proof of Lemma 2

C-A * Deterministic equivalent of $\frac{1}{K}\hat{\textbf{h}}_{k,g}^{H}\hat{\textbf{C}}_{[k,g]}^{-1}\hat{\textbf{h}}_{k,g}$:*

It has been shown in ([27], Appendix E) that $\textbf{Z}_{g}$ in (15) satisfies $\text{lim sup}_{N}||\textbf{Z}_{g}\textbf{Z}_{g}^{H}||<\infty$ almost surely. We can therefore assume $\textbf{Z}_{g}$ to be deterministic as a starting point and use trace lemma and rank-one lemma (Lemma 3 and Lemma 5 from [33]) to have,

[TABLE]

where $\boldsymbol{\Phi}_{g}$ is given by (17). Utilizing the result in (55) for $\textbf{U}=\boldsymbol{\Phi}_{g}$ yields $\frac{1}{K}\text{tr }\boldsymbol{\Phi}_{g}\hat{\textbf{C}}^{-1}-e_{g}\xrightarrow[N\rightarrow\infty]{a.s.}0$. The result can be extended to random $\textbf{Z}_{g}$s using Fubini Theorem and the deterministic equivalent of $e_{g}$ is obtained as $m_{g}$ in Appendix B.

C-B * Deterministic equivalent of $\frac{1}{K}{\textbf{h}}_{k,g}^{H}\hat{\textbf{C}}_{[k,g]}^{-1}\hat{\textbf{h}}_{k,g}$:*

Utilizing the expressions of the double scattering channel in (15) and its MMSE estimate in (6) and the independence of $\tilde{\textbf{w}}_{k,g}$ and $\textbf{n}_{k,g}^{tr}$, we have,

[TABLE]

Next exploit trace lemma and rank-1 lemma (Lemmas 3 and 5 from [33]) to get,

[TABLE]

where $\hat{\textbf{C}}$ is given by (54) under the assumption of deterministic $\textbf{Z}_{g}$ and $\tilde{\bar{\boldsymbol{\Phi}}}_{g}=\frac{K_{g}}{K}\bar{\boldsymbol{\Phi}}_{g}$, where $\bar{\boldsymbol{\Phi}}_{g}=d_{g}\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{Z}_{g}\textbf{Z}_{g}^{H}$. Utilizing the result in (55) for $\textbf{U}=\bar{\boldsymbol{\Phi}}_{g}$, we have,

[TABLE]

where $(e_{g},\bar{e}_{g})$ are given as a unique solution to the set of implicit equations in (56), (57), Denote $\frac{1}{K_{g}}\text{tr }\tilde{\bar{\boldsymbol{\Phi}}}_{g}(\sum_{i=1}^{G}\bar{e}_{i}\tilde{\boldsymbol{\Phi}}_{i}+\alpha\textbf{I}_{N})^{-1}$ in (82) as $T_{3}$ and let $\textbf{Z}_{g}$’s be random now. A similar development as done for the term $T_{1}$ in Appendix B will yield,

[TABLE]

Using the convergence in (70) and $f_{g}(\textbf{R}_{g})$ defined in (72), we obtain,

[TABLE]

Replace $\bar{e}_{g}$ and $f_{g}(\textbf{R}_{g})$ in (84) by their deterministic approximations $\bar{m}_{g}$ and $\delta_{g}(\textbf{R}_{g})$ respectively from Appendix B. This yields the deterministic equivalent of $\frac{1}{K}\textbf{h}_{k,g}^{H}\hat{\textbf{C}}_{[k,g]}^{-1}\hat{\textbf{h}}_{k,g}$ denoted as $h_{g}$.

Appendix D Proof of Theorem 2

We are interested in finding the deterministic equivalent of:

[TABLE]

First we need to control the variance of ${\chi}_{g}$ and prove that $\text{var}(\chi_{g})$ converges to zero. This can be done using standard tools from RMT (see [26]) and will imply that,

[TABLE]

where $\chi^{o}_{g}=\mathbb{E}[\chi_{g}]$. This allows us to focus directly on $\mathbb{E}[\chi_{g}]$.

Using the expression of $\boldsymbol{\Phi}_{g}$ in (17) we have,

[TABLE]

where,

[TABLE]

and $\kappa^{o}_{g}$ and $\xi^{o}$ denote their respective expectations. We therefore start the proof by deriving the deterministic equivalents of $\kappa_{g}(\textbf{A},\textbf{B})$ and $\xi(\textbf{A},\textbf{B})$, where A and B are deterministic matrices of uniformly bounded spectral norm. Plugging these results in (88) will complete the proof of Theorem 2.

D-A * Deterministic equivalent of $\kappa_{g}(\textbf{A},\textbf{B})=\frac{1}{K^{2}}\text{tr }\textbf{A}\textbf{Z}_{g}\textbf{Z}_{g}^{H}\textbf{B}\hat{\textbf{C}}^{-1}\hat{\textbf{H}}^{H}\textbf{P}\hat{\textbf{H}}\hat{\textbf{C}}^{-1}$:*

The aim of this section is to derive a deterministic equivalent for the random quantity $\kappa_{g}(\textbf{A},\textbf{B})$. We can again control the variance of $\kappa_{g}$ and show that it goes to zero implying that,

[TABLE]

where $\kappa^{o}_{g}(\textbf{A},\textbf{B})=\mathbb{E}[\kappa_{g}(\textbf{A},\textbf{B})]$. Using the result from last section in (55) that $\frac{1}{K}\text{tr }\hat{\textbf{C}}^{-1}-\frac{1}{K}\text{tr }\textbf{T}\xrightarrow[N\rightarrow\infty]{a.s.}0$, where $\textbf{T}=(\sum_{i=1}^{G}\frac{\tilde{\boldsymbol{\Phi}}_{i}}{1+m_{i}}+\alpha\textbf{I}_{N})^{-1}$ and the resolvent identity given as $\hat{\textbf{C}}^{-1}-\textbf{T}=\textbf{T}(\textbf{T}^{-1}-\hat{\textbf{C}})\hat{\textbf{C}}^{-1}=\textbf{T}\left(\sum_{i=1}^{G}\frac{\tilde{\boldsymbol{\Phi}}_{i}}{1+m_{i}}-\frac{1}{K}\hat{\textbf{H}}_{i}^{H}\hat{\textbf{H}}_{i}\right)\hat{\textbf{C}}^{-1}$, we decompose $\kappa_{g}$ as,

[TABLE]

We start by applying matrix inversion lemma (Lemma 1 in [2]) on $\mathbb{E}[Z_{1}]$ and re-writing the expression to get,

[TABLE]

Assuming $\textbf{Z}_{g}$ to be deterministic and using the fact that $\text{lim sup}_{N}||\textbf{Z}_{g}\textbf{Z}_{g}^{H}||<\infty$ almost surely, we can use the expression of $\hat{\textbf{h}}_{k,g}$ in (16) and show with the help of trace lemma (Lemma 3 from [33]) that the first term on the right side of the above equation is negligible. Now using the rank-one perturbation lemma (Lemma 5 from [33]) and Lemma 2 on the second term yields,

[TABLE]

Using trace lemma (Lemma 3 from [33]) with the expression of $\hat{\textbf{h}}_{k,g}$ in (16) we have,

[TABLE]

Now using the rank-one lemma (Lemma 5 from [33]) and Corollary 1 from [28] we have,

[TABLE]

Note that the analysis in the last two steps can be extended to random $\textbf{Z}_{g}$s using Fubini Theorem and we need the deterministic equivalent of $\bar{\beta}_{g,i}(\textbf{A},\textbf{B})=\frac{1}{K_{i}}\mathbb{E}[\text{tr }\tilde{\boldsymbol{\Phi}}_{i}\textbf{T}\textbf{A}\textbf{Z}_{g}\textbf{Z}_{g}^{H}\textbf{B}\textbf{T}]$ under the actual random $\textbf{Z}_{g}$s. Using the expression of $\boldsymbol{\Phi}_{i}$ from (17) results in $\bar{\beta}_{g,i}(\textbf{A},\textbf{B})=$

[TABLE]

We analyze the two terms separately. First we use matrix inversion lemma (Lemma 1 from [2]) on $T_{4}$ to remove the dependence of T on $\textbf{z}_{i,j}$ resulting in $T_{4}=$

[TABLE]

where \tilde{\textbf{T}}_{i}=\Big{(}\sum_{l=1}^{G}\frac{K_{l}}{K}d_{l}^{2}\bar{m}_{l}\textbf{R}_{BS_{l}}\textbf{Q}_{l}\left(\textbf{Z}_{l}\textbf{Z}_{l}^{H}+\frac{1}{\rho_{tr}}\textbf{I}_{N}\right)\textbf{Q}^{H}_{l}\textbf{R}^{H}_{BS_{l}}-\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}\textbf{R}_{BS_{i}}\textbf{Q}_{i}\textbf{z}_{i,j}\textbf{z}_{i,j}^{H}\textbf{Q}^{H}_{i}\textbf{R}^{H}_{BS_{i}}\\ +\alpha\textbf{I}_{N}\Big{)}^{-1}. Next we use trace lemma and rank-one lemma (Lemmas 3 and 5 from [33]) as,

[TABLE]

Using matrix inversion lemma to remove the dependence of T on $\textbf{z}_{g,j}$ and using $\frac{1}{S_{i}}\text{tr }\bar{\textbf{D}}_{i}\textbf{T}-\delta_{i}(\bar{\textbf{D}}_{i})\xrightarrow[N\rightarrow\infty]{a.s.}0$ from Appendix B, where $\bar{\textbf{D}}_{i}=\textbf{R}_{BS_{i}}\textbf{Q}_{i}\textbf{R}_{BS_{i}}\textbf{Q}^{H}_{i}\textbf{R}^{H}_{BS_{i}}$, yields,

[TABLE]

Finally we use trace lemma and rank-one perturbation lemma from [33] to get,

[TABLE]

Similar steps would yield the following deterministic equivalent for $T_{5}$,

[TABLE]

where $\tilde{\textbf{D}}_{i}=\textbf{R}_{BS_{i}}\textbf{Q}_{i}\textbf{Q}^{H}_{i}\textbf{R}^{H}_{BS_{i}}$.

In order to complete the calculation of the deterministic equivalent of $\bar{\beta}_{g,i}(\textbf{A},\textbf{B})$, we need the deterministic equivalent of $u^{\prime}_{g}(\textbf{R}_{g},\textbf{L})=\mathbb{E}[\frac{1}{S_{g}}\text{tr }\textbf{R}_{g}\textbf{T}\textbf{L}\textbf{T}]$ which is stated in the following Lemma.

Lemma 3: Define $u^{\prime}_{g}(\textbf{R}_{g},\textbf{L})=\mathbb{E}[\frac{1}{S_{g}}\text{tr }\textbf{R}_{g}\textbf{T}\textbf{L}\textbf{T}]$ where $\textbf{R}_{g}$ and L, $g=1,\dots,G$ are deterministic matrices with uniformly bounded spectral norm. Then under the setting of assumptions A-1, A-2 and A-3 and for $\alpha>0$,

[TABLE]

where $\textbf{W}_{i}=\left(\textbf{I}_{S_{i}}+\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}\delta_{i}(\bar{\textbf{D}}_{i})\bar{\textbf{S}}_{i}\right)^{-1}$ and $\textbf{u}^{\prime}(\bar{\textbf{D}},\textbf{L})=[u^{\prime}_{1}(\bar{\textbf{D}}_{1},\textbf{L}),u^{\prime}_{2}(\bar{\textbf{D}}_{2},\textbf{L}),\dots,u^{\prime}_{G}(\bar{\textbf{D}}_{G},\textbf{L})]^{T}$, which can be expressed as a system of linear equations as follows:

[TABLE]

for $g,i=1,\dots,G$. The proof of Lemma 3 can be found in Appendix E.

Using Lemma 3, we have the deterministic equivalents of $T_{4}$ and $T_{5}$ (denoted as $\bar{\beta}^{1}_{g,i}(\textbf{A},\textbf{B})$ and $\bar{\beta}^{2}_{g,i}(\textbf{A},\textbf{B}$) in Theorem 2) and hence $\bar{\beta}_{g,i}(\textbf{A},\textbf{B})$. This yields the expression of $\mathbb{E}[Z_{1}]$ as,

[TABLE]

We now study $Z_{3}$. Using matrix inversion lemma and common inverses of resolvents we have,

[TABLE]

We sequentially deal with terms $Z_{31}$ and $Z_{32}$. Using Lemma 2 and some manipulation we have,

[TABLE]

Assuming $\textbf{Z}_{g}$ to be deterministic, we can use Lemmas 6, 3 and 5 from [33] and Corollary 1 from [28] sequentially to obtain,

[TABLE]

Extending the analysis to random $\textbf{Z}_{g}$ based on the Fubini Theorem and using Lemma 2 we have,

[TABLE]

The term $\Psi_{1}$ is compensated by $Z_{2}$. To see this, observe that the first order of the term does not change if we substitute $\hat{\textbf{H}}_{[l,i]}$ by $\hat{\textbf{H}}$ and $\textbf{P}_{[l,i]}$ by P and then apply trace lemma, rank-one perturbation lemma (Lemmas 3 and 5 from [33]) and Fubini theorem.

Finally, we deal with $Z_{32}$. Using Lemma 2 and some manipulation we have,

[TABLE]

Analogously to before, $\mathbb{E}[Z_{32}]$ can be simplified as,

[TABLE]

[TABLE]

where $\chi^{o}_{i}$ is defined in (88). Combining (103), (108) and (109), we obtain,

[TABLE]

D-B * Deterministic equivalents of $\xi(\textbf{A},\textbf{B})$ and $\chi_{g}$:*

Using similar argument for $\xi(\textbf{A},\textbf{B})=\frac{1}{K^{2}}\text{tr }\textbf{A}\textbf{B}\hat{\textbf{C}}^{-1}\hat{\textbf{H}}^{H}\textbf{P}\hat{\textbf{H}}\hat{\textbf{C}}^{-1}$ as done for $\kappa_{g}(\textbf{A},\textbf{B})$ in (91), we have $\xi(\textbf{A},\textbf{B})-\xi^{o}(\textbf{A},\textbf{B})\xrightarrow[N\rightarrow\infty]{a.s.}0$, where $\xi^{o}(\textbf{A},\textbf{B})=\mathbb{E}[\xi(\textbf{A},\textbf{B})]$.

Repeating the same steps as done for $\kappa^{o}_{g}(\textbf{A},\textbf{B})$, we obtain,

[TABLE]

where $\tilde{\beta}_{i}(\textbf{A},\textbf{B})=\frac{d_{i}^{2}}{K}\left(\sum_{j=1}^{S_{i}}\frac{\bar{s}_{i,j}u^{\prime}_{i}(\bar{\textbf{D}}_{i},\textbf{A}\textbf{B})}{(1+\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}\bar{s}_{i,j}\delta_{i}(\bar{\textbf{D}}_{i}))^{2}}+\frac{S_{i}}{\rho_{tr}}u^{\prime}_{i}(\tilde{\textbf{D}}_{i},\textbf{A}\textbf{B})\right)$, where $u^{\prime}_{i}(\textbf{R}_{i},\textbf{L})$ has been defined in Lemma 3.

Plugging (110) and (111) into (88) yields the deterministic equivalent of $\chi^{o}_{g}$ as,

[TABLE]

Now $\boldsymbol{\chi}^{o}=[\chi^{o}_{1},\chi^{o}_{2},\dots,\chi^{o}_{g}]^{T}$ can be expressed as a system of linear equations as,

[TABLE]

Also, plugging (112) into (110) and (111) completes the deterministic equivalents of $\kappa^{o}_{g}(\textbf{A},\textbf{B})$ and $\xi^{o}(\textbf{A},\textbf{B})$ respectively. These terms will later be needed in the proof of Theorem 3.

Appendix E Proof of Lemma 3

To derive the deterministic equivalent of $u^{\prime}_{g}(\textbf{R}_{g},\textbf{L})=\mathbb{E}[\frac{1}{S_{g}}\text{tr }\textbf{R}_{g}\textbf{T}\textbf{L}\textbf{T}]$, where $\textbf{R}_{g}$ and L, $g=1,\dots,G$ are deterministic matrices, we use the results from Appendix B in (70) and (73) along with $||f_{g}(\textbf{R}_{g})-\delta_{g}(\textbf{R}_{g})||\xrightarrow[N\rightarrow\infty]{a.s.}$ to have $\frac{1}{S_{g}}\text{tr }\textbf{R}_{g}\textbf{T}-\delta_{g}(\textbf{R}_{g})\xrightarrow[N\rightarrow\infty]{a.s.}0$, where $\textbf{T}=(\sum_{i=1}^{G}\bar{m}_{i}\tilde{\boldsymbol{\Phi}}_{i}+\alpha\textbf{I}_{N})^{-1}$ and $\delta_{g}(\textbf{R}_{g})=\frac{1}{S_{g}}\textbf{R}_{g}\bar{\textbf{T}}$, where $\bar{\textbf{T}}$ is given by (19). To this end note that,

[TABLE]

The deterministic equivalent of $\frac{1}{S_{g}}\text{tr }\textbf{R}_{g}(\textbf{T}^{-1}-l\textbf{L})^{-1}$ is obtained using Corollary 1 from [28] as,

[TABLE]

where $F_{g}({\textbf{R}}_{g},\textbf{L}),\bar{F}_{g}$ are defined as the unique solution to,

[TABLE]

such that $(F_{g}({\textbf{R}}_{g},\textbf{L}),\bar{F}_{g})\geq 0$. Substituting $\bar{F}_{g}$ in $F_{g}({\textbf{R}}_{g},\textbf{L})$, we have,

[TABLE]

Note that $F_{g}({\textbf{R}}_{g},\textbf{L})$ reduces to $\delta_{g}(\textbf{R}_{g})$ for $l=0$.

The expression of $u^{\prime}_{g}({\textbf{R}}_{g},\textbf{L})$ can be obtained by using (116) and (E) to get $u^{\prime}_{g}({\textbf{R}}_{g},\textbf{L})=\frac{d}{dl}F_{g}({\textbf{R}}_{g},\textbf{L})|_{l=0}+o(1)$ resulting in,

[TABLE]

where $\textbf{W}_{i}=\left(\textbf{I}_{S_{i}}+\frac{K_{i}}{K}d_{i}^{2}\bar{m}_{i}\delta_{i}(\bar{\textbf{D}}_{i})\bar{\textbf{S}}_{i}\right)^{-1}$.

Note that $u^{\prime}_{g}(\textbf{R}_{g},\textbf{L})$ depends on the values of $u^{\prime}_{z}(\bar{\textbf{D}}_{z},\textbf{L})$. The latter can be expressed as a system of linear equations by solving (121) for $\textbf{R}_{g}=\bar{\textbf{D}}_{g}$. This system is represented by (100).

Appendix F Proof of Theorem 3

The deterministic equivalents of the energy term $|\mathbb{E}[\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g}]|^{2}$, the term $\Theta$ of the power normalization, the interference term $\mathbb{E}[\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{V}}\textbf{h}_{k,g}]$ and the variance term

$\text{var }(\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g})$ are worked out separately to yield the deterministic equivalent of the SINR.

F-A * Deterministic equivalent of $|\mathbb{E}[\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g}]|^{2}$:*

Note that $\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g}$ can be written as $\frac{1}{K}\textbf{h}_{k,g}^{H}\hat{\textbf{C}}^{-1}\hat{\textbf{h}}_{k,g}$, where $\hat{\textbf{C}}=\frac{1}{K}\hat{\textbf{H}}^{H}\hat{\textbf{H}}+\alpha\textbf{I}_{N}$. In order to remove the dependency of $\hat{\textbf{C}}$ on $\hat{\textbf{h}}_{k,g}$, we use common inverses of resolvents lemma (Lemma 2 from [33]) to get,

[TABLE]

The deterministic equivalents of both $\frac{1}{K}\hat{\textbf{h}}_{k,g}^{H}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{h}}_{k,g}$ and $\frac{1}{K}{\textbf{h}}_{k,g}^{H}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{h}}_{k,g}$ have been derived in Lemma 2 as $m_{g}$ and $h_{g}$ respectively. Note that Lemma 2 not only implies almost sure convergence but also convergence in mean. Therefore, by the dominated convergence theorem and the continuous mapping theorem, we have

[TABLE]

F-B *Deterministic equivalent of $\Theta=\mathbb{E}[\text{tr }\textbf{P}\hat{\textbf{H}}\hat{\textbf{V}}^{2}\hat{\textbf{H}}^{H}]$: *

[TABLE]

where $\xi^{o}(\textbf{A},\textbf{B})$ was derived in Appendix D and is given by (111).

F-C Deterministic equivalent of $\mathbb{E}[\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{V}}\textbf{h}_{k,g}]$:

Denote $\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{V}}\textbf{h}_{k,g}$ as $\Upsilon_{k,g}$ and use Lemma 2 from [33] to decompose it as,

[TABLE]

Let us begin by treating $X_{1}$. Using trace lemma and results from Appendix D we have,

[TABLE]

Therefore,

[TABLE]

where $\kappa^{o}_{g}(\textbf{A},\textbf{B})$ was derived in the last section and is given by (110).

Next note that,

[TABLE]

Using the definition of $\textbf{h}_{k,g}$ and $\hat{\textbf{h}}_{k,g}$ from (15) and (6) respectively, we have

[TABLE]

Using the trace lemma and the independence of $\tilde{\textbf{w}}_{k,g}$ and $\textbf{n}^{tr}_{k,g}$ we have,

[TABLE]

Note that $\frac{1}{K^{2}}\text{tr }\textbf{R}_{BS_{g}}\textbf{Q}_{g}\textbf{Z}_{g}\textbf{Z}_{g}^{H}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{C}}^{-1}_{[k,g]}=\kappa_{g}(\textbf{R}_{BS_{g}}\textbf{Q}_{g},\textbf{I}_{N})+o(1)$. Therefore,

[TABLE]

where $\kappa^{o}_{g}(\textbf{A},\textbf{B})$ was derived in the last section and is given by (110). Combining this with Lemma 2 yields the deterministic equivalent of $X_{2}$ as,

[TABLE]

Similar analysis yields the deterministic equivalent of $X_{3}$ as,

[TABLE]

Next note that,

[TABLE]

where $Y_{4}=\frac{1}{K^{2}}\hat{\textbf{h}}_{k,g}^{H}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{h}}_{k,g}$. Using Fubini theorem and trace lemma we have $Y_{4}-\frac{1}{K^{2}}\boldsymbol{\Phi}_{g}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{C}}^{-1}_{[k,g]}\xrightarrow[N\rightarrow\infty]{a.s.}0$. Notice that $\frac{1}{K^{2}}\boldsymbol{\Phi}_{g}\hat{\textbf{C}}^{-1}_{[k,g]}\hat{\textbf{H}}_{[k,g]}^{H}\textbf{P}_{[k,g]}\hat{\textbf{H}}_{[k,g]}\hat{\textbf{C}}^{-1}_{[k,g]}=\chi_{g}+o(1)$, where $\chi_{g}$ is defined in (85). This yields the deterministic equivalent of $X_{4}$ as,

[TABLE]

Combining (128), (134), (135) and (137) yields the deterministic equivalent of $\Upsilon_{k,g}$. The lemmas utilized here not only imply almost sure convergence but also convergence in mean. Therefore,

[TABLE]

F-D Deterministic Equivalent of $\text{var }(\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g})$

Define the quantities $a=\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g}$ and $\bar{a}=\mathbb{E}[{\textbf{h}}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g}]$. By matrix inversion lemma we have $0\leq a,\bar{a}\leq 1$. Thus $\text{var }(\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g})=\mathbb{E}[|a-\bar{a}|^{2}]\leq 2\mathbb{E}[|a-\bar{a}|]$. We have already shown that $a-\frac{h_{g}}{1+m_{g}}\xrightarrow[N\rightarrow\infty]{a.s.}0$. Since $a$ and $\bar{a}$ are bounded, so by dominated convergence theorem we have $\mathbb{E}|a-\bar{a}|\xrightarrow[N\rightarrow\infty]{a.s.}0$. Therefore $\text{var }(\textbf{h}_{k,g}^{H}\hat{\textbf{V}}\hat{\textbf{h}}_{k,g})\xrightarrow[N\rightarrow\infty]{a.s.}0$.

Combining the results of the four subsections completes the proof of Theorem 3.