# Efficient Numerical Methods for Analysis of Square Ratio ofκ-μ and η-μ Random Processes with Their Applications in Telecommunications

Authors: Gradimir V. Milovanović; Suad Suljović; Stefan R. Panić; Ismet Kalčo; Mihajlo Č. Stefanović
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

We will provide statistical analysis of the square ratio ofκ-μ and η-μ random processes and its application in the signal-to-interference ratio (SIR) based performance analysis of wireless transmission subjected to the influence of multipath fading, modelled by κ-μ fading model, and undesired occurrence of co-channel interference (CCI), distributed as η-μ random process. First contribution of the paper is deriving exact closed expressions for the probability density function (PDF) and cumulative distribution function (CDF) of square ratio of κ-μ and η-μ random processes. Further, a verification of accuracy of these PDF and CDF expressions was given by comparison with the corresponding approximations obtained by the high-precision quadrature formulas of Gaussian type with respect to the weight functions on (0,+∞). The computational procedure of such quadrature rules is provided by using the constructive theory of orthogonal polynomials and theMathematica package OrthogonalPolynomials created by Cvetković and Milovanović (2004). Capitalizing on obtained expression, important wireless performance criteria, namely, outage probability (OP), have been obtained, as functions of transmission parameters. Also, possible performance improvement is observed through a glance at SC (selection combining) reception employment based on obtained expressions.


---

## Body

## 1. Introduction

LetX1i and X2i be mutually independent Gaussian random processes, with variances V(X1i)=V(X2i)=σ2. Modelling mean values of Gaussians, with pi and qi, and assuming μ initially to be an integer, the resulting random process x is obtained as the following function of Gaussians: (1)x=∑i=1μX1i+pi2+∑i=1μX2i+qi2which follows κ-μ distribution, whose PDF is given in the following form (cf. [1]): (2)fXx=2μ1+κμ+1/2xμκμ-1/2eμκΩμ+1/2exp⁡-μ1+κx2ΩIμ-12μκ1+κx2Ω,with (3)p2=∑i=1μpi2,q2=∑i=1μqi2,κ=p2+q22μσ2,Ω=2μσ21+κ,and I0(x) is the zero-the order modified Bessel function of the first kind [2, Eq. 8.445].As a general probability distribution model, this model includes some classical distribution models as its particular cases, for example, Rician, and Nakagami-m distribution models as special cases (as the One-Sided Gaussian and the Rayleigh distributions since they also represent special cases of Nakagami-m).LetY1i and Y2i, i=1,…,μ, be Gaussian random processes, with E(Y1i)=E(Y2i)=0. Assuming μ initially to be an integer, resulting random process y obtained as a function of Gaussians (4)Y=∑i=1μY1i2+Y2i2follows η-μ distribution, whose PDF is given in the following form (cf. [1]):(5)fYy=4πμμ+1/2hμy2μΓμHμ-1/2Ωμ+1/2exp⁡-2μhy2ΩIμ-1/22μHy2Ω,with Ω defined as Ω=E(y2), the ratio of arbitrary variances of independent Gaussians defined as η=E(Y1i2)/E(Y2i2), and distribution parameters H and h defined as H=(η-1-η)/4, h=(2+η-1+η)/4,  η≥0. As a general probability distribution model, this model also includes some classical distribution models as its particular cases, for example, Nakagami-q (Hoyt), one-sided Gaussian, Rayleigh, and Nakagami-m.Multipath fading is physical phenomena that occurs as randomly delayed components of desired signal combine in constructive or destructive manner at the reception [1]. In the case when a line of sight (LOS) component between the transmitter and the receiver is present, the most general fading model which describes the short-term signal variation is kappa-mu (κ-μ) fading model [3, 4]. Fading κ-μ model, defined as the function of parameter κ, related to the dominant/scattered components powers quotient, and parameter μ, related to the number of propagating clusters, easily reduces to other fading models, by setting corresponding values for parameters κ and μ. When μ=1, the κ-μ model transforms to the Rice model with arbitrary K factor, K=κ [4]. By assigning κ=0, the observed κ-μ model transforms to the Nakagami-m model [4]. Co-channel interference (CCI) signal which is transmitted at the same frequency as the carrier signal is also exposed to multipath fading phenomenon [5]. Due to complicated propagation of interfering signal, CCI at the reception is usually modelled with random process that describes small scale signal variations in general non-line-of-sight condition. Recently proposed [6] eta-mu (η-μ) fading distribution meets those conditions. As mentioned in [7], η-μ model is presented in function of the parameter η, related to the in-phase/in-quadrature components scattered wave powers quotient of each cluster, and parameter μ, related to the number of multipath clusters in the environment. As a general distribution, this model spans through other some well-known fading models as its particular cases, that is, Hoyt model and Nakagami-m model. Namely, Nakagami-m distribution can be obtained in an exact manner by assigning parameter values μ=m and η→0 (cf. [6]).Observing interference-limited fading environment, where corruptive effects of thermal noise can be ignored, it is necessary to determine properties of signal-to-interference ratio (SIR) in order to carry out effective performance analysis of observed system. Indeed, SIR has been measured in base and mobile stations by using SIR estimators [8]. In order to determine behavior of instantaneous SIR random process, we must determine square ratio of κ-μ and η-μ random processes, λ=X2/Y2, and analyze its properties. Here, exact closed expressions for the probability density function (PDF) and cumulative distribution function (CDF) of square ratio of κ-μ and η-μ random processes are shown. Further, a verification of accuracy of these PDF and CDF expressions was given by comparing with the corresponding approximations obtained by the high-precision quadrature formulas of Gaussian type with respect to the weight functions on (0,+∞). The computational procedure of such quadrature rules is provided by using Gautschi’s constructive theory of orthogonal polynomials (see [9, 10]) and the Mathematica package OrthogonalPolynomials created by Cvetković and Milovanović (see [11, 12]). Based on obtained expressions, rapidly converging infinite-series expression for wireless communication systems performance measure, outage probability (OP) is also presented, and OP is considered for different values of transmission parameters. At the end, possible performance improvement will be considered through a glance at SC (selection combining) reception employment.The paper is organized as follows. System model and closed-form of the PDF and CDF expressions are presented in Section2. Section 3 is devoted to an alternative approach based on the constructive theory of orthogonal polynomials and the corresponding quadrature rules of Gaussian type. System performances and numerical results are presented in Sections 4 and 5, respectively. Finally, concluding comments are given in Section 6.


## 2. System Model and Closed-Form of the PDF and CDF Expressions

The desired signal followsκ-μ distributed random process modelled as (cf. [3]) (6)pRR=CRRμdexp⁡-μd1+κdR2ΩdIμd-12μdκd1+κdR2Ωd,where (7)CR=2μd1+κdμd+1/2κdμd-1/2eμdκdΩdμd+1/2,with Ωd=E[R2], denoting desired signal average power, while Iν(x) being the νth order modified Bessel function of the first kind [2, Eq. 8.445]. Due to physical separation between the co-channel interferer and target receiver, which is very large in practice, it is very unlikely that direct LOS condition (s) exist, for CCI signal. The most general fading model, which does not include LOS condition is η-μ fading model, is chosen here for modelling random envelope fluctuations of CCI signal with corresponding PDF [6]: (8)prr=4πμcμc+1/2hcμcr2μcΓμcHcμc-1/2Ωcμc+1/2exp⁡-2μchcr2ΩcIμc-1/22μcHcr2Ωc,where Ωc=E[R2] stands for the average power of CCI, while Γ(a) denotes Gamma function [2, Eq. 8.310.1]. Here Hc and hc are interfering signal parameters, written in the function of parameter ηc as follows [6]: (9)Hc=ηc-1-ηc4,hc=2+ηc-1+ηc4.Now, we will observe random envelope statistics in observed interference-limited system. Here the novelty and significance of this model must be pointed out again. Namely, for the first time scenario is observed when desired signal envelope variations are caused byκ-μ fading occurrence (most general fading model which includes LOS component existence), while CCI signal envelope variations are caused by η-μ fading occurrence (most general fading model that does not include LOS component). Previously, various other scenarios have been observed in [1, 13–15]. In [1, 14] LOS component existence in CCI signal has been assumed, which is not so common case, because of large physical distance between the interferer and receiver. In [13] non-LOS model was observed for desired signal, while Gamma approximation of η-μ random variables sum has been observed. Finally, neither in general nor in comprehensive study, such as [15], has this model been considered. In the case when CCI is much stronger than Gaussian noise in same channel, PDF for the instantaneous SIR, λ=R2/r2, can be determined according to [16] in the form (10)pλt=12t∫0+∞pRrtprrrdr.After substituting (6) and (8) into (10), we can derive PDF of SIR in the following closed-form:(11)pλλ=2πhcμcexp⁡μdκdΓμc∑p=0∞∑l=0∞Aplμd2p+μd1+κdp+μdκdpμc2l+2μcHc2l×S/λ2l+2μc+1Sμd1+κd+2S/λμchcp+2l+2μc+μd,where (12)Apl=Γp+2l+2μc+μdΓμd+pΓμc+l+1/2p!l!and S=Ωd/Ωc stands for the average SIR.Namely, with respect to well-known series representation of modified Bessel function (cf. [2, Eq. 8.445]) (13)Iνz=∑k=0+∞1k!Γk+ν+1z22k+ν,substituting (6) and (8) in (10) results in (14)pλλ=2πhcμcexp⁡μdκdΓμc∑p=0∞∑l=0∞μd2p+μd1+κdp+μdκdpμc2l+2μcHc2lΓμd+pΓμc+l+1/2p!l!×S/λ2l+2μc+1Sμd1+κd+2S/λμchcp+2l+2μc+μd×∫0∞tp+2l+2μc+μd-1exp⁡-tdt.By using the well-known definition of the Gamma function [2, Eq. 8.310.1]: (15)Γa=∫0∞ta-1exp⁡-tdt,(10) can be reduced into (11). A similar procedure has been used in [7]. By taking finite numbers of terms p=l=m in the last expression for pλ(λ) we obtain its approximation in a closed-form: (16)pλλm=2πhcμcexp⁡μdκdΓμc∑p=0m∑l=0mApl×μd2p+μd1+κdp+μdκdpμc2l+2μcHc2lS/λ2l+2μc+1Sμd1+κd+2S/λμchcp+2l+2μc+μd.The closed-form expression (16) converges rapidly, since only m  (10–15) terms are needed to be summed in each sum in order to reach an accuracy of five significant decimal digits, for observed set of system parameters. Otherwise, in a case of slowly convergent series, we need some methods for accelerating their convergence (cf. [17–19]).Following [16], the cumulative distribution function (CDF) can be obtained as (17)Fλt=∫0tpλxdx.After substituting (11) into (17), by applying the same mathematical transformations given in Appendix of [7], we can obtain the following closed-form expression of rapid convergence: (18)Fλt=πeμdkdΓμc∑p=0∞∑l=0∞μdpκdpHc2l22l+2μc-1hc2l+μcQp,l,where (19)Qp,l=Γp+2l+2μc+μdBzp+μd,2l+2μcΓμd+pΓμd+l+1/2p!l!,z=λμd1+κdλμd1+κd+2Sμchc,and Bz(a,b) is the incomplete Beta function [2, Eq. 8.38] (for some extended special functions see [20]). Detailed convergence analysis of this expression is provided in Table 1.Table 1
Number of terms that should be added in expression (18) in order to reach accuracy at 5th significant digit.
μd μ c κ d κ c S = - 10 dB S = 0 dB S = 10 dB 1.5 1.5 2.0 1.5 22 22 21 2.5 2.5 2.0 1.5 21 21 20 2.5 2.5 1.5 1.5 21 20 19 2.5 2.5 1.5 0.5 20 19 19


## 3. Method of Quadrature Sums

In this section we give an alternative approach for getting very high-precision approximation of the probability density function (PDF), as well as one of the cumulative distribution function (CDF). Our goal is to obtain fast convergence of quadrature sums,Qn(f;w)=∑k=1nAk(n)f(xk(n)), for integrals on half-infinite ranges or finite intervals for functions of the form F(x)=w(x)f(x), where w(x)(≥0) is the so-called weight function and f(x) is a smooth factor if it is possible. We use here a computational procedure for constructing quadratures of Gaussian type with respect to the weight functions w(x) on (0,+∞), as well as ones on finite intervals, supposing the polynomial behavior of the function f(x) (for nonpolynomial cases see [21, 22], and for some nonstandard Gaussian quadratures see [23–25]). Our computational procedure provides the ability to verify the accuracy of expressions obtained in the previous section. As we will see, these quadrature processes give much better accuracy of the considered integrals.As already mentioned PDF of this random process could be efficiently expressed by using Gaussian quadrature formulas, namely, starting from (10), which can be rewritten in the following form: (20)pλλ=A∫0+∞xμc+μd/2e-xIμc-1/2αcxIμd-1αdxdx,with parameters αc, αd, and A depending on (being functions of) normalized SIR (z=S/λ) expressed as (21)αcz=2μcHczμd1+κd+2μchcz,αdz=2μdκd1+κdμd1+κd+2μchcz,Az=2πμcμc+1/2hcμcμd1+κdμd+1/2zμc+3/2Hcμc-1/2Γμcκdμd-1/2eκdμdμd1+κd+2μchczμd/2+μc+1S.For calculating the weighted integral (20) we use the n-point generalized Gauss-Laguerre quadrature formula (cf. [26, pp. 324–328]) (22)∫0+∞xse-xfxdx=∑k=1nAknfxkn+Rnf,with respect to the generalized Laguerre weight function w(x)=w(s)(x)=xse-x on (0,+∞), with a special selection of the parameter s. The quadrature parameters, the nodes xk(n)≡xk, and the weights (Christoffel numbers) Ak(n)≡Ak can easily be constructed by the well-known Golub-Welsch procedure [27] implemented in the Mathematica package OrthogonalPolynomials (see [11, 12]). This package is freely downloadable from the website http://www.mi.sanu.ac.rs/~gvm/.In fact, the nodesxk are eigenvalues of the following symmetric tridiagonal Jacobi matrix (cf. [26, pp. 325–328]): (23)Jn=α0β1Oβ1α1β2β2α2⋱⋱⋱βn-1Oβn-1αn-1,and the weight coefficients Ak are given by Ak=β0vk,12, k=1,…,n, where vk,1 is the first component of the eigenvector vk=[vk,1⋯vk,n]T corresponding to the eigenvalue xk and normalized such that vkTvk=1.Remark 1.
In some special cases for calculating the weight coefficientsAk it is better to use the complete eigenvectors, instead of their first components. An analysis of such cases is given in [28] (see also [29, p. 339] for an application).Elements in the Jacobi matrix (23) are recursion coefficients in the three-term recurrence relation (24)πk+1t=t-αkπkt-βkπk-1tk=0,1,…,π0t=1,π-1t=0,for polynomials orthogonal with respect to the generalized Laguerre weight function w(x)=w(s)(x) on (0,+∞): (25)αk=2k+s+1k=0,1,…,βk=kk+sk=1,2,….The coefficient β0 can be arbitrary, but it is convenient to put (26)β0=∫0+∞wsxdx=Γs+1,in order to hold certain compact formulas, including the previously mentioned ones for Christoffel numbers Ak, k=1,…,n.Remark 2.
The recursion coefficients in the three-term recurrence relation (24) are known in explicit form only for some narrow classes of orthogonal polynomials, among which are the most famousclassical orthogonal polynomials (Jacobi, the generalized Laguerre, and Hermite polynomials) (cf. [26, pp. 121–146], [30], [31, pp. 19–24]). In the case of the so-calledstrongly nonclassical polynomials, the recursion coefficients must be constructed numerically using moment information on a given weight function (cf. [32–34]). For details onGautschi’s constructive theory of orthogonal polynomials see [9, 10, 35], [26, pp. 160–166].Knowing the firstN recursion coefficients in (24), we can obtain the nodes xk and the weights Ak in (22) for any number of nodes n≤N.It is known that if the integrand in the weighted integral (20) is smooth, then we can expect Gaussian quadrature applied to (20) to converge rapidly. According to the series representation (13), we see that x↦x-νIν(x) is an entire function, and, therefore, the function (27)x⟼gx;z=AzIμc-1/2αczxxμc-1/2·Iμd-1αdzxxμd-1/2is also entire. This means that, for calculating (20), this integral should be written in the following form: (28)pλλ=∫0+∞x2μc+μd-1e-xgx;zdx,where z=S/λ and g(x;z) is defined by (27), and then quadrature formulas are applied with respect to the generalized Gauss-Laguerre weight with the parameter s=2μc+μd-1.Then the sequence of quadrature sums{pλ(n)(λ)}n, defined by (29)pλnλ=∑k=1nAkngxkn;z,shows a rapid convergence.For example, for some typical parameters(30)Ωd=1,μd=2,κd=1,ηc=12,ρc=1,μc=2.01,Hc=ηc-1-ηc4,hc=2+ηc-1+ηc4,S=ΩdΩc,and λ∈(0,3), we present the graphic of λ↦pλ(20)(λ) in Figure 1 as an approximation of the function λ↦pλ(λ).Figure 1
Approximation of PDF,λ↦pλ(λ)≈pλ(n)(λ), using n=6 nodes in (29).Relative errors(31)errnλ=pλnλ-pλλpλλare presented in log-scale in Figure 2, for five different numbers of nodes in the quadrature sums (n=4, 6, 8, 10, and 12). As the exact PDF we take the approximation pλ(λ)≈pλ(20)(λ). For example, we can see that the quadrature sum with only six points gives an approximation of PDF for each λ>0, with the relative error less than 10-6; that is, pλ(6)(λ) gives the result with at least six decimal digits, but the 10-point quadrature sum approximates the PDF with at least eleven digits.Figure 2
Relative errors inpλ(n)(λ) obtained by the quadrature formula with n=4(2)12 nodes.According to (17) the cumulative distribution function (CDF) can be transformed as (32)Fλt=t∫01pλtτdτ.For numerically calculating this integral we use the Gauss-Legendre quadrature rule (transformed to the interval [0,1]), (33)∫01fτdτ=∑j=1NBjNfτjN+RNf,whose parameters, the nodes and weights, τj(N) and Bj(N), j=1,…,N, can be obtained from the eigenvalue problem for the corresponding Jacobi matrix of the form (23), with recursion coefficients (34)αj=12j=0,1,…,β0=1,βj=j244j2-1j=1,2,….In this symmetric case (with respect to τ=1/2), it is clear that τN+1-j(N)=1-τj(N) and BN+1-j(N)=Bj(N), j=1,2,…,N. For a sufficiently smooth function f∈C2N[0,1], the remaining term in (33) has the form [26, p. 333 & p. 97] (35)RNf=β0β1⋯β2N2N!f2Nζ0<ζ<1.As before, for approximatingpλ(x), we use (29), so that we have (36)FλN,nt=t∑j=1NBjNpλntτjN,that is, (37)FλN,nt=t∑j=1NBjN∑k=1nAkngxkn;StτjN;where g(x;z) is defined by (27).The corresponding relative errors(38)errN,nt=FλN,nt-FλtFλtare presented in log-scale in Figure 3, for N=n=8 and 12 in the quadrature sums (36), that is, (37). As the exact CDF we used the approximation Fλ(t)≈Fλ(20,20)(t).Figure 3
Relative errors in the quadrature sumsFλ(N,n)(t) for N=n=8 and 12 nodes.As we can see, in the case whent is large, this formula needs a bigger N, and therefore the following approach is much more convenient. Namely, we use a decomposition of the interval [0,t] into several subintervals for larger values of t, for example, (39)0,t=0,1∪1,2∪⋯∪ν-1,ν∪ν,tν≤t,and then apply the Gauss-Legendre rule (33) to each of the subintervals. Thus, (40)Fλt=∑i=1ν∫i-1ipλxdx+∫νtpλxdx=∑i=1ν∫01pλξ+i-1dξ+t-ν∫01pλν+t-νξdξ=∫01∑i=1νpλξ+i-1+t-νpλν+t-νξdξ,where the sums on the right-hand side are empty for t<1(ν=0). The corresponding quadrature approximation, analogous to (36), becomes (41)F^λN,nt=∑j=1NBjN∑i=1νpλnτjN+i-1+t-νpλnν+t-ντjN.Takingν+1 subintervals (ν=0,1,2,3), the relative errors in the quadrature sums F^λ(N,n)(t) for N=n=6 are presented in Figure 4 for 0≤t≤10. We see that the last case, with four subintervals [0,1], [1,2], [2,3], and [3,t](t≥3), gives the relative errors less than 10-6 for each t∈[0,10]. The corresponding cumulative distribution function is presented in Figure 5.Figure 4
Relative errors in the quadrature sumsF^λ(N,n)(t) for 0≤t≤10, taking ν+1 subintervals (ν=0,1,2,3) and N=n=6.Figure 5
Approximation of CDF,t↦Fλ(t)≈F^λ(N,n)(t), for N=n=6 with four subintervals (ν=3).


## 4. System Performances

One of standard performance criteria, used for managing CCI level in faded channel, and for meeting grade-of-service and quality-of-service requirements, is OP [36]. When considering cellular systems, OP is defined as the probability that received SIR falls below a given outage threshold qth also known as a protection ratio (42)Pout=PRσ<qth=∫0qthpσtdt=Fσqth.Diversity combining indicates a concept of upgrading reception system performances by combining replicas of information carrying signal [37]. By using separated antennas on reception diversity terminal, independent copies of desired signal are accomplished, without increase in signal bandwidth or transmitted signal power. Simplest space diversity combining method is selection combining (SC). In communication environments where the level of the CCI is significantly higher than the level of channel AWGN, SC outputs the antenna with the highest SIR. When reception system has uncorrelated branches (sufficient spacing between antennas), the PDF of the SIR at the output of the receiver is given by [38]: (43)pλt=∑j=1Npλjt∏i=1,i≠jNFλkt,where N denotes for the diversity reception order.


## 5. Numerical Results

At Figure6 OP in the function of normalized SIR in dB (decibels, 10log⁡(S/λ)) is presented. It is visible how smaller values of OP are obtained when parameter κd arises, since dominant component of desired signal then straightens. Also better performance (smaller values of OP) is obtained when fading, so that influences desired and CCI signal are less severe, namely, when parameters μd and μc take higher values. With SC reception performance improvement reached is also shown. Considerably lower OP values are obtained, comparing to the case with no-diversity reception. Finally it can be seen how balanced reception leads to achieving better quality of service. For example, it can be seen from the figure that when SC diversity has been used, for reaching the same level of OP (i.e., 10-4), it is necessary to apply about 4 dB lower level power at the transmission, comparing to the level of transmitting power necessary to reach given OP level for the case of single branch reception when same values of propagation parameters are considered. Also, it can be seen from the figure that when the same level of transmitted signal power exists (i.e., 12 dB) higher level of outage (i.e., 10-3) is obtained for single branch reception compared to case when SC diversity has been applied (i.e., 10-4), when same values of propagation parameters are considered. Presented efficient numerical methods in such way could enable wireless system designers to perform trade-off studies among the required values of system performances at the reception (quality of service) in order to determine the optimal choice of transmission parameters in the presence of their available constraints.Figure 6
Outage probability versus normalized SIR.


## 6. Conclusion and Future Work

In this paper the exact PDF and CDF expressions of the square ratio ofκ-μ and η-μ random processes were derived. Also, a verification of accuracy of these PDF and CDF expressions was given by comparing with the corresponding approximations obtained by the high-precision quadrature formulas of Gaussian type with respect to the weight functions on (0,+∞). The computational procedure of such quadrature rules is provided by using the contemporary constructive theory of orthogonal polynomials and the corresponding software. Our efficient numerical methods provide necessary control of level of error for the obtained approximation expressions. Also, the convergence of the proposed expressions is secured. Based on the provided analysis of this novel random process, interference-limited wireless communication subjected to the influence of κ-μ fading was observed when undesired occurrence of CCI has been modelled with η-μ random process. Rapid converging infinite-series expressions are presented for received SIR’s PDF and CDF. Major communication performance measure, OP, has been calculated in terms of transmission parameters. OP upgrade accomplished by SC reception usage has also been analyzed.Contribution presented in this paper provides potentials for applying proposed efficient numerical methods in carrying out analysis of other important wireless telecommunication performance measures for observed propagation scenario, such as channel capacity (CC) and average bit error probability (ABER). Closed-form expressions for ABER values could be derived when both coherent and noncoherent modulation techniques are applied for wireless propagation overη-μ interference-limited channels. In a similar manner CC values could be efficiently evaluated for proposed propagation scenario. By using provided high-precession approximation techniques that are based on presented quadrature formulas of Gaussian type, a verification of accuracy of these ABER and CC expressions could also be provided. Finally, presented approximations could find their application in many fields of telecommunication theory, since proposed solutions can be applied for performance analysis of large number of wireless transmission scenarios.



---

*Source: 4967613-2018-03-11.xml*