Interference Modeling for Cellular Networks under Beamforming Transmission

TL;DR

This paper develops and validates analytical interference models for MIMO beamforming cellular networks in different scattering environments, using heavy tail distributions and mixture models to accurately predict interference power.

Contribution

It introduces a three-parameter mixture model combining Inverse Gaussian and Inverse Weibull distributions for interference power, with parameter estimation methods and goodness-of-fit validation.

Findings

The mixture model fits simulated interference data remarkably well.

The model accurately predicts network capacity under beamforming.

Goodness-of-fit is validated using relative entropy metrics.

Abstract

We propose analytical models for the interference power distribution in a cellular system employing MIMO beamforming in rich and limited scattering environments, which capture non line-of-sight signal propagation in the microwave and mmWave bands, respectively. Two candidate models are considered: the Inverse Gaussian and the Inverse Weibull, both are two-parameter heavy tail distributions. We further propose a mixture of these two distributions as a model with three parameters. To estimate the parameters of these distributions, three approaches are used: moment matching, individual distribution maximum likelihood estimation (MLE), and mixture distribution MLE with a designed expectation maximization algorithm. We then introduce simple fitted functions for the mixture model parameters as polynomials of the channel path loss exponent and shadowing variance. To measure the goodness of these models, the information-theoretic metric relative entropy is used to capture the distance from the model distribution to a reference one. The interference models are tested against data obtained by simulating a cellular network based on stochastic geometry. The results show that the three-parameter mixture model offers remarkably good fit to simulated interference power. The mixture model is further used to analyze the capacity of a cellular network employing joint transmit and receive beamforming and confirms a good fit with simulation.

Tables

Table 1. Algorithm 1. MLE Models EM Algorithm

1.	Initialization:	Initialize $w_j^{\left(0\right)}$, ${\lambda }^{\left(0\right)}$, and $c^{\left(0\right)}$, $j\in \{1,2\}$, and compute the initial log-likelihood: $L^{\left(0\right)}$ $=$ $\frac{1}{n}\sum _{i=1}^n\log \left(w_1^{\left(0\right)}{f}_{{\gamma }_{IG}}\left(y_i|{\lambda }^{\left(0\right)}\right)+w_2^{\left(0\right)}{f}_{{\gamma }_{IW}}\left(y_i|c^{\left(0\right)}\right)\right).$ (44)
2.	E-step:	Compute ${\gamma }_{ij}^{\left(m\right)}$ as given in Eq. (39) and $n_j^{\left(m\right)}={\sum }_{i=1}^n{\gamma }_{ij}^{\left(m\right)}$.
3.	M-step:	Compute the new estimate for $w_j^{\left(m+1\right)}$, ${\lambda }^{\left(m+1\right)}$, and $c^{\left(m+1\right)}$, $j\in \{1,2\}$ as given in Eqs. (42), and (43), respectively.
4.	Convergence check:	Compute the new log-likelihood function $L^{\left(m+1\right)}$ $=$ $\frac{1}{n}\sum _{i=1}^n\log \left(w_1^{\left(m+1\right)}{f}_{{\gamma }_{IG}}\left(y_i|{\lambda }^{\left(m+1\right)}\right)+w_2^{\left(m+1\right)}{f}_{{\gamma }_{IW}}\left(y_i|c^{\left(m+1\right)}\right)\right).$ (45)
	Return to 2 if $\left|L^{\left(m+1\right)}-L^{\left(m\right)}\right|\ge \delta$ for a preset threshold $\delta$; Otherwise end the algorithm.

Table 2. (a) IG weight parameter w 1 ( × 10 3 ) w_{1}(\times 10^{3}) .

	$a_3$	$a_2$	$a_1$	$a_0$
2	-111.99	-589.09	-1021.1	-873.55
3	-41.663	-214.46	-443.27	-433.55
3.5	13.4287	-43.972	-230.91	-243.78
4	16.8309	-17.117	-137.71	-149.15
5	5.7399	-1.4256	-40.481	-122.15

Equations

$$\mathbf{H}=\sqrt{l(r)}\tilde{\mathbf{H}}$$

$$l(r)$$

$$\tilde{\mathbf{H}}=\sum _{k_1=1}^K{a}_{k_1}{\mathbf{u}}_{rx}\left({\theta }_{k_1}^{rx}\right){\mathbf{u}}_{tx}^{\ast }\left({\theta }_{k_1}^{tx}\right)$$

$${\mathbf{H}}_0$$

$${\text{y}}_0$$

$$=$$

$${\mathbf{h}}_k$$

$${\text{v}}_0$$

$${\bar{v}}_0={\bar{\text{w}}}_0^{\ast }{\text{v}}_0$$

$$G_k=\{\begin{array}{ll}MM,& \text{mbox}withprobability{p}_1=K(\frac{{\theta }_M}{2\pi }{)}^2\\ mM,& \text{mbox}withprobability{p}_2=2K(\frac{2\pi -K{\theta }_M}{2\pi })(\frac{{\theta }_M}{2\pi })+(K^2-K)(\frac{{\theta }_M}{2\pi }{)}^2\\ mm,& \text{mbox}withprobability{p}_3=(\frac{2\pi -K{\theta }_M}{2\pi }{)}^2\end{array}$$

$${\text{R}}_0={\mathbf{\Sigma }}_0+{\sigma }^2\mathbf{I},$$

$${\tilde{\text{y}}}_0={\text{s}}_d^{\ast }{\text{y}}_0={\text{s}}_d^{\ast }{\text{H}}_0{\text{w}}_0\sqrt{P}{U}_0+{\text{s}}_d^{\ast }{\text{R}}_0{\tilde{\text{z}}}_0,$$

$$C$$

$$C_{IU}={\log }_2\left(1+\frac{|{\lambda }_{max}({\text{H}}_0){|}^2}{{\sigma }^2+{\mathbf{v}}_1^{\ast }{\mathbf{\Sigma }}_0{\mathbf{v}}_1}P\right),$$

$$C_{IA}={\log }_2\left(1+{\Vert {\tilde{\text{s}}}_d^{\star }\Vert }^{2}P\right).$$

$$\displaystyle\!\!\!\!\!\!C={\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\log_{2}}{\left(1+\frac{M^{2}|a_{0}|^{2}\mathit{l}\left({\|\mathbb{p}_{0}\|}_{2}\right)P}{\sigma^{2}+\sigma_{0}^{2}}\right)},$$

$$p_o(R_T)=\mathbb{P}\{C<R_T\},$$

$$D_{KL}(P||Q)={\int }_{-\infty }^{\infty }p(x)\log \frac{p(x)}{q(x)}dx$$

$$f_{{\gamma }_{IG}}(t|\mu ,\lambda )=\sqrt{\frac{\lambda }{2\pi t^3}}\exp \left\{\frac{-\lambda (t-\mu {)}^2}{2{\mu }^{2}t}\right\},t>0,\lambda ,\mu >0.$$

$$\mathbb{E}[{\gamma }_{IG}]=\mu ,\text{var}[{\gamma }_{IG}]=\frac{{\mu }^3}{\lambda }.$$

$$f_{{\gamma }_{IW}}(t|b,c)$$

$$\mathbb{E}[{\gamma }_{IW}]$$

$$q_0$$

$${\sigma }_0^2={\bar{q}}_0=\sum _{k>0:z_k\in {\Phi }_2}{G}_{k}l\left({\Vert z_k\Vert }_2\right)g_k{P}_k={\bar{q}}_0^{(1)}+{\bar{q}}_0^{(2)}+{\bar{q}}_0^{(3)},$$

$$\frac{{\bar{q}}_0^{(i)}}{G^{(i)}}=\sum _{k>0:z_k\in {\mathcal{N}}_i}l\left({\Vert z_k\Vert }_2\right)g_k{P}_k,$$

$$\mathbb{E}\left[q_0\right]$$

$${\zeta }_1$$

$$\nu$$

$$\mathbb{E}\left[{\tilde{q}}_0\right]$$

$$C_{xy}=\frac{{\tilde{q}}_0}{\mathbb{E}\left[q_0\right]},$$

$$\text{var}\left[C_{xy}\right]$$

Full text

Keywords: mmWave cellular; interference model; stochastic geometry; moment matching; maximum likelihood estimation; mixture distribution.

I Introduction

The current scarcity of wireless spectrum coupled with the predicted exponential increase in capacity demand has focused attention on denser network deployment and the underutilized millimeter wave (mmWave) bands (30-300 GHz) [1, 2, 3]. Fifth generation (5G) wireless systems aiming to cater for the capacity demand are expected to provide a minimum of 1 Gb/s data rate anywhere with up to 5 Gb/s for high mobility users and 50 Gb/s data rates for pedestrians [2]. Both dense deployments and mmWave make promising candidates for the 5G systems through the resources reuse over smaller areas and the huge amount of available spectrum [4, 2]. Modeling and characterization of wireless interference under these scenarios is essential for cellular system analysis and design. Due to the different characteristics of signal propagation in conventional microwave and mmWave bands, we use the terminology rich and limited scattering environments to refer to the signal propagation in microwave and mmWave bands.

I-A Background and Related Works

We argue for the important role that interference characterization plays in evaluating and predicting the network performance in linecolor=red,backgroundcolor=red!25,bordercolor=red,], including both microwave and mmWave bands. Traditionally, mmWave bands are considered for backhaul in cellular systems and for high-volume consumer electronics such as personal area and local area networks, but not for cellular access due to concerns about short-range and non-line-of-sight coverage issues [4, 2]. MmWave has, however, recently been shown to be suitable for cellular communications, provided short cell radius of the order of 100-200 meters and sufficient beamforming gain between communicating nodes [2]. Reducing the cell radius leads to dense base station (BS) deployments. Even under beamforming, these high BS and user densities can drive cellular networks to be more interference rather than noise limited. While large adaptive arrays with narrow beams can boost the received signal power and hence reduce the impact of out-of-cell interference [5, 3], this interference remains an important performance-limiting factor in dense mmWave networks [6].

In next generation of wireless networks with a large number of subscribers and dense BSs deployment, interference modeling is an important step towards network realization. Existing stochastic geometry based work either uses the interference Laplace transform to evaluate simple transmission schemes [7, 8], or models the interference using moment matching gamma distribution [9, 10]. Gaussian distribution is another approximation that is considered for the centered and normalized aggregate wireless interference power [11, 12]. The central limit theorem which justifies the Gaussian approximation, however, does not apply when some of the interferers are dominant. Also, the Gaussian distribution does not model the interference very well at low density of interferers or when the exclusion region, the region with no interferers, is relatively small, as the cell sizes shrink. The distribution of the interference power at relatively small exclusion regions has a heavy tail which can not be captured by the Gaussian distribution [12]. Here we propose analytical distribution models characterized by only a few parameters, which can be fitted to simple polynomial functions of channel path loss exponent and shadowing variance, and test their fitness against network simulation based on stochastic geometry. A closed form distribution of the interference helps in designing a cellular system by allowing the analysis of system capacity and comparing different transmission techniques, which may not be directly achievable with a characteristic function of the interference using the Laplace transform. Further, having the interference distribution can significantly speed up system capacity evaluation by generating the interference directly to use in capacity numerical computation, instead of running time-consuming simulations.

To evaluate an interference model, we apply it to a cellular network and analyze the system performance. For cellular network analysis, stochastic geometry is shown to capture the main performance trends with analytical tractability. Stochastic geometry is used to develop a tractable model for downlink heterogeneous cellular networks [8], which is applied to analyze coordinated multipoint beamforming [7]. These stochastic geometry based networks are useful in verifying that the performance of a cellular network matches experimental trends, and can also be used to verify the accuracy of an interference model.

I-B Approaches to Interference Modeling

The underlying question is "What parameterized distribution can best model out-of-cell interference and how can we estimate the parameters of this distribution?". We focus on parameterized distributions with as few parameters as possible to make the model the simplest while having a good fit. We test the goodness of the proposed interference models against simulation of a cellular network based on stochastic geometry. We note that testing against actual measurement data is also feasible and is desirable when such data are available.

To estimate the parameters of the proposed interference models, we use both moment matching and maximum likelihood estimation (MLE) techniques. Specifically, we consider three approaches: analytical moment matching, individual distribution MLE, and mixture distribution MLE. We exploit the iterative expectation maximization (EM) algorithm to estimate the parameters of the mixture MLE model. To evaluate the goodness of each model, we introduce the use of the information-theoretic relative entropy or Kullback-Leibler (KL) divergence as a measure for the relative distance between the modeled distributions and simulated data.

In our interference model development, we consider the interference at a BS in a MIMO uplink scenario, but the results are also applicable to downlink. The developed interference models apply to all transmit beamforming techniques, as long as the beamforming vector is independent of the interference channels. We verify the modeled interference power distribution against network simulation for a wide range of realistic channel propagation parameters.

I-C Main Results and Contributions

We develop out-of-cell interference models that represent the interference power as a random variable with a known, parametrized probability density function. This representation is important for cellular network performance evaluation, prediction and design. We further apply the developed models in analyzing the capacity and outage performance of a MIMO cellular network employing joint transmit and receive beamforming.

The main contributions and novelties of this paper are a new analytical model for interference power distribution and methods for estimating its parameters for both rich scattering and limited scattering environments, which are summarized as:

We propose the use of two distributions to model interference power, each distribution characterized by two parameters: the inverse Gaussian (IG) as a light-to-heavy tailed distribution and the inverse Weibull (IW) as a heavy tailed one. Further, we propose a novel model as a mixture of these two distributions with remarkably good fit to simulation data while having only three parameters. 2. 2)

We apply three approaches to estimate the interference models parameters: moment matching (MM), individual distribution MLE, and mixture distribution MLE. In the MM approach, a simple matching of the first two interference power moments is used. In the other two approaches, a combination of MM and MLE techniques is used in designing an iterative EM algorithm which maximizes a log likelihood function of the interference data. 3. 3)

We propose the use of the relative entropy or Kullback-Leibler distance from information theory to measure the goodness of each model. This metric measures the relative distance between a modeled distribution and reference (simulated) data, which gives a good indication of how far the proposed interference model is from the referenced interference. 4. 4)

We provide simple polynomial functions with fitted coefficients to express the mixture MLE model parameters in terms of channel characteristics, including the path loss exponent and shadowing standard deviation. These polynomials can be used as a simple representation of network interference in complex system-level simulations. 5. 5)

We apply the interference models to evaluate the performance of a cellular network with joint transmit and receive linecolor=red,backgroundcolor=red!25,bordercolor=red,]dominant mode beamforming. User’s outage probabilities using the MM and mixture MLE interference models are compared to stochastic geometry based system simulation. 6. 6)

Our proposed mixture model shows excellent fit in both interference distribution and network performance evaluation for a wide range of channel propagation parameters, including path loss exponent from 2 to 5 and shadowing standard deviation from 0 dB to 9 dB.

This work can be extended to model interference in more complicated propagation environments, including those with a probabilistically parametrized LOS and NLOS channel model, and to more heterogeneous network deployment such as those involving relaying nodes [13, 14, 15].

The importance of our developed mixture model comes from the fact that we can fit the data set of out-of-cell interference into a simple mathematical model, a distribution with known probability density function. This mathematical model can be used to study the performance of a cellular network where model parameters can be specifically tailored for each network setting. These parameters can be derived by applying the MLE method to interference data obtained either via simulation (such as those based on stochastic geometry networks, including point processes and random shapes), or via actual measurement campaigns. The flexibility of adapting the parameters of this mixture distribution makes it a versatile tool in modeling interference.

II Models for network and channel propagation

II-A Network model based on stochastic geometry

We consider a cellular system consisting of multiple cells with an average cell radius $R_{0}$. Each cell has a single BS that is equipped with $N_{\text{BS}}$ antennas and serves multiple user equipments (UEs). Each UE is equipped with $N_{\text{UE}}$ antennas and uses a distinct resource block in each cell, hence there is no intra-cell interference. However, each UE suffers from out-of-cell interference due to frequency reuse in all other cells. In this paper, we denote all transmissions on the same resource block from all cells other than the cell under study as the out-of-cell interference.

Due to the irregular structure of current and future cellular networks, we employ stochastic geometry to describe the network. We consider uplink transmissions in this paper, although results are also applicable to downlink. We model the active UEs in different cells contending for the same resource block and causing interference to each other as being distributed on a two-dimensional plane according to a homogeneous and stationary Poisson point process (PPP) $\Phi_{1}$ with intensity $\tilde{\lambda}_{1}=\eta\lambda_{1}$, where $\eta$ represents the user density factor and can be varied. Furthermore, under the assumption that each BS serves a single mobile in a given resource block, we adopt the model in [16] which places each BS uniformly in the Voronoi cell of its served UE. According to this model, the distance, $r_{c}$, between each BS and its served UE is Rayleigh distributed and has an average value of $R_{c}=1/({2\sqrt{\eta\lambda_{1}}})$. Throughout this paper, we assume that $\lambda_{1}$ is fixed and that $R_{c}$ changes with the parameter $\eta$.

In order to develop a model for the out-of-cell interference power in a stochastic geometry based network, we consider a cell under study with a typical radius $R_{c}$ in a field of active UE interferers with intensity $\tilde{\lambda}_{1}$. As such, our network model consists of a typical cell under study centered at its BS and surrounded by Voronoi cells of other active uplink UEs who interfere with the considered BS at the origin. Our goal is to model this out-of-cell interference to the considered BS receiver, taking into account channel propagation features as discussed next.

II-B Channel propagation model

We consider two channel models for two different environments, the rich and the limited scattering environments, which are typical for signal propagation in the microwave and mmWave bands, respectively.

For the rich scattering environment, we consider a complete channel model with shadowing, path loss and small scale fading. As such, we express a typical MIMO channel with Tx-Rx distance $r$ in the following form:

[TABLE]

where $\mathbf{\tilde{H}}$ is a random matrix that captures the effects of small scale fading and which can be modeled as i.i.d. $\mathcal{CN}({0},1)$; and function $\mathit{l}(r)$ captures the large scale fading which includes both path loss and shadowing. The large scale fading function $\mathit{l}(r)$ is modeled as a log-normal random variable multiplied with the pathloss as a function of the distance $r$ as follows:

[TABLE]

where $\mathit{L}_{s}\sim log(0,\sigma_{SF})$ is a log-normal random variable with standard deviation $\sigma_{SF}$ dB; $\alpha$ is the pathloss exponent; and $\beta$ is the intercept of the pathloss formula. The intercept represents the reference attenuation point that determines the tilt of the path loss model [17].

For the limited scattering environment, as typical in mmWave, the signal propagation characteristics differ linecolor=red,backgroundcolor=red!25,bordercolor=red,]differs considerably. The path loss in mmWave systems is severe with distance-dependent LOS and NLOS propagation, often modeled as a probabilistic function of mixed LOS and NLOS [4, 13]. Further, mmWave signals can be severely vulnerable to shadowing and blockage effects [18]. For limited scattering, however, we only focus on the NLOS interference component in this paper, leaving the composite interference from both NLOS and LOS components as a future work. As such, the large scale fading can be modeled as having a single slope as in (2).

The main feature, due to limited scattering, is that the channel becomes highly directional. Specifically, for the small scale fading, a directional $\mathbf{\tilde{H}}$ can be defined for $K$ scatterers as

[TABLE]

where $K$ denotes the number of channel path clusters; $a_{k_{1}}$ is the complex channel coefficient of the $k_{1}^{th}$ path cluster, assuming (2-D) beamforming; $\theta_{k_{1}}^{rx}$ and $\theta_{k_{1}}^{tx}$ respectively are the horizontal angles of arrival (AoA) and departure (AoD) of the $k_{1}^{th}$ path cluster; and $\mathbf{u}_{rx}(\cdot)\in\mathbb{C}^{N_{\text{BS}}}$ and $\mathbf{u}_{tx}\left(\cdot\right)\in\mathbb{C}^{N_{\text{UE}}}$ are the vector response functions for the BS and UE antenna arrays.

Based on these channel models, we can express the channels of the direct link between the considered active UE and its BS and the interfering links from other active UEs as follows:

[TABLE]

​​​​ where $\mathbb{p}_{0}$ and $\mathbb{z}_{k}$ are vectors representing the 2-$D$ location of the considered active UE and the $k^{th}$ interfering UE in $\Phi_{1}$ with respect to the origin (i.e. the considered BS). Here $\mathbf{H}_{0}$ is the direct channel from the considered active UE and $\mathbf{H}_{k}$ is the channel from the $k^{th}$ interfering UE in $\Phi_{1}$ to the considered BS.

Through a number of recent mmWave channel measurement campaigns, channel parameters for the model in (2) have been identified [19, 18, 20, 21]. For channel propagation at 28 GHz, typical values of the path loss exponent range from 1.68 to 2.55 in LOS and from 2.92 to 5.76 in NLOS environments. Shadowing standard deviation ranges from 0.2 dB to 8.66 dB in LOS and from 8.7 dB to 9.7 dB in NLOS. In this paper, we develop an interference model that works well for these ranges of channel parameters and beyond, as shown later in the numerical section.

III Interference Formulation under MIMO Beamforming

In this section, we present the MIMO beamforming signal model and formulate the out-of-cell interference at the considered BS under transmit beamforming from all active UEs in rich and limited scattering environments. Beamforming has been adopted as an essential technique for mmWave communications to overcome the huge propagation loss at high mmWave carrier frequencies [2, 22]. We then formulate the per-user achievable rate and outage probability under an example of dominant mode joint transmit and receive beamforming in the rich scattering environment, and directional analog beamforming in the limited scattering environment.

III-A MIMO Beamforming Signal Model

We now describe the signal model for single-stream MIMO beamforming, i.e., no spatial multiplexing. We model the received signal $\{\textbf{y}_{0}\in\mathbb{C}^{N_{\text{BS}}\times 1}\}$ at the considered BS as

[TABLE]

where $\textbf{x}_{0}\in\mathbb{C}^{N_{\text{UE}}\times 1}$ is the transmitted signal vector from the considered UE; $\textbf{z}_{0}\in\mathbb{C}^{N_{\text{BS}}\times 1}$ is an i.i.d. linecolor=red,backgroundcolor=red!25,bordercolor=red,]i.i.d noise vector distributed as $\mathcal{CN}(\textbf{0},\sigma^{2}I)$; $\textbf{H}_{0}\in\mathbb{C}^{N_{\text{BS}}\times N_{\text{UE}}}$ is the considered UE-to-BS channel matrix; and $\textbf{v}_{0}\in\mathbb{C}^{N_{\text{BS}}\times 1}$ represents the interference vector received at the considered BS from interfering UEs in all other cells.

The transmit signal vector under beamforming from each user can be generally described as

[TABLE]

where $U$ is a standard Gaussian signal with zero mean and unit variance; $P$ represents the total power allocated to the active UE within a single transmission period; and $\textbf{w}\in\mathbb{C}^{N_{\text{UE}}\times 1}$ is the unit norm beaforming vector at the active UE. We will use this signal model with the appropriate subscript to denote the transmit vector from an active UE, either intended or interfering one.

III-B Interference Signal Formulation

For the purpose of modeling network-wide interference, we emphasize that the distribution of interference is independent of the beamforming scheme employed. The only condition is that the transmit and receive beamforming vectors, which are designed for the intended channel, are independent of the interference channels. This condition is realistic in most practical scenarios.

For the interfering active UE in the $k^{th}$ cell, denote the interference channel as ${\mathbf{{H}}_{k}}$ and the unit-norm interfering beamforming vectors as $\textbf{w}_{k}\in\mathbb{C}^{N_{\text{UE}}\times 1}$, which can represent either a baseband digital precoder or analog RF beamforming or a hybrid of both of these forms.

III-B1 Rich Scattering Environment

Given that $\textbf{w}_{k}$ depends on the direct channel between the $k^{th}$ active UE and its associated BS, it is independent of the interference channel ${\mathbf{{H}}_{k}}$. We can denote the effective interference vector from the $k^{th}$ interfering active UE as

[TABLE]

where ${{\mathbf{\breve{h}}}_{k}}$ is a random vector with i.i.d. elements as $\mathcal{CN}(0,1)$. This follows from the fact that the inner product of a vector of i.i.d. complex Gaussian variables (the column of $\mathbf{\tilde{H}}_{k}$) with an arbitrary unit norm vector is a complex Gaussian variable of the same distribution as each element in the original vector.

The interference vector at the receiving antenna array as presented in (5) can be expressed as

[TABLE]

where the summation is over all interfering users.

III-B2 Limited Scattering Environment

In this environment, we assume directional analog beamforming only and consider the composite effect of both transmit beamforming vector, $\textbf{w}_{k}$, and receive combining vector, $\bar{\textbf{w}}_{0}\in\mathbb{C}^{N_{\text{BS}}\times 1}$, on the received interference power. We further approximate the BSs and UEs antenna array patterns by a sectored antenna model with $M$ and $m$ representing the main and back lobes gains, and $\theta_{M}$ representing the beamwidth of the main lobe. Then, given the independence between $\textbf{w}_{k}$, $\bar{\textbf{w}}_{0}$, and ${\mathbf{{H}}_{k}}$, the received interference signal after both transmit and receive beamforming is

[TABLE]

where $a_{k}$ is the effective channel coefficient and $G_{k}$ is the effective antenna gain from the $k^{th}$ interfering active UE and can be modeled, assuming channel clusters are well separated, as

[TABLE]

where $K$ is the number of clusters (the special case of a single cluster channel model, $K=1$, is discussed in [13, 23]). Note the difference between the formulations in (8) and (9) in that (8) is an interference vector at the receiving antenna array, before receiver processing, whereas the interference in (9) applies after receive beamforming.

III-B3 Modeling the Interference

For each spatial and channel fading realization, based on the large number of interferers, the out-of-cell interference can be modeled as a complex Gaussian random vector for (8) or a complex Gaussian scalar for (9) with zero mean and covariance $\boldsymbol{\Sigma}_{0}$ or variance ${\sigma}_{0}^{2}$. These covariance and variance, however, are random and dependent on user locations and channel fading. The matrix $\boldsymbol{\Sigma}_{0}$ is symmetric with diagonal elements representing the interference power at each receiving antenna element, and the off-diagonal elements representing the correlation among interference signals at different antenna elements. The variance ${\sigma}_{0}^{2}$ on the other hand represents the total interference power after receiver beamforming. We analyze the off-diagonal and diagonal elements of the covariance matrix separately in Sec. V-A. We show later in the numerical analysis section that the correlation among antenna elements is weak and negligible; hence we focus on modeling the interference power elements in $\boldsymbol{\Sigma}_{0}$ and ${\sigma}_{0}^{2}$.

III-C Achievable Rate and Outage Probability Formulation

The main goal of this paper is to model the interference as formulated above. In order to evaluate this model, we apply it in capacity analysis. In this section, we formulate the capacity and show how it is affected by interference.

To establish the per-user capacity in the rich scattering environment, we assume the transmit beamforming vector $\mathbf{w}_{0}$ as the right singular vector corresponding to the dominant mode of $\textbf{H}_{0}$, and consider either an interference-aware (IA) or interference-unaware (IU) receive combining vectors. The noise plus interference in the received signal in (5) can be treated as a Gaussian random vector with a random covariance matrix $\textbf{R}_{0}$, which is dependent on interfering nodes locations and their interference channels. This covariance matrix can be expressed as

[TABLE]

where $\boldsymbol{\Sigma}_{0}$ is the interference covariance matrix discussed in Sec. V-A.

Applying a receive combining vector $\textbf{s}_{d}\in\mathbb{C}^{N_{\text{BS}}\times 1}$ to the received signal ${\textbf{y}}_{0}$ in (5) to get

[TABLE]

where $\tilde{\textbf{z}}_{0}$ is a vector of i.i.d. zero mean and unit variance complex normal entries, we can then write the per-user capacity for the given combining vector $\textbf{s}_{d}$ as

[TABLE]

For interference-unaware (IU) combining, we design the receive combining vector without knowledge of the interference statistics – particularly its covariance. As such, the IU combining vector is chosen as the left singular vector corresponding to the dominant mode of $\textbf{H}_{0}$, $\textbf{s}_{d}=\mathbf{v}_{1}(\textbf{H}_{0})$, resulting in the following per-user capacity:

[TABLE]

where $\lambda_{max}(\textbf{H}_{0})$ is the maximum singular value of $\textbf{H}_{0}$.

As seen in capacity formula (13), with IU combining, the effective post-combining interference is $\sigma_{I}^{2}=\mathbf{v}_{1}^{\ast}\boldsymbol{\Sigma}_{0}\mathbf{v}_{1}$. Since $\boldsymbol{\Sigma}_{0}$ is a random covariance, this post-combining interference is also random. Its distribution, however, can be characterized in the same way as the pre-combining interference power at each antenna element, i.e., the diagonal elements of $\boldsymbol{\Sigma}_{0}$, as stated in Lemma 1. As a result, the IU capacity depends only on the interference power at each antenna element but not on the correlation between interference at different antenna elements.

Lemma 1** (Interference Power Distribution for IU Combining).**

The post-combining interference power $\!\mathbf{v}_{1}^{\ast}\boldsymbol{\Sigma}_{0}\mathbf{v}_{1}$ is a random variable that has the same distribution as a diagonal element of $\boldsymbol{\Sigma}_{0}$.

Proof.

See Appendix AAppendix A: Proof of Lemma 1 for detail. ∎

For interference-aware (IA) combining, the receive combining vector can be designed as a function of the interference statistics. Define $\tilde{\textbf{s}}_{d}=\mathbf{R}_{0}^{1/2}\textbf{s}_{d}$; then replacing $\textbf{s}_{d}=\mathbf{R}_{0}^{-1/2}\tilde{\textbf{s}}_{d}$ into (12), it is straightforward to see that the vector which maximizes (12) is $\tilde{\textbf{s}}_{d}^{\star}={\left(\mathbf{R}_{0}\right)^{-1/2}\textbf{H}_{0}\textbf{w}_{0}}$. The interference-aware per-user capacity can then be obtained as

[TABLE]

This IA capacity depends on the complete interference covariance $\boldsymbol{\Sigma}_{0}$; it is therefore dependent on both the interference power at each antenna element and the correlation between interference at different antenna elements.

In the limited scattering environment, we can express the per-user capacity, assuming perfect beam alignment between the UE and its serving BS [13], as

[TABLE]

where $M$ is the antenna array main lobe gain; $a_{0}$ is the effective complex channel coefficient; $\mathbb{p}_{0}$ is the location of the active UE with respect to the BS; and $\sigma_{0}^{2}$ is the total interference power. Achieving the capacity in (15) does not require knowledge of the interference covariance matrix.

The rates in Eqs. (13) – (15) are random quantities due to channel fading and interference. We investigate the outage probability at a given target rate as a performance metric, defined as

[TABLE]

where $R_{T}$ is the target rate at which the active UE under consideration transmits.

IV Interference Modeling Approaches and Fitness Metric

In wireless cellular networks at system level simulations, the most difficult part is to simulate interference from other cells to the cell under study. This is mainly due to the processing and memory power required. These simulations can be substantially simplified if we can characterize the interference by a simple parametrized distribution that represents a good fit to the actual interference. In this section, we introduce the approaches we use to model the interference power terms on the diagonal of $\boldsymbol{\Sigma}_{0}$ and ${\sigma}_{0}^{2}$. Then, we introduce the information-theoretic based metric used to measure the fitness of each of those models. Finally, we discuss candidate distributions considered for the interference model, the Inverse Gaussian and Inverse Weibull distributions.

IV-A Interference Modeling Approaches

We introduce three different approaches for modeling the distribution of the interference power at the receive antenna. These approaches are Moment Matching (MM), individual distribution MLE (individual MLE), and mixture distribution MLE (mixture MLE).

IV-A1 MM Approach

In this approach, we leverage tools from stochastic geometry to analytically derive the first two moments of each interference power term. Then, we match these two moments with those of the two-parameter candidate distributions to estimate these two parameters and hence specify the distribution. This approach is the simplest and does not require complex computations to estimate the distributions parameters. Thus, it is often used in the literature to produce analytically tractable models [24, 10, 25]. We will show, however, that this approach may not produce a good fit for the model.

IV-A2 Individual MLE Approach

In this approach, to simplify parameter estimation but improve the model fitness, we use a combination of MM and MLE, a powerful estimation technique which is often used in signal detection and estimation [26]. Here we only match the first moment using the analysis results derived based on stochastic geometry, which helps in either determining one of the parameters for the candidate distribution, or getting one parameter in terms of the other. Then, we use MLE to determine the other parameter of the candidate distribution.

IV-A3 Mixture Approach

This approach is the most complex but also most accurate, in which we model the interference probability density function as a mixture of the two candidate distributions. We again use a combination of MM and MLE techniques to simplify parameter estimation. In this approach, we do not determine the optimal values of the distributions parameters in closed form, but design an efficient Expectation Maximization (EM) iterative algorithm to identify these paramters, as discussed in Sec. VI-C.

IV-B KL Divergence Fitness Metric

In order to measure the developed models goodness or fitness, we use a concept from information theory as the relative entropy [27]. In this paper, we measure the relative entropy between a proposed model and the data obtained by simulating the network based on stochastic geometry as discussed in Sec. II-A, but it should be noted that other data such as those obtained from measurement campaigns can also be used to test our models via relative entropy. The relative entropy of a distribution $P$ with respect to another distribution $Q$, also called Kullback-Leibler divergence $D_{KL}(P||Q)$, reflects the difference, or distance, between these two probability distributions. The distribution $P$ usually represents the true distribution of data or observations, while $Q$ is typically used to represent a modeled distribution or an approximation to the true distribution $P$. This relative entropy is defined for continuous probability distributions as

[TABLE]

where $p(x)$ and $q(x)$ denote the probability density functions of $P$ and $Q$, respectively. This definition can also be described as the expectation of the logarithmic difference between the probabilities $P$ and $Q$. Kullback-Leibler divergence metric has been used to measure the accuracy of modeling composite fading and shadowing wireless channels such as the Rayleigh-Inverse Gaussian and mixture Gamma distributions [28, 29].

IV-C Candidate Distributions for the Interference Model

We consider two candidate distributions for interference modeling: the inverse Gaussian (IG) and the inverse Weibull (IW), both are characterized by two parameters. The IG distribution is also known as the Wald distribution and is usually used to model nonnegative positively skewed data. For example, a composite Rayleigh-IG distribution is used to approximate the composite Rayleigh-Lognormal distribution in [28]. The name "Inverse Gaussian" comes from the inverse relationship between the cumulant generating functions of these distributions and those of Gaussian distributions. It is a two-parameter family of continuous probability distributions which is specified by a shape parameter ${{\lambda}}$ and a scale parameter ${{\mu}}$. Given an IG random variable $\gamma_{IG}({{\mu}},{{\lambda}})$, its probability density function is defined as [30, 31, 32]

[TABLE]

The mean and variance of $\gamma_{IG}({{\mu}},{{\lambda}})$ can be written in the following form [32]

[TABLE]

The IW distribution is also known as the complementary Weibull distribution and is usually used to model nonnegative positively skewed data that exhibits a long right tail. It is also specified by two parameters, a shape parameter ${{c}}$ and a scale parameter ${{b}}$ (or $\lambda_{IW}$ in some text). Given an IW random variable $\gamma_{IW}({{b}},{{c}})$, its probability density function is defined as [33, 34]

[TABLE]

linecolor=red,backgroundcolor=red!25,bordercolor=red,]changed comma to full stop after Eq. (20)The mean and variance of $\gamma_{IW}({{b}},{{c}})$ can be written in the following form [33]

[TABLE]

where the Gamma function $\Gamma(t)$ is defined as $\Gamma(t)=\int_{0}^{\infty}x^{t-1}e^{-x}dx$.

In a recent and independent work, the Gaussian, Gamma, Inverse Gamma, and Inverse Gaussian distributions are considered as models for the interference distribution in a stochastic geometry based spatial network [35]. Different point processes are considered in this work, including PPP, Strauss process, and Poisson cluster process. The inverse Gaussian distribution is shown to be a reasonable model for the interference created by BSs distributed according to the Strauss point and Poisson cluster processes. Further, both the Gamma and IG were shown to be suitable for the PPP network geometry. This prior work, however, considered a bounded path loss model $l(r)=(1+r^{\alpha})^{-1}$ and Rayleigh fading with no shadowing.

In our work, we find through simulations that the Gaussian, Gamma and Inverse Gamma distributions are all poor models for the interference power when shadowing is present. Shadowing introduces large variation in the interference power and causes a medium to heavy tail, which none of these distributions are capable of capturing. Thus we need to find other distributions that can more closely model heavy tail interference. As such, we examine the Inverse Gaussian and also introduce the Inverse Weibull as a good candidate for parameterized heavy tail distributions that can capture shadowing with large standard deviation.

V Analytical Moments of Interference and Moment Matching Models

In this section, we analytically derive the first two moments of the diagonal and off-diagonal elements of the interference covariance matrix. We then use the derived moments of the interference power at each antenna element in rich scattering, and the post receive beamforming interference power $\sigma_{0}^{2}$ in the limited scattering environment, to develop the MM interference models.

V-A Analytical Moments of Interference

V-A1 Interference Power

In the rich scattering environment, the diagonal elements of the interference covariance matrix $\boldsymbol{\Sigma}_{0}$ represent the interference power at each antenna element and they are uncorrelated with the off-diagonal elements. A correlation between any two of these diagonal elements exists and can be determined similarly to the second moment of the off-diagonal elements (Lemma 3). However, since this correlation is not strong as shown in the numerical analysis in Sec. VII-A, we model the diagonal elements of $\boldsymbol{\Sigma}_{0}$ as a random vector of independent elements. Each of these diagonal elements can be generally expressed as

[TABLE]

where ${g}_{k}$ are all $\text{i.i.d. }\mathcal{\exp}(1)$ and represent the power gain of the small-scale Rayleigh channel fading from the $k^{th}$ interferer.

Similarly, in the limited scattering environment, the interference power can be expressed as

[TABLE]

where we assume that the interferers are distributed according to PPP $\Phi_{2}$ with intensity $\lambda_{2}$. Using the thinning property of random point processes [23], from $\Phi_{2}$, we define three independent PPPs $\mathcal{N}_{i},i\in\{1,2,3\}$ with intensities $\tau_{i}=p_{i}\lambda_{2}$. Then, we can define the normalized interference power terms as

[TABLE]

where $G^{(i)}\in\{MM,Mm,mm\}$. Comparing Eqs. (22) and (24), we note that the distribution of the normalized interference power for each thinned process in (24) can be characterized as that of the interference power in (22) by matching the intensity $\tau_{i}$ to the intensity of the interferers in the rich scattering environment, $\tilde{\lambda}_{1}$. Hence, we focus on the characterization of the rich scattering environment interference power in (22), based on which the three separate interference power terms for the limited scattering environment in (23) can be derived.

Lemma 2 characterizes the first two moments of the out-of-cell interference power $\mathit{q}_{0}$ at each receiving antenna.

Lemma 2** (Interference Power Moments).**

*For network-wide deployment of MIMO transmit beamforming, the out-of-cell interference power at each antenna of the interfered receiver has the following first and second moments: *

[TABLE]

[TABLE]

Proof.

See Appendix BAppendix B: Proof of Lemmas 2 and 3 for details. ∎

V-A2 Interference Correlation among Different Antennas

The first moment of the off-diagonal elements of the covariance matrix $\boldsymbol{\Sigma}_{0}$ represents the direct correlation among interference signals at different receiving antenna elements. The second moment, however, represents linecolor=red,backgroundcolor=red!25,bordercolor=red,]represent the correlation between received interference powers. These off-diagonal elements are expressed as in (47) for $i\neq j$ in Appendix A. Assuming antenna elements are spaced sufficiently far apart such that they experience independent fading, Lemma 3 characterizes their first and second moments.

Lemma 3** (Interference Correlation Moments).**

*For network-wide deployment of MIMO transmit beamforming, the correlation among interference at different receiving antenna elements of the BS has the following first and second moments: *

[TABLE]

*where $\text{var}\left[\mathit{q}_{0}\right]$ is defined in Lemma 2. *

Proof.

See Appendix BAppendix B: Proof of Lemmas 2 and 3 for details. ∎

The interference correlation coefficient $C_{xy}$ can then be written as

[TABLE]

Based on Lemmas 2 and 3, we can derive the variance of this correlation coefficient as

[TABLE]

This variance increases with increasing path loss exponent, $\alpha$, and shadow fading standard deviation, $\sigma_{SF}$. We perform numerical analysis of this correlation coefficient in Sec. VII-A and show that the correlation among interference at different antenna elements is in general weak.

V-B Moment Matching Interference Models

Given the first two moments of the interference power in Lemma 2, we can use moment matching as the simplest and most straightforward method to estimate the parameters of the candidate distributions. This method, however, may not produce the best fit. Specifically, we use the two analytically derived moments to estimate the $({{\mu}},{{\lambda}})$ and $({{b}},{{c}})$ parameters of the IG and IW distributions. The fitted distribution can then be used to represent the marginal distribution of the interference power at each receiving antenna element, as stated in Lemmas 4 and 5 next.

Lemma 4** (Estimation of IG Parameters using MM).**

Given the derived moments in Lemma 2, the shape and scale parameters ${{\lambda}}$ and ${{\mu}}$ of an IG model for the interference power $\gamma_{IG}[{{\mu}},{{\lambda}}]$ can be estimated as

[TABLE]

Proof.

Follows directly from the definition of $\lambda$ and $\mu$ in (19). ∎

Lemma 5** (Estimation of IW Parameters using MM).**

Given the derived moments in Lemma 2, the shape and scale parameters ${{c}}$ and $b$ of an IW model for the interference power $\gamma_{IW}[{{b}},{{c}}]$ can be estimated, assuming ${{c}}>2$, by solving the following equation:

[TABLE]

If the above equation gives a value of $c<2$, then set $c=2.01$ and compute $b$ accordingly.

Proof.

Follows directly from the equations for $b$ and $c$ in (21). ∎

The equation in (31) can be solved numerically to obtain ${{c}}$, which is then used to obtain ${{b}}$.

Remark 1**.**

For the expressions in (21), we note that the mean is infinite for ${{c}}<1$ and that the variance is infinite for ${{c}}<2$. This is a consequence of the heavy right tail of the IW distribution. Hence, for the IW moment matching in Lemma 5, we assume that ${{c}}>2$.

Remark 2**.**

For the path loss exponent value $\alpha=2$, we use numerically evaluated mean and variance to perform moment matching since the expressions in Lemma 2 can not be evaluated analytically at this value of $\alpha$.

VI Individual and Mixture MLE Interference Models

In this section, we apply maximum likelihood estimation to compute the model parameters. We discuss both the individual MLE and mixture interference models and develop an EM algorithm to estimate the model parameters. Maximum likelihood is one of the dominant methods of estimation in statistics. A typical MLE problem is to estimate the parameter set $\theta$ of a candidate distribution $p(y|\theta)$ from a given data set $y$ and can be expressed as an optimization as follows:

[TABLE]

where $\log p(y|\theta)$ is the log likelihood of the observed data $y$, and $\Theta$ is a closed convex set.

VI-A Individual MLE Interference Models

In individual MLE interference model, we use only a single distribution from the candidates to model the interference power at each antenna. We use the maximum likelihood optimization problem in (32) to estimate the parameters of each distribution. However, to simplify the estimation, we use moment matching first. In moment matching first, we match the mean, $\mu_{Y}$, of the observed data set to that of each of the candidate distributions. This matching results in the scale parameters ${{\mu}}$ and ${{b}}$ for the IG and IW distributions as

[TABLE]

Next, we replace these expressions of the scale parameters into the log-likelihood function in (32) for each candidate distribution. In the case of IG distribution, this substitution results in a single-parameter log-likelihood function as follows

[TABLE]

Hence, the MLE problem in (32) becomes very simple and has a closed form solution for ${{\lambda}}$. Similarly, in the case of IW, the substitution results in a single-parameter log-likelihood function

[TABLE]

The derivation and final expressions of the optimal parameters, ${{\lambda}}$ and ${{c}}$, that optimize the fit of each candidate distribution are presented in Sec. VI-C.

VI-B Mixture MLE Interference Model

In our mixture model, the probability density function (pdf) of the interference model is now a mixture of both IG and IW distributions as defined in the next theorem.

Theorem 1**.**

Given a data set $Y$, the pdf of the mixture interference model can be written as

[TABLE]

where $\{w_{1},w_{2}:\sum_{i=1}^{2}w_{i}=1\}$ are the weight parameters for mixing the IG and IW distributions and $f_{\gamma_{IG}}(y|{{\lambda}})$ and $f_{\gamma_{IW}}(y|{{c}})$ are as given in Eqs. (34) and (35) with the scale parameters obtained through moment matching first as in (33). This mixture model is a 3-parameter distribution with the parameter set as $\theta=\{w_{1},\lambda,c\}$.

Given the pdf of the mixture model in Theorem 1, we note that the individual MLE models are special cases of this mixture model. By setting $w_{1}=1$ , we have the IG individual MLE model and by setting $w_{2}=1$, we have the IW individual model. Subsequently we focus on the mixture model parameters estimation, which can be applied to the individual models as well.

VI-C MLE Models Parameters Estimation

We have three parameters to be optimized in the mixture model which is significantly more than the single parameter in the individual MLE model. This mixture leads to a more complex MLE optimization problem and the EM algorithm [36] becomes an appealing approach. Here, we identify the complete data $X$ as the observed data $Y$ plus some hidden data $Z$, i.e. $X=(Y,Z)$, where $Y$ is the observed set of points that we model by a weighted IG and IW distribution and $Z\in\{1,2\}$ is a discrete random variable representing the assignment of each data point to the two candidate distributions. Then, we can define the pdf of the complete data $X\in\mathbb{R}^{+}$ as

[TABLE]

where $\phi_{1}(y|\theta_{1})=f_{\gamma_{IG}}(y|{{\lambda}})$ and $\phi_{2}(y|\theta_{2})=f_{\gamma_{IW}}(y|{{c}})$.

Next, we present the algorithm for estimating the parameters of the MLE interference models. We focus on the representative case of the mixture MLE model, which includes the individual MLE models as special cases. Given $n$ observations, in order to develop the EM algorithm to estimate the parameters of our mixture model, we use Proposition 1 below [37] which can be easily verified using the independence assumption, the fact that the observed data $Y=f(X)$ is a deterministic function of the complete data $X$ for some function $f$, and Bayes’ rule.

Proposition 1**.**

[37]* Let the complete data $X$ consist of n i.i.d. samples: $X_{1},X_{2},...,X_{n}$, which satisfies $f(X|\theta)=\prod_{i=1}^{n}f(X_{i}|\theta)$ for all $\theta\in\Theta$, and let $y_{i}=f(x_{i})$, $i=1,2,...,n$, then*

[TABLE]

We also define a responsibility function $\gamma_{ij}^{(m)}$, which represents our guess at the $m^{th}$ iteration of the probability that the $i^{th}$ sample belongs to the $j^{th}$ distribution component, as

[TABLE]

where $\theta_{1}$ and $\theta_{2}$ represent the shape parameters $\lambda$ and $c$, respectively.

An EM algorithm includes two steps: the E-step to calculate the conditional expectation of the log likelihood of the complete data, and the M-step to maximize this conditional expectation function. Lemma 6 provides the final expression for the EM algorithm $Q$-function in the E-step.

Lemma 6**.**

*Let the complete data $X$ consist of n i.i.d. samples: $X_{1},X_{2},...,X_{n}$, which satisfies $f(X|\theta)=\prod_{i=1}^{n}f(X_{i}|\theta)$ for all $\theta\in\Theta$, and let $y_{i}=T(x_{i})$, $i=1,2,...,n$, then we have *

[TABLE]

where $n_{j}^{(m)}=\sum_{i=1}^{n}\gamma_{ij}^{(m)},j\in\{1,2\}$; and ${{{\theta}}}=\{w_{1},\lambda,c\}$.

Proof.

See Appendix C for details. ∎

In the M-step, to update our estimate of the parameter ${{{\theta}}}$, we solve the following optimization problem for the optimal ${{{{\theta}}}}^{*}$:

[TABLE]

The optimal ${{{{\theta}}}}^{*}$ is then found as in Theorem 2.

Theorem 2**.**

The optimal new estimate of ${{{\theta}}}^{(m+1)}=\left\{w_{1}^{(m+1)},{{{\lambda}}}^{(m+1)},c^{(m+1)}\right\}$ that maximizes the Q-function in (40) at the $m^{th}$ iteration are determined as follows:

[TABLE]

and ${{c}}^{(m+1)}$ is obtained by solving the following equation numerically:

[TABLE]

Proof.

See Appendix D for details. ∎

Algorithm 1 summarizes the EM algorithm for the mixture distribution interference model.

In the case of IG individual MLE model, we can set $w_{1}=1$ and the resulting $Q$-function in (52) represents the log-likelihood function of the IG distribution. Hence, problem (41) results in the MLE of a single parameter, ${{\lambda}}$. The optimal value of ${{\lambda}}$ can be obtained from (42) by setting $n_{1}^{(m)}=n$ and $\gamma_{i1}^{(m)}=1$. Similarly, the optimal value of ${{c}}$ for the IW individual MLE model can be obtained from the numerical solution of (43) after setting $w_{2}=1,n_{2}^{(m)}=n$ and $\gamma_{i2}^{(m)}=1$.

VI-D Functional Fitting of MLE Model Parameters to Channel Propagation Characteristics

By performing the MLE parameters estimation, we observe that these parameters can be fitted directly to rather simple functions of the channel features, including the shadowing variance and path loss exponent. These functions present excellent fit with simulated data for the range of mmWave propagation parameters as measured in recent campaigns [19, 18, 20, 21].

Corollary 1 presents the polynomial function form that we use to fit the MLE parameters to the channel shadowing variance. Note the same functional form applies to all three parameters of the mixture MLE model, each parameter with a different set of fitting coefficients.

Corollary 1**.**

For a given path loss exponent, the parameters of the mixture MLE algorithm can be estimated directly as a function of the shadowing standard deviation $\sigma_{SF}$ as

[TABLE]

where $\theta$ represents any of the parameters $\{w_{1},{{c}},{{\lambda}}\}$. The function coefficients for each parameter are obtained via numerical fitting based on the EM estimated parameter values. Parameters ${{b}}$ and ${{\mu}}$ can then be obtained analytically based on Eq. (33) and Lemma 2.

In Appendix E, we provide representative values of the functional coefficients fitted against simulation data based on stochastic geometry.

Similarly, the MLE parameters can be fitted to functions of the channel path loss exponent using a simple least square optimization procedure. Corollary 2 presents the polynomial function forms that we use to fit the MLE parameters to the path loss exponent.

Corollary 2**.**

At a fixed value of the shadowing variance, the IG weight $w_{1}$ and the IW shape parameter $c$ of the mixture MLE model can be estimated directly as a $4^{th}$ order polynomial function of the path loss exponent $\alpha$. The logarithm of the IG shape parameter $\log_{10}(\lambda)$ can be estimated as a polynomial function of $\alpha$ of the $3^{rd}$ degree. The coefficients for each function can be obtained via least square fitting based on the EM estimated parameter values.

Using these functional fits, we can bypass the MLE iterative algorithm in Algorithm 1 and directly use the table of fitted function coefficients. These functional fits therefore provide a simple and powerful way of modeling the interference and its direct dependence on the channel propagation features.

Note that the coefficients examples in Appendix E are verified against simulated data using stochastic geometry based network settings and need further verification against actual measurement data, for which the coefficient values may change slightly. We expect, however, that the identified function forms in (46) remain a good fit. To fit the function coefficients, apply the MLE algorithm to the measured data to identify the interference model parameters. Then fit these model parameters to the function forms in (46) to identify the coefficient values. Once the fit is verified, the function with the identified coefficients can be used directly to model the interference in a new environment with given shadowing variance and path loss exponent, without the need to measure new data or perform the MLE algorithm again.

VII Numerical Performance Results

In this section, we use simulation data obtained from the stochastic geometry based network setting in Sec. II-A to numerically verify the validity of the analytical results and proposed interference models. We first verify the validity of our moment matching interference models. We then discuss the fit of the MLE individual and mixture interference models. Last, we apply these interference models to evaluate the outage performance of a cellular system employing MIMO joint transmit-receive dominant mode beamforming.

Since interference in the limited scattering environment can be derived as the sum of interference components from independent, thinned processes, each of which exhibits similar characteristics as the interference in a rich scattering environment but with a different interferer density as in [23], we focus on the rich scattering environment in our simulation. For the simulations setting, we consider a typical cell of radius $R_{c}=150$ m, which is related to the active UEs density $\lambda_{1}$ as $\eta\lambda_{1}=0.25R_{c}^{-2}$, where $\eta=1$ in the typical case. Later we vary $\eta$ to examine the effect of user density on performance.

VII-A Interference Correlation among Different Antenna Elements

This part aims to examine the correlation between interference at different antenna elements to see if this correlation is negligible. We use a kernel distribution to fit the actual interference correlation. A kernel distribution is a nonparametric representation of the probability density function of a random variable. It can be used when a parametric distribution cannot properly describe the simulation data, or when assumptions about the distribution of the data are to be avoided, as is the case for testing the fit of our proposed parameterized interference distributions.

We examine the interference correlation via the correlation coefficient defined in Sec. V-A2. In Fig. 1, we plot the probability that the correlation coefficient is smaller than 0.3, and show that this probability is significant for a wide range of the shadowing standard deviation and the path loss exponent. This figure demonstrates that the interference correlation among different antenna elements is in general weak. As shown later in the system performance evaluation (Sec. VII-E), the impact of interference correlation on transmission rate differs depending on the receive combining vector, as discussed in Sec. III-C. For IA combining, interference correlation can have an impact on capacity; whereas for IU combining, the interference correlation has no impact and can be ignored. In this paper we do not model the interference correlation, hence IU combining is more appropriate.

VII-B Moment Matching Interference Models Evaluation

We now use the relative entropy metric described in Sec. IV-A to measure the goodness of the MM interference models in Section V-B. In Fig. 2(a), we plot the relative entropy versus shadowing standard deviation $\sigma_{SF}$ at a sample path loss exponent value of $\alpha=3.5$. The result shows that IG approximates the simulated interference distribution well at low values of shadowing standard deviation $\sigma_{SF}$ and diverges at higher values.

The IW distribution, on the other hand, does not show a good fit at the simulated path loss value. Our extensive simulation and testing suggests that IW can be a good approximation of the simulated data distribution only for low values of $\alpha$ (lower than shown in Fig. 2). The moment matching IW distribution diverges significantly as the path loss exponent increases. This divergence is caused by the fact that the optimal shape parameter ${{c}}$, as we will see later in this section, obtained at high values of $\alpha$ is close to or less than $2$. At this range of ${{c}}$ values, the variance of the IW distribution is infinite and the moment matching model is not applicable.

In Fig. 2(b), we show a sample interference power distribution at one antenna element for $\alpha=3.5$ and $\sigma_{SF}=4$ dB to visually confirm the match, or rather the mismatch, of the moment-matched IW interference distribution. The moment-matched IG distribution also starts to diverge from the simulated data at this value of $\sigma_{SF}$, as indicated in Fig. 2(a).

VII-C Individual and Mixture MLE Interference Models Evaluation

Next, we evaluate the fit of the individual and mixture MLE interference models. Fig. 3 shows the relative entropy measures for these 3 models against simulated data as well as representative sample distributions. In Fig. 3(a), we see that our proposed mixture model outperforms both the individual IG and IW models. This result aligns with our expectation since the pdf of the mixture model is a combination of both the individual IG and IW distributions. Hence, performance of the mixture model should be at least as good as that of the individual models. The fact that the mixture model performs better than either individual one is because the relative entropy measure is a non-linear function of the distribution densities.

Comparing between Figs. 2 and 3 also reveals that all MLE models perform much better than the MM models: at the worst $\sigma_{SF}=9$ dB, the IW distribution under individual MLE has a maximum relative entropy of about merely $0.23$ compared to a maximum of $1.8$ under MM. Further, we see that the proposed mixture model follows the simulated data distribution almost exactly even at large shadowing: the maximum relative entropy is barely $0.03$.

Consistently in all our extensive simulations and testing, the mixture MLE interference model provides the best fit. This model requires an iterative algorithm to fit the parameters initially, before fitting the functional coefficients as in Corollary 1. The MM models are attractive from the simplicity point of view since the model parameters can be determined analytically without optimization. However, our simulation results suggest that the IG MM model should only be used for the range of shadowing standard deviation $\sigma_{SF}<3$ dB, at which the relative entropy of the IG MM model is approximately within $1\%$ of the mixture model.

In Fig. 5, we plot the number of iterations required by Algorithm 1 to obtain the optimal values of the mixture MLE model parameters versus $\sigma_{SF}$ at two samples of the path loss exponent $\alpha=\{2.5,4\}$, where the preset threshold value is chosen to be $\delta=10^{-6}$, the model weights are initialized as $w_{1}^{(0)}=w_{2}^{(0)}=0$, and the model shape parameters $\lambda^{(0)}$ and $c^{(0)}$ are initialized to the optimal values of the corresponding individual MLE models. The number of iterations required by the EM algorithm to obtain the optimal mixture MLE model parameters is in general less than 200. In Fig. 5, we study the sensitivity of Algorithm 1 to the initial values of the parameters. We consider the sample data obtained at path loss exponent $\alpha=2$, shadowing value $\sigma_{SF}=4$ dB, and density factor $\eta=1$. In this figure, we plot the IG weight parameter $w_{1}^{(m)}$ at each iteration $(m)$ for different initial values $w_{1}^{(0)}$, which shows that the weight parameter converges to the same value after at most 170 iterations. We also note from simulations that the likelihood function is concave in both of the shape parameters $c$ and $\lambda$ at a fixed value of the weight parameter, therefore it has a unique maximum. Hence, we can conclude based on numerical evidence that Algorithm 1 converges to the same optimal parameter values irrespective of their initializations.

VII-D Mixture Model Parameters and Function Coefficients Fitting

In this section, we use the EM algorithm discussed in Sec. VI-C. along with Theorem 2 to obtain the mixture model parameters. We then fit the estimated parameters $w_{1},{{c}},{{\lambda}},$ into a function in terms of $\sigma_{SF}$ as in Corollary 1. After that, Eq. (33) can be used to determine the parameter ${{b}}$. We can also use (25) to determine the parameter ${{\mu}}$ directly, but instead we provide a function fitted to the simulated data for this parameter, since (25) is invalid at $\alpha=2$. The coefficients of the fitted functions are provided for representative values in Appendix E.

In Fig. 7, we plot the mixture model weight factor obtained from both simulation (cross markers) and functional fit (lines). The estimated parameters and fitted functions show an excellent fit to the simulated data. The results confirm that larger shadowing results in an increase in the IW component and a decrease in the IG component in our mixture interference model. The increase in shadowing standard deviation $\sigma_{SF}$ implies a higher probability of high interference power, leading to a heavier interference distribution tail. As the path loss exponent $\alpha$ increases, however, the IW distribution starts to diverge from the simulated data and hence reduces in weight, thus the IG weight $w_{1}$ increases as shown. Note also that as $\alpha$ increases, shadowing has less effect on the mixture model weights, especially at low active UEs density. As such, the weight distribution at high path loss is less sensitive to shadowing.

In Fig. 7, we plot the IW and IG shape parameters $\{{c},\lambda\}$ versus $\sigma_{SF}$ at $\alpha=3$. These curves match the simulated data perfectly, verifying the functional fits in Corollary 1. At higher path loss exponents, however, we noted some irregularities in the values of simulated ${{c}}$ specifically at $\alpha\geq 4$ and $\sigma_{SF}\leq 4$ (not shown in the figure). This irregularity is caused by the fact that the IG distribution has a much higher weight than IW in this range, hence the IW distribution can assume a looser range of values for its parameters without affecting the goodness of the fit. The functional fit for $c$ in (46) thus can still be applied in this high path loss range.

These plots demonstrate that the fitted functions for parameters of the mixture MLE model are accurate for the practical ranges of channel path loss and shadowing. These functions therefore provide an extremely simple yet accurate way to estimate interference distribution parameters.

VII-E Rate CDF and Outage Performance

We now apply the proposed interference models to evaluate the cellular system performance. In particular, we evaluate the transmission rate cumulative distribution function (CDF), which is related to the outage probability, of a cell edge user at a distance of $145$ m from the base station with $\eta=1$ (note $R_{c}=150$ m and $\eta\lambda_{1}=0.25R_{c}^{-2}$). In this simulation scenario, we do not consider uplink scheduling; instead, we consider a random user with a random channel to the considered BS. We further normalize the simulation results by the channel bandwidth and show rates in bps/Hz. The uplink maximum transmission power is $30$ dBm. We use the dominant mode transmit beamforming discussed in Sec. III-C and assume a noise variance of $\sigma^{2}=-124$ dBm/Hz and a path loss intercept $\beta=-72.3$ dB.

In Fig. 8, we plot the interference-aware capacity in (14) and compare the user transmission rate CDFs based on IG MM and mixture interference models to the simulated data. Both models show good accuracy for the interference power as in the left plot at a low $\sigma_{SF}$ value. At larger shadowing, however, as in the right plot with $\sigma_{SF}=9$ dB, only the mixture model remains accurate while the IG MM model diverges significantly. This result also suggests that interference has a strong impact on user transmission rate CDF even at the somewhat high noise variance of $\sigma^{2}=-124$ dBm/Hz. We also note that the slight mismatch between the mixture model and simulation is due to the interference correlation, which we did not model but does have an impact on capacity in an interference-aware combining scheme.

To further examine the accuracy of the mixture interference power model, we consider the interference-unaware combining scheme with capacity in (13). IU capacity does not depend on interference correlation as a result of Lemma 1. Results in Fig. 9 shows that the capacity using interference power model matches with simulation almost perfectly in all considered antenna configurations, even at a high shadowing variance of $\sigma_{SF}=9$ dB. This result confirms the validity of our interference power distribution in capacity evaluation.

Our interference distribution models have been validated against simulation based on stochastic geometry. The ideal next step would be to evaluate the models against actual measurement data, but at the point of writing this paper, such data are not available. As stochastic geometry has been shown to reasonably capture cellular performance compared to actual data [8], we expect that the mixture model with adjusted functional coefficients will also present a good fit to actual measurement data. The mixture MLE model can provide a valuable tool for cellular system performance evaluation and prediction.

VIII Conclusion

We introduced new interference models based on moment matching and maximum likelihood estimation techniques for both the rich scattering and limited scattering NLOS propagations. A novel model as a mixture between the Inverse Gaussian and Inverse Weibull presents remarkably good fit with simulation, capturing the heavy-tail characteristics of interference especially in high shadowing environments. We designed an expectation maximization algorithm to estimate the parameters of this mixture model. We also fitted the mixture model parameters to simple polynomial functions of the channel propagation features, which provide a simple way to model the interference without any optimization. The fitted mixture model can then be integrated into a more complex system level simulation to evaluate and predict cellular performance, or used to aid system design. The next step would be to test the mixture model with appropriately adjusted parameters and functional fitting coefficients against data from measurement campaigns to evaluate the model in an actual system setting.

Appendix A: Proof of Lemma 1

Given the interference vector in (8), the element in the $i^{th}$ row and $j^{th}$ column of the interference covariance matrix $\boldsymbol{\Sigma}_{0}$ can be expressed as

[TABLE]

where $\breve{h}_{(k,i)}$, $\breve{h}_{(k,j)}$ are i.i.d. $\mathcal{CN}(0,1)$ representing the Rayleigh small-scale fading channels from the $k^{th}$ UE to the $i^{th}$ and $j^{th}$ antenna element of the considered BS, respectively. The post-combining interference power term $\mathbf{v}_{1}^{\ast}\boldsymbol{\Sigma}_{0}\mathbf{v}_{1}$ can then be written as

[TABLE]

Since the direct channel small-scaling fading component is assumed to be independent of interference channels small-scale fading, $\mathbf{v}_{1}$ is independent of all $\mathbf{\breve{h}}_{k}$. Thus the distribution of $\sum_{i}[\mathbf{v}_{1}]_{i}\breve{h}_{(k,i)}$ is $\mathcal{CN}(0,1)$ as $\mathbf{v}_{1}$ is unit-norm, which renders the distribution of $\mathbf{v}_{1}^{\ast}\boldsymbol{\Sigma}_{0}\mathbf{v}_{1}$ in (48) the same as that of the diagonal elements of $\boldsymbol{\Sigma}_{0}$ in (47).

Appendix B: Proof of Lemmas 2 and 3

We first derive the $Laplace$ transform of the interference power term, $\mathit{{q}}_{0}$, which is then used to characterize its moments. We develop the $Laplace$ transform of the interference power at the destination as follows:

[TABLE]

where the last equality follows from the $Laplace$ functional expression for PPP using polar coordinates and assuming the field of interferers outside a cell of fixed radius $R_{c}$; $\mathcal{L}_{T}(s)$ is the Laplace transform of the composite shadowing and small scale fading channel, where $T\buildrel d\over{=}\mathit{L}_{s}G$ with $G\sim\exp(1)$ as an exponential random variable, and $\mathcal{L}_{\mathcal{J}_{0}}\left(s,{\|\mathbb{z}_{k}\|}_{2}\right)$ is expressed as

[TABLE]

Then, we obtain the results in Lemma 2 based on the interference power $Laplace$ transform and evaluating the following formulas:

[TABLE]

The proof to Lemma 3 follows the same procedure with the only difference in that we replace $\mathit{{q}}_{0}$ in Eqs. (Appendix B: Proof of Lemmas 2 and 3) and (51) with $\mathit{\tilde{q}}_{0}$, and the Laplace transform, $\mathcal{L}_{T}(s)$, in eq. (50) with the Laplace transform $\mathcal{L}_{T^{*}}(s)$ of the composite random variable $T^{*}\buildrel d\over{=}\mathit{L}_{s}\tilde{H}_{1}\tilde{H}_{2}$ where $\tilde{H}_{1},\tilde{H}_{2}\sim\mathcal{CN}(0,1)$ are independent and identically distributed complex normal random variables.

Appendix C: Proof of Lemma 6

We start by writing the expression for the EM algorithm Q-function using Proposition 1 as

[TABLE]

Using the probability density functions of IG and IW distributions in Eqs. (18) and (20) along with moment matching first, we can write $\log\phi_{1}(y_{i}|\lambda)$ and $\log\phi_{2}(y_{i}|c)$ as in Eqs. (34) and (35), respectively. Then, substituting into Eq. (52), replacing $\sum_{i=1}^{n}\gamma_{ij}^{(m)}$ by $n_{j}^{(m)}$, and dropping constant terms without affecting the EM algorithm, we get the final expression in Eq. (40).

Appendix D: Proof of Theorem 2

We first solve for $w_{j}$ using the method of Lagrange multipliers as follows.

[TABLE]

Substituting into the constraint $\sum_{j=1}^{2}w_{j}=1$, we get the weight $w_{j}$ update as

[TABLE]

Then, we can find the update of $\lambda_{IW}$ by simply taking the partial derivative of the $Q$-function w.r.t. $\lambda_{IW}$ as

[TABLE]

Solving for ${{\lambda}}^{(m+1)}$ results in Eq. (42).

Similarly, taking the partial derivative of the $Q$-function w.r.t. ${{c}}$ to get

[TABLE]

Equating to zero, we obtain the equation in (43).

Appendix E: Fitted Functional Coefficients for the Mixture MLE Model

The fitted coefficient values in Table I are obtained for the functional fits in (46) by comparing our Mixture model against simulation data. These coefficients are estimated for a maximum transmission power of $P_{s}=30$ dBm. To estimate the mixture distribution parameters at a different transmission power $P_{new}$ dBm, we only need to update the $b,\lambda$ and $\mu$ parameters by scaling each of them by the corresponding transmit power scaling factor, i.e. $b_{new}=10^{0.1(P_{new}-P_{s})}\times b$. Parameter $c$ stay unchanged with respect to the transmission power.