Student Member, IEEE Kezhi Wang
email: kezhi.wang@warwick.ac.uk

Senior Member, IEEE Yunfei Chen
email: yunfei.chen@warwick.ac.uk

Member, IEEE Marco Di Renzo
email: marco.direnzo@lss.supelec.fr

Outage Probability of Dual-Hop Selective AF With Randomly Distributed and Fixed Interferers

Keywords: Amplify-and-forward, interference, outage probability, Poisson point process, relay selection

The outage probability performance of a dual-hop amplify-and-forward selective relaying system with global relay selection is analyzed for Nakagami-m fading channels in the presence of multiple interferers at both the relays and the destination. Two different cases are considered. In the first case, the interferers are assumed to have random number and locations. Outage probability using the generalized Gamma approximation (GGA) in the form of one-dimensional integral is derived. In the second case, the interferers are assumed to have fixed number and locations. Exact outage probability in the form of one-dimensional integral is derived. For both cases, closed-form expressions of lower bounds and asymptotic expressions for high signal-to-interference-plus-noise ratio are also provided.

Simplified closed-form expressions of outage probability for special cases (e.g., dominant interferences, i.i.d. interferers, Rayleigh distributed signals) are studied. Numerical results are presented to show the accuracy of our analysis by examining the effects of the number and locations of interferers on the outage performances of both AF systems with random and fixed interferers.

I. INTRODUCTION

Wireless relaying can extend the network coverage by using idle nodes as relays in the network.

It can also provide diversity gain by using idle nodes as "virtual" antennas [START_REF] Laneman | Energy-efficient antenna sharing and relaying for wireless networks[END_REF]. Consequently, a huge amount of works have been conducted on its application in future wireless networks. Among all the relaying strategies, amplify-and-forward (AF) and decode-and-forward (DF) are perhaps the most widely used ones. In AF relaying, the source broadcasts its information to the relays in the first phase and then the relays simply amplify the received signals from the source and forward them to the destination in the second phase, while in DF relaying, the source broadcasts its information in the first phase but the relays have to decode the received signals from the source and then re-encode the signals before forwarding them to the destination. Due to its lower complexity, AF relaying is more attractive for some applications [START_REF] Hasna | A performance study of dual-hop transmissions with fixed gain relays[END_REF]. On the other June 12, 2018 DRAFT hand, in practical systems, it is often the case that more than one idle nodes are available at the same time such that multiple relays can be used by the source. If all the idle nodes are used in wireless relaying, orthogonal channels between relays have to be used in the second relaying phase such that the relayed signals will not interfere with each other. This will reduce the usage efficiency of the system resources considerably. To solve this problem, relay selection can be used by choosing only one node out of all available idle nodes in the second phase [START_REF] Chen | Novel partial selection schemes for af relaying in Nakagamim fading channels[END_REF]. It can be shown that relay selection can achieve the same diversity gain as the scheme that uses all available idle nodes, with proper designs [START_REF] Jing | Single and multiple relay selection schemes and their achievable diversity orders[END_REF]. Thus, this paper focuses on relay selection using AF as an effective technology to achieve reliable communications.

Several researchers have studied the performance of relay selection using AF. In [START_REF] Zhao | Improving amplify-and-forward relay networks: optimal power allocation versus selection[END_REF], the optimal relay selection criterion was proposed by selecting the relay with the largest instantaneous end-to-end or global signal-to-noise ratio (SNR) for forwarding. The performance of this criterion was analyzed in [START_REF]Symbol error rate of selection amplify-and-forward relay systems[END_REF]. In [START_REF] Bletsas | A simple cooperative diversity method based on network path selection[END_REF], two suboptimal relay selection schemes based on two upper bounds to the instantaneous global SNR were proposed and analyzed. Reference [START_REF] Krikidis | Amplify-and-forward with partial relay selection[END_REF] proposed partial relay selection scheme where relay selection is based on only the instantaneous SNR of the first hop. In [START_REF] Chen | Novel partial selection schemes for af relaying in Nakagamim fading channels[END_REF], both the optimal selection scheme and the partial selection scheme were analyzed for Nakagami-m fading channels. Two new partial relay selection schemes were also proposed in [START_REF] Chen | Novel partial selection schemes for af relaying in Nakagamim fading channels[END_REF]. However, none of these works considered the interferences from other transmitting sources in the network. In a multiple-access system or a frequency-reused cell, interferences from other transmitting sources, such as interferers, may cause performance degradation and therefore, cannot be ignored. Moreover, the positions of the nodes may not be optimized such that interferers may be randomly distributed. In this case, the spatially averaged (over the distributions of the positions) performance metrics may be of more practical use for system design and optimization by considering random locations of interferers. Reference [START_REF] Guidotti | Outage and symbol error probabilities of dual-hop af relaying in a poisson field of interferers[END_REF] provided the closed-form expression of the outage probability of dual-hop AF relaying in the presence of interference at the destination over Rayleigh fading channels. Reference [START_REF] Suraweera | Performance analysis of fixed gain relay systems with a single interferer in Nakagami-m fading channels[END_REF] analyzed the outage probability of a dual-hop AF relaying system where both relay and destination are interfered by a single source in Nakagami-m fading channel. However these two works either considered interferers only at the destination over Rayleigh fading channels or only a signal interferer of fixed location.

In this paper, we provide a comprehensive analytical framework to derive the outage probability performance of an AF relay selection system where the relays suffer from path loss, independent but non-identically distributed Nakagami-m fading as well as Nakagami-m interferences. In the first case, the interferers have random number and locations. This is the case for multipleaccess systems with mobile nodes. In the second case, the interferers have fixed number and locations. This is the case for fixed-access wireless systems where wireless interconnections are mainly provided to replace wires with considerably low or little mobility. The optimal criterion that selects the relay with the largest instantaneous global signal-to-noise-plus-interference ratio (SINR) is studied. The exact outage probability in terms of either one-dimensional integrals or closed-form approximations are derived. Also, lower bounds to outage probability are given.

Finally, asymptotic expression of outage probability in Rayleigh fading are studied for large SINR values. Numerical results are presented to show the accuracy of our analysis and therefore to examine the effects of interferences on relay selection using AF.

The remainder of this paper is organized as follows. Section II introduces the system model.

Section III considers the case when the interferers have random number and locations, while

Section IV studies the case when the interferers have fixed number and locations. Numerical results are presented in Section V, followed by concluding remarks in Section VI.

II. SYSTEM MODEL

Consider a wireless relaying system with one source S, one destination D and multiple relays

J of R j , j = 1, 2, • • • , J.
There is no direct link between the source and the destination, which is the case when relays are used to extend network coverage and is the focus of this paper. All nodes have a single antenna and are in half-duplex mode. Assume that there are I sj interferers of

I sj i , i = 1, 2, • • • , I sj , and I jd interferers of I jd v , v = 1, 2, • • • , I jd ,
that are transmitting at the same time as the source to the j-th relay and the j-th relay to the destination, causing interferences to June 12, 2018 DRAFT the j-th relays and the destination, respectively. Assume that the distance between the source S and the j-th relay R j , the j-th relay R j and the destination D, the i-th interferer I sj i and the j-th relay R j and the v-th interferer I jd v and the destination D are l sj , l jd , l ij and l vj , respectively.

Also, assume that the path loss between S and R j , R j and D, I sj i and R j , I jd v and D are η(l sj ), η(l jd ), η(l ij ) and η(l vj ), respectively. As the singular path loss model leads to impractical power condition in the network when l < 1, we assume the non-singular model for the path loss as

η(l) = 1 l β + 1 ( 1
)
where β is the path loss exponent.

In the case when the interferers have random number and locations, we assume the numbers I sj , I jd and the distances l ij , l vj are random. We assume a Poisson point process (PPP) with density λ I for the spatial distribution of the interferers. Then, the probability density function (PDF) of the number of interferers is given as

P r {i = I} = (λ I A I ) I I! e -λ I A I , I = 0, 1, • • • (2)
where A I is the distribution area of interferers. Also, we assume that the distance l follows a general distribution with a PDF of f l (x) which can be specified for different applications considered. In the case when the interferers have fixed locations, both the number of interferers and their locations are fixed such that I sj , I jd , l ij and l vj are deterministic values.

The received signal from the source S to the j-th relay R j is given by

y sj = Ω sj h sj x + I sj i=1 Ω ij h ij x ij + n sj (3)
and it can be further amplified and forwarded such that the received signal at the destination is given by

y jd = Ω jd h jd • G j • y sj + I jd v=1 Ω vj h vj x vj + n jd (4)
where the source S and the j-th relay R j , the j-th relay R j and the destination D, the i-th interferer I sj i and the j-th relay R j , the v-th interferer I jd v and the destination D, respectively, P sj , P jd , P ij and P vj are the transmitted power of S, R j , I sj i and I jd v , respectively, K sj , K jd , K ij and K vj are constants that take other power factors, such as antenna gains and the average signal power factors, into account, h sj , h jd , h ij and h vj are the fading gains with unit average power between S and R j , R j and D, I sj i and R j , I jd v and D, respectively, x, x ij and x vj are the transmitted symbol of S, I sj i and I jd v , respectively, n sj , n jd are the additive white Gaussian noise (AWGN) in the channel between S and R j , R j and D, respectively, and G j is the amplification factor.

Ω sj = K sj P sj |h sj | 2 η(l sj ), Ω jd = K jd P jd |h jd | 2 η(l jd ), Ω ij = K ij P ij |h ij | 2 η(l ij
In the above, assume enough distances between relays and between source and relays and between relays and destination such that |h sj |, j = 1, 2, • • • , J, are independent of each other, we assume that n sj are independent of each other for different j and n jd are independent of each other for different j and n sj are independent of n jd . Also, assume enough distances between interferers at the relay and between interferers at the destination such that

|h jd |, j = 1, 2
I sj i , |h ij |, i = 1, 2, • • • , I sj ,
at R j are independent of each other, respectively, and

I jd v , |h vj |, v = 1, 2
= E{Y sj } σ 2 sj , γ IN R jd = E{Y jd } σ 2 jd
, respectively, and the SINR between S and R j , R j and

D as γ SIN R sj = E{Ω sj |h sj | 2 } σ 2 sj +E{Y sj } , γ SIN R jd = E{Ω jd |h jd | 2 } σ 2
jd +E{Y jd } , respectively. Using these assumptions, for variable-gain relaying, the amplification factor is given by

G j = 1 Ω sj |h sj | 2 + σ 2 sj + I sj i=1 Ω ij |h ij | 2 . ( 5
)
Using ( 5) in ( 4), the instantaneous end-to-end SINR of the j-th relaying link can be derived as

Γ j = Γ sj Γ jd Γ sj + Γ jd + 1 (6)
where

Γ sj = Ω sj |h sj | 2 σ 2 sj + Y sj (7)
and

Γ jd = Ω jd |h jd | 2 σ 2 jd + Y jd . (8)
In relay selection, the relay with the largest end-to-end instantaneous SINR is selected. Thus, the outage probability for a threshold of γ th is given by

P o (γ th ) = P r{max{Γ j } < γ th }, j = 1, 2, • • • , J (9)
where γ th = 2 2R -1 and R is the transmission rate. In the next sections, this outage probability is derived in different cases.

III. RANDOM INTERFERERS

We first consider the case when the interferers have random number and locations. In this case, the randomness comes from the Nakagami-m fading powers, the number of interferers I sj and I jd as well as the distances l ij and l vj . In the first subsection, the PDFs of Y sj and Y jd are derived using the generalized Gamma approximation (GGA), as their closed-form expressions are

     B d sj +1 p sj , d sj +2 p sj = E{Y sj }E{Y 2 sj } E{Y 3 sj } B d sj p sj , d sj +3 p sj B d sj +1 p sj , d sj +1 p sj = E{Y sj } 2 E{Y 2 sj } B d sj p sj , d sj +2 p sj (11)
where E{Y sj }, E{Y 2 sj } and E{Y 3 sj } can be found in ( 44), ( 45) and ( 46), respectively, and then one can calculate a sj by inserting the solved values of p sj and d sj into one of the equations in (47). Also, when all the interferers are i.i.d., with the help of

e λ = ∞ I=0 λ I I! , e λ = 1 λ ∞ I=1 I λ I I! or e λ = 1 λ 2 ∞ I=1 I(I -1) λ I I! [13], E{Y sj }, E{Y 2 sj } and E{Y 3 sj } can be simplified as                  E{Y sj } = I sj K ij P ij η 1,ij E{Y 2 sj } = I sj K 2 ij P 2 ij η 2,ij m ij +1 m ij + (I sj ) 2 K i 1 j P i 1 j K i 2 j P i 2 j η 2 1,ij E{Y 3 sj } = I sj K 3 ij P 3 ij (m ij +1)(m ij +2) m 2 ij η 3,ij + (I sj ) 2 K 2 i 1 j P 2 i 1 j K i 2 j P i 2 j m ij +1 m ij η 2,ij η 1,ij +(I sj ) 3 K i 1 j P i 1 j K i 2 j P i 2 j K i 3 j P i 3 j η 3 1,ij . (12)
The PDF expression f Y jd (x; a jd , d jd , p jd ) of Y jd can be obtained using the same method as above. They are not listed here due to the limited space.

B. PDF and CDF of Γ sj and Γ jd

Denote W sj = Ω sj |h sj | 2 . Since |h sj | 2 is a Gamma random variable with shape parameter m sj and scale parameter 1/m sj , W sj is also a Gamma random variable with PDF

f W sj (x) = m sj Ω sj m sj x m sj -1 Γ(m sj ) e - m sj Ω sj x , x > 0. (13)
Also, denote

Z sj = Y sj + σ 2 sj , where the PDF of Z sj is determined by f Z sj (x) = f Y sj (x -σ 2 sj ).
Thus, one has

Γ sj = W sj
Z sj and the PDF of Γ sj is given by

f Γ sj (u) = ∞ -∞ |x| f W sj (xu) f Y sj x -σ 2 sj dx. (14)
Using ( 10) and ( 13) in ( 14) and defining p sj = l sj /k sj such that gcd(l sj , k sj ) = 1 where gcd(•, •)

is the great common devisor operator [START_REF] Gradshteyn | Table of Integrals, Series, and Products[END_REF], one has

f Γ sj (u) = m sj r 1 =0 µ 1,sj,r 1 e - m sj σ 2 sj u Ω sj u -d sj +m sj -r 1 -1 G k sj ,l sj l sj ,k sj    µ 0,sj u -l sj | I(l sj , 1 -d sj -r 1 ) I(k sj , 0)    , (15)
DRAFT June 12, 2018 where µ 1,sj,r 1 =

l d sj +r 1 -1 2 sj σ 2(m sj -r 1 ) sj m sj Ω sj -d sj +m sj -r 1 √ k sj p sj m sj !a -d sj sj (2π) - k sj 2 - l sj 2 +1 r 1 !(m sj -r 1 )!Γ(m sj )Γ d sj p sj , µ 0,sj = k -k sj sj l l sj sj a -k sj p sj sj m sj Ω sj -l sj , I(n, ξ) = (ξ/n, (ξ + 1)/n, • • • , (ξ + n -1)/n) and G c,d a,b (•) denotes
the Meijer'G-function [START_REF] Gradshteyn | Table of Integrals, Series, and Products[END_REF] which is available as a built-in function in many mathematical software packages, such as MATLAB, MATHEMATICA and MAPLE.

Proof : See Appendix B.
Also, one can get the CDF of Γ sj as

F Γ sj (u) = 1 - m sj -1 r 2 =0 r 2 r 3 =0 µ 2,sj,r 2 ,r 3 e - m sj σ 2 sj u Ω sj u -d sj +r 2 -r 3 G k sj ,l sj l sj ,k sj   µ0,sju -l sj | I(l sj , 1 -d sj -r 3 ) I(k sj , 0)    (16)
where

µ 2,sj,r 2 ,r 3 = l d sj +r 3 -1 2 sj σ 2(r 2 -r 3 ) sj m sj Ω sj -d sj +r 2 -r 3 √ k sj p sj a -d sj sj (2π) - k sj 2 - l sj 2 +1 r 3 !(r 2 -r 3 )!Γ d sj p sj . Proof : See Appendix C.
The PDF and CDF of Γ sj in ( 15) and ( 16) in terms of the Meijer'G-function are computationally complex. Therefore, we provide the high SINR approximations next that have simpler forms. Using the following Taylor's series expansion

e x = N n=0 x n n! + o(x N ), as x → 0, (17)
where o(x) denotes the higher-order term of an arbitrary function a(x), one can get the PDF and CDF of Γ sj as ( 18) and ( 19), respectively.

f Γ sj (u) = N 1 n 1 =0 N 3 n 3 =0 m sj r 1 =0 µ 5,sj,r 1 ,n 1 ,n 3 u m sj +n 1 +n 3 -1 + o (u/Ω sj ) N 1 + o (σ 2 sj u/Ω sj ) N 3 , (18)
where

µ 5,sj,r 1 ,n 1 ,n 3 = (-1) n 1 +n 3 a n 1 +r 1 sj m m sj +n 1 +n 3 +1 sj Ω -m sj -n 1 -n 3 sj σ 2(m sj -r 1 +n 3 ) sj Γ d sj +n 1 +r 1 p sj n 1 !n 3 !Γ(r 1 +1)Γ(m sj -r 1 +1)Γ d sj p sj , F Γ sj (u) = 1 - N 2 n 2 =0 N 4 n 4 =0 m sj -1 r 2 =0 r 2 r 3 =0 µ 6,sj,r 2 ,r 3 ,n 2 ,n 4 u n 2 +r 2 +n 4 + o (u/Ω sj ) N 2 + o (σ 2 sj u/Ω sj ) N 4 , (19)
where µ 6,sj,r 2 ,r 3 ,n 2 ,n 4 = In high SINR conditions when Ω sj → ∞ or σ 2 sj → 0 or in the low outage regime when u → 0, [START_REF] Al-Qahtani | Exact outage analysis of dual-hop fixed-gain af relaying with cci under dissimilar nakagami-m fading[END_REF] and ( 19) to obtain corresponding approximations, respectively.

one can get rid of o (u/Ω sj ) N 1 , o (u/Ω sj ) N 2 , o (σ 2 sj u/Ω sj ) N 3 and o (σ 2 sj u/Ω sj ) N 4 in
The PDF and CDF expressions f Γ jd (u) and F Γ jd (u) and their high SINR approximations can be obtained using the same method as above. They are not listed here to make the paper compact.

C. Outage probability

Using the derived PDF and CDFs of Γ sj and Γ jd in the previous subsection, the CDF of the instantaneous end-to-end SINR in ( 6) can be derived as [START_REF] Chen | Novel partial selection schemes for af relaying in Nakagamim fading channels[END_REF]]

F E Γ j (x) = ∞ 0 P r Γ sj Γ jd Γ sj + Γ jd + 1 ≤ x|t f Γ sj (t)dt = F Γ sj (x) + ∞ 0 F Γ jd x 2 + x + xt t f Γ sj (t + x)dt. (20)
One can see that [START_REF] Ikki | Performance evaluation and optimization of dual-hop communication over nakagami-m fading channels in the presence of co-channel interferences[END_REF] only has one-dimensional integral, which can be calculated numerically using mathematical software. Also, by using the lower bound as

F Γ j (x) 1 -1 -F Γ sj (x) 1 -F Γ jd (x) = F LB Γ j (x), (21)
one can get the lower bound to the CDF as

F LB Γ j (x) = 1 - m sj -1 r 2 =0 r 2 r 3 =0 m jd -1 r ′ 2 =0 r ′ 2 r ′ 3 =0 µ 2,sj,r 2 ,r 3 µ 2,jd,r ′ 2 ,r ′ 3 e - m sj σ 2 sj x Ω sj - m jd σ 2 jd x Ω jd x -d sj +r 2 -r 3 -d jd +r ′ 2 -r ′ 3 G k sj ,l sj l sj ,k sj    µ 0,sj x -l sj | I(l sj , 1 -d sj -r 3 ) I(k sj , 0)    G k jd ,l jd l jd ,k jd    µ 0,jd x -l jd | I(l jd , 1 -d jd -r ′ 3 ) I(k jd , 0)    , (22)
where (43), one can get the asymptotic expression of the CDF as

µ 2,jd,r ′ 2 ,r ′ 3 = l d jd +r ′ 3 -1 2 jd σ 2(r ′ 2 -r ′ 3 ) jd m jd Ω jd -d jd +r ′ 2 -r ′ 3 √ k jd p jd a -d jd jd (2π) - k jd 2 - l jd 2 +1 r ′ 3 !(r ′ 2 -r ′ 3 )!Γ
F ∞ Γ j (x) = 1 - N 2 n 2 =0 N 4 n 4 =0 m sj -1 r 2 =0 r 2 r 3 =0 N ′ 2 n ′ 2 =0 N ′ 4 n ′ 4 =0 m jd -1 r ′ 2 =0 r ′ 2 r ′ 3 =0
x n 2 +r 2 +n 4 +n ′

2 +r ′ 2 +n ′ 4 (-1) n 2 +n 4 +n ′ 2 +n ′ 4 σ 2(r 2 -r 3 +n 4 ) sj σ 2(r ′ 2 -r ′ 3 +n ′ 4 ) jd m sj Ω sj n 2 +r 2 +n 4 m jd Ω jd n ′ 2 +r ′ 2 +n ′ 4 E Y n 2 +r 3 sj E Y n ′ 2 +r ′ 3 jd n 2 !n 4 !r 3 !(r 2 -r 3 )!n ′ 2 !n ′ 4 !r ′ 3 !(r ′ 2 -r ′ 3 )! (23)
where

E Y n 2 +r 3 sj = a n 2 +r 3 sj Γ d sj +n 2 +r 3 p sj Γ d sj p sj , E Y n ′ 2 +r ′ 3 jd = a n ′ 2 +r ′ 3 jd Γ d jd +n ′ 2 +r ′ 3 p jd Γ d jd p jd .
Then, the outage probability for AF relay selection is given by

P Ψ o (γ th ) = J j=1 F Ψ Γ j (γ th ) (24)
where Ψ in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF] can be E using F E Γ j (x) to get the exact outage probability, LB using F LB Γ j (x) to get lower bound and ∞ using F ∞ Γ j (γ th ) to get the asymptotic expression for the outage probability.

Using the simple form in [START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF], several insights can be obtained. For example, in the low outage regime when γ th → 0, the above result becomes exact. Also, one can see from ( 23) that when the Nakagami-m parameters of the interference m ij or m vj (≫ 1) are large and increase, the outage probability remains almost unchanged. This is because m ij or m vj only have an influence on the order of interference power, as can be seen from ( 44

m ij +1
m ij in (45) and

(m ij +1)(m ij +2) m 2 ij in (46) cannot be ignored if m ij ≈ 1
(They approach 1 when m ij → 0). This phenomenon will be shown in Fig. 1 or m jd = 1), the outage probability remains almost unchanged if one changes INR but keeps SINR fixed. These explanations will be verified in Fig. 2 and Fig. 3 in Section V.

Also, since the possible boundary of the interferers and the pass loss between the interferer and the signal have influence on the order of interference power (referring to (44), ( 45) and

(46) in Appendix A), changing the possible boundary and the pass loss also changes the outage probability, even if one keeps the SINR and INR fixed. Similar to before, for the Rayleigh case, the outage probability changes little if one keeps SINR and INR fixed. These discussions will be verified in Fig. 5 and Fig. 6 in Section V.

Another insight that can be obtained from ( 23) is that, when m sj or m jd increase, the outage probability decreases, as m sj or m jd are in the upper limits of the summations in [START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF], which will also be examined in Fig. 2 -Fig. 6 in Section V.

D. Special cases

For the sake of simplicity, we now focus on the case when the interferences are dominant at both relay and destination. By setting σ 2 sj ≈ 0 and σ 2 jd ≈ 0 in [START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF], one can get the asymptotic CDF in this case as

F ∞ Γ j (x) = 1 - N 2 n 2 =0 m sj -1 r 2 =0 N ′ 2 n ′ 2 =0 m jd -1 r ′ 2 =0 (-1) n 2 +n ′ 2 m sj Ω sj n 2 +r 2 m jd Ω jd n ′ 2 +r ′ 2 E Y n 2 +r 2 sj E Y n ′ 2 +r ′ 2 jd n 2 !r 2 !n ′ 2 !r ′ 2 ! x n 2 +r 2 +n ′ 2 +r ′ 2 , (25)
where E Y 1 and m jd = 1, (23) is specialized to

F ∞ Γ j (x) = 1 - N 2 n 2 =0 N 4 n 4 =0 N ′ 2 n ′ 2 =0 N ′ 4 n ′ 4 =0 (-1) n 2 +n 4 +n ′ 2 +n ′ 4 Ω -n 2 -n 4 sj Ω -n ′ 2 -n ′ 4 jd σ 2n 4 sj σ 2n ′ 4 jd E Y n 2 sj E Y n ′ 2 jd n 2 !n 4 !n ′ 2 !n ′ 4 ! x n 2 +n 4 +n ′ 2 +n ′ 4 , (26)
where

E Y n 2 sj = a n 2 sj Γ d sj +n 2 p sj Γ d sj p sj , E Y n ′ 2 jd = a n ′ 2 jd Γ d jd +n ′ 2 p jd Γ d jd p jd .
Then, using [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF], the outage probability for the two special cases above can be obtained.

IV. FIXED INTERFERERS

In this section, we consider the case when the interferers have fixed number and locations.

In this case, the numbers of interferers I sj and I jd as well as the distances l ij and l vj are deterministic such that they can be all treated as constants. The only randomness comes from the Nakagami-m fading. Thus, Γ sj is a function of only the random channel gains. Note that similar derivations have also been conducted for dual-hop AF relaying without relay selection in the literature [START_REF] Suraweera | Performance analysis of fixed gain relay systems with a single interferer in Nakagami-m fading channels[END_REF], [START_REF] Zhong | Dual-hop systems with noisy relay and interference-limited destination[END_REF]- [START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF]. However, they either consider interferences at only one of the relay and the destination [START_REF] Zhong | Dual-hop systems with noisy relay and interference-limited destination[END_REF]- [START_REF] Benevides Da Costa | Interference-limited relaying transmissions in dual-hop cooperative networks over nakagami-m fading[END_REF], [START_REF] Lin | Outage performance of relaying with directional antennas in the presence of co-channel interferences at relays[END_REF], [START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF], only a single interferer [START_REF] Suraweera | Performance analysis of fixed gain relay systems with a single interferer in Nakagami-m fading channels[END_REF], for fixed-gain relaying [START_REF] Al-Qahtani | Exact outage analysis of dual-hop fixed-gain af relaying with cci under dissimilar nakagami-m fading[END_REF], for performance upper bounds [START_REF] Ikki | Performance evaluation and optimization of dual-hop communication over nakagami-m fading channels in the presence of co-channel interferences[END_REF], for interference-limited case with identical Nakagami-m channels [START_REF] Da Costa | Outage performance of two hop af relaying systems with co-channel interferers over Nakagami-m fading[END_REF], or for Rayleigh channels [START_REF] Lee | Outage probability for dual-hop relaying systems with multiple interferers over rayleigh fading channels[END_REF]. In the following, we will derive the exact performance for the case when both the relay and the destination suffer from multiple non-identically distributed Nakagami-m interferers and we also consider relay selection in our derivation.

A. PDF and CDF of Γ sj and Γ jd

Similarly, we derive the PDF and CDF of Γ sj first. Since

Y sj = I sj i=1 Ω ij |h ij | 2 and |h ij | 2 are
independent Gamma random variables, by proper scaling, Y sj is actually a sum of independent Gamma random variables. A closed-form expression for the PDF of this sum was derived in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF]. However, this expression uses an infinite series in order to consider the general case of June 12, 2018 DRAFT arbitrary shape parameters and scale parameters. To avoid this infinite series in our result, we use the PDF in [START_REF] Mathai | Storage capacity of a dam with gamma type inputs[END_REF] for the case when the Nakagami-m parameters could be different but are integers. In this case, the PDF of Y sj is given by [25, eq. ( 4)]

f Y sj (x) =   I sj i * =1 - Ω i * j m i * j -m i * j   I sj i=1 m ij r=1 (-1) r b ir (r -1)! x r-1 e - m ij Ω ij x , x > 0 ( 27
)
where b ir is a constant given by [25, eq. ( 5)], Ω i * j , m i * j are the same as Ω ij , m ij , respectively for the same i and j and all other symbols are defined as before. It is derived in Appendix E that the PDF of Γ sj can be written as

f Γ sj (u) = I sj i=1 m ij r=1 m sj f =0 ϕ 1,sj,ij,i,r,f u m sj -1 e - m sj Ω sj σ 2 sj u ( m sj Ω sj u + m ij Ω ij ) f +r , (28)
where ϕ 1,sj,ij,i,r,f =

I sj i * =1 - Ω i * j m i * j -m i * j (-1) r b ir m sj Ω sj m sj ( m sj f )(σ 2 sj ) m sj -f Γ(f +r) Γ(m sj )(r-1)!
, and the CDF of Γ sj can be derived as

F Γ sj (u) = 1 - I sj i=1 m ij r=1 m sj -1 f =0 f h=0 ϕ 2,sj,ij,i,r,f,h u f e - m sj Ω sj σ 2 sj u ( m sj Ω sj u + m ij Ω ij ) h+r , (29)
where ϕ 2,sj,ij,i,r,f,h =

I sj i * =1 - Ω i * j m i * j -m i * j (-1) r b ir m sj Ω sj f ( f h )(σ 2 sj ) f -h Γ(h+r) (r-1)!f ! .
Using [START_REF] Benevides Da Costa | Interference-limited relaying transmissions in dual-hop cooperative networks over nakagami-m fading[END_REF] and the following Taylor's series expansion

(1 + x) -n = N i=0 -n i x i + o(x N ), as x → 0, (30)
where n and N are positive integers, one can get the high SINR approximations for PDF and CDF of Γ sj as (31) and (32), respectively.

f Γ sj (u) = I sj i=1 N 5 n 5 =0 N 6 n 6 =0 m ij r=1 m sj f =0 ϕ 3,sj,ij,i,r,f,n 5 ,n 6 u m sj +n 5 +n 6 -1 + o (σ 2 sj u/Ω sj ) N 5 + o (uΩ ij /Ω sj ) N 6 (31) DRAFT June 12, 2018
where ϕ 3,sj,ij,i,r,f,n 5 ,n 6 =

I sj i * =1 - Ω i * j m i * j -m i * j (-1) r+n 5 b ir m sj Ω sj m sj +n 5 +n 6 m ij Ω ij -f -r-n 6 ( m sj f )( -f -r n 6 )(σ 2 sj ) m sj -f +n 5 Γ(f +r) Γ(m sj )(r-1)!n 5 ! , F Γ sj (u) = 1 - I sj i=1 N 7 n 7 =0 N 8 n 8 =0 m ij r=1 m sj -1 f =0 f h=0 ϕ 4,sj,ij,i,r,f,h,n 7 ,n 8 u f +n 7 +n 8 + o (σ 2 sj u/Ω sj ) N 7 + o (uΩ ij /Ω sj ) N 8 , (32)
where ϕ 4,sj,ij,i,r,f,h,n 7 ,n 8 =

I sj i * =1 - Ω i * j m i * j -m i * j (-1) r+n 7 b ir m sj Ω sj f +n 7 +n 8 m ij Ω ij -h-r-n 8 ( f h )( -h-r n 8 )(σ 2 sj ) f -h+n 7 Γ(h+r) (r-1)!f !n 7 ! . In high SINR conditions, o (σ 2 sj u/Ω sj ) N 5 , o (uΩ ij /Ω sj ) N 6 , o (σ 2 sj u/Ω sj ) N 7 and o (uΩ ij /Ω sj ) N 8
in the above equations can be removed to obtain corresponding approximations.

The PDF and CDF expressions f Γ jd (u) and F Γ jd (u) and their high SINR approximations can be also obtained using the same methods as above.

B. Outage probability

Using the derived exact PDF and CDFs of Γ sj and Γ jd into (20), one can get the CDF of the instantaneous end-to-end SINR in one-dimensional integral, which can be calculated numerically using mathematical software.

Then, following the same process as in Section III, one can get the lower bound of the CDF as (33), when the exact CDFs of Γ sj and Γ jd are used,

F Γ j (x) 1 - I sj i=1 m ij r=1 m sj -1 f =0 f h=0 I jd v=1 m vj r ′ =1 m jd -1 f ′ =0 f ′ h ′ =0 u f e - m sj Ω sj σ 2 sj u ( m sj Ω sj u + m ij Ω ij ) h+r u f ′ e - m jd Ω jd σ 2 jd u ( m jd Ω jd u + m vj Ω vj ) h ′ +r ′ ϕ 2,sj,ij,i,r,f,h ϕ 2,jd,vd,v,r ′ ,f ′ ,h ′ = F LB Γ j (x), (33)
where ϕ 2,jd,vd,v,r ′ ,f ′ ,h ′ =

I jd v * =1 - Ω v * j m v * j -m v * j (-1) r ′ b vr ′ m jd Ω jd f ′ f ′ h ′ (σ 2 jd ) f ′ -h ′ Γ(h ′ +r ′ ) (r ′ -1)!f ′ ! .
Also, using high SINR approximations of F Γ sj and F Γ jd , one can get the asymptotic expression of the CDF as

F ∞ Γ j (x) = 1 - I sj i=1 N 7 n 7 =0 N 8 n 8 =0 m ij r=1 m sj -1 f =0 f h=0 I jd v=1 N ′ 7 n ′ 7 =0 N ′ 8 n ′ 8 =0 m vj r ′ =1 m jd -1 f ′ =0 f ′ h ′ =0 u f +n 7 +n 8 +f ′ +n ′ 7 +n ′ 8 ϕ 4,sj,ij,i,r,f,h,n 7 ,n 8 ϕ 4,jd,vd,v,r ′ ,f ′ ,h ′ ,n ′ 7 ,n ′ 8 , (34)
where

ϕ 4,jd,vd,v,r ′ ,f ′ ,h ′ ,n ′ 7 ,n ′ 8 = I jd v * =1 - Ω v * j m v * j -m v * j (-1) r ′ +n ′ 7 b vr ′ m jd Ω jd f ′ +n ′ 7 +n ′ 8 m vj Ω vj -h ′ -r ′ -n ′ 8 f ′ h ′ -h ′ -r ′ n ′ 8 (σ 2 jd ) f ′ -h ′ +n ′ 7 Γ(h ′ +r ′ ) (r ′ -1)!f ′ !n ′ 7 !
.

Then, using [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF], the outage probability can be obtained.

C. When the interferences are dominant

In the case when the interferences are dominant at both relay and destination such that σ 2 sj ≈ 0 and σ 2 jd ≈ 0, exact CDF expression (20) can be solved in closed-form as

F Γ j (x) = F Γ sj (x) + I sj i=1 m ij r 1 =1 ϕ 7 ϕ 9 (x) - I sj i=1 I jd v=1 m ij r 1 =1 m vj r 2 =1 m jd -1 f =0 m sj -1 j 1 =0 f j 2 =0 ϕ 8 ϕ 10 (x) (35)
where

ϕ 7 = m sj Ω sj -r 1 I sj i * =1 - Ω i * j m i * j -m i * j (-1) r 1 b ir 1 Γ(m sj +r 1 )B(1,r 1 ) Γ(m sj )(r 1 -1)! , ϕ 8 = m sj Ω sj m sj I sj i * =1 - Ω i * j m i * j -m i * j I jd v * =1 - Ω v * j m v * j -m v * j f j 2 m sj -1 j 1 (-1) r 1 +r 2 b ir 1 bvr 2 Γ(m sj +r 1 )Γ(f +r 2 )Ω r 2 jd (Ω ij Ω sj ) m sj +r 1 m j 1 +j 2 +1 jd Ω j 1 +j 2 +r 2 +1 vj (r 2 -1)!f !Γ(m sj )(r 1 -1)! B(j 1 + j 2 + r 2 + 1, m sj + r 1 - j 1 -j 2 + f -1), ϕ 9 (x) = 2 F 1 1, 1 -m sj ; r 1 + 1; - m ij Ω sj m sj Ω ij x m ij Ω sj m sj Ω ij + x -m sj -r 1 +1
x m sj -1 ,

ϕ 10 (x) = (x+1) j 1 +1 x j 2 +m sj (m ij Ω sj +m sj Ω ij x) -m sj -r 1 (m jd Ω vj x+m vj Ω jd ) -j 1 -j 2 -r 2 -1 2 F 1 (m sj + r 1 , j 1 + j 2 + r 2 + 1; m sj + r 1 + f + r 2 ; 1 - m jd Ω vj m sj Ω ij (x 2 +x) (m ij Ω sj +m sj Ω ij x)(m vj Ω jd +m jd Ω vj x) ) and 2 F 1 (•, •; •; •) is the hypergeometric function. Proof : See Appendix F.
Note that (35) is a very good closed-form approximation to the exact CDF when INR is large in the case of fixed interferences.

D. When the interferences are i.i.d.

In this subsection, we focus on the case when all the interferences are i.i.d. One can see that the PDF f Y sj (x) in ( 27) is only suitable when the Nakagami-m parameters of interferences are DRAFT June 12, 2018 different. When all the interferences are i.i.d., [START_REF] Prudnikov | More Special Functions[END_REF] simply becomes

f Y sj (x) = m ij Ω ij I sj m ij x I sj m ij -1 Γ(I sj m ij ) e - m ij Ω ij x , x > 0 (36)
where b ir does not exist any more and all other symbols are defined as before. Using ( 13), ( 14), ( 21) and (36), the lower bound can be derived as

F Γ j (x) 1 -ϕ 12 m sj -1 v 1 =0 m jd -1 v 2 =0 v 1 s 1 =0 v 2 s 2 =0 ϕ 13 ϕ 14 (x) = F LB Γ j (x) (37)
where

ϕ 12 = Ω vj m vj -I jd m vj Ω ij m ij -I sj m ij e - m jd σ 2 jd x Ω jd - m sj σ 2 sj x Ω sj Γ(I jd m vj )Γ(I sj m ij ) , ϕ 13 = Ω jd m jd -v 2 Ω sj m sj -v 1 (σ 2 jd ) v 2 -s 2 (σ 2 sj ) v 1 -s 1 Γ(I jd m vj +s 2 )Γ(I sj m ij +s 1 ) s 1 !s 2 !(v 1 -s 1 )!(v 2 -s 2 )! , ϕ 14 (x) = x v 1 +v 2 m jd x Ω jd + m vj Ω vj -I jd m vj -s 2 m ij Ω ij + m sj x Ω sj -I sj m ij -s 1 .
Similar, the asymptotic outage probability in the low outage regime can be derived as

F ∞ Γ j (x) = 1 - m sj -1 v 1 =0 m jd -1 v 2 =0 v 1 s 1 =0 v 2 s 2 =0 N 9 n 9 =0 N 10 n 10 =0 N 11 n 11 =0 N 12 n 12 =0
ϕ 15 ϕ 16 x n 10 +n 11 +n 12 +n 9 +v 1 +v 2 , (

where ϕ 15 = (-1) n 10 +n 9 Γ(I jd m vj +s 2 )Γ(I sj m ij +s 1 )

-I jd m vj -s 2 n 12 -I sj m ij -s 1 n 11 n 10 !n 9 !s 1 !s 2 !(v 1 -s 1 )!(v 2 -s 2 )!Γ(I jd m vj )Γ(I sj m ij )
,

ϕ 16 = Ω ij m ij n 11 +s 1 Ω vj m vj n 12 +s 2 (σ 2 jd ) n 10 -s 2 +v 2 (σ 2 sj ) n 9 -s 1 +v 1 Ω jd m jd -n 10 -n 12 -v 2 Ω sj m sj -n 11 -n 9 -v 1 .
Then, the lower bound and asymptotic expression for the outage probability can be derived by using (37) and ( 38) in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF]. Several insights can be obtained from (38) for i.i.d. fixed interferers.

For example, one can see that with the increase of the Nakagami-m parameters of the interference m ij or m vj , the outage probability remains almost the same, which will be examined in Fig. 7.

Similar to the analysis for random interferers, with the increase of SINR, the outage probability for fixed interferer decreases accordingly, which will be examined in Fig. 9 and Fig. 10 in Section V. However, different from random interferers, changing INR will not have a great influence on the outage probability for fixed interferers if the SINR is fixed. This is because SINR dominates the outage probability for fixed interferes and changing INR just change the ratio between the noise power and the interference power but has negligible influence on the overall outage probability, which will be shown in Fig. 8 and Fig. 9 in Section V. Another insight that can be June 12, 2018 DRAFT obtained from (38) is that, with the increase of m sj or m jd , the outage probability decreases, as m sj or m jd are in the upper limits of the summations in (38), which will be examined in Fig. 8 -Fig. 10 via simulation in Section V.

In the i.i.d. case when the interferences are dominant at both relay and destination, one can further get the closed-form expression of the CDF as

F Γ j (x) = 1 - m sj -1 f 1 =0 ϕ 17 x f 1 m ij Ω ij + m sj x Ω sj -f 1 -I sj m ij + ϕ 18 ϕ 20 (x) - m jd -1 f 2 =0 f 2 s 1 =0 m sj -1 s 2 =0 ϕ 19 ϕ 21 (x) (39)
where

ϕ 17 = Γ(f 1 +I sj m ij ) f 1 ! Ω sj m sj f 1 Ω ij m ij I sj m ij Γ(I sj m ij ) , ϕ 18 = m -I sj m ij sj Ω ij m ij -I sj m ij Ω I sj m ij sj B(1,I sj m ij ) B(I sj m ij ,m sj ) , ϕ 19 = Γ(f 2 +I jd m vj )m -I jd m vj jd Ω I jd m vj jd Ω vj m vj -I jd m vj m -I sj m ij sj Ω ij m ij -I sj m ij Ω I sj m ij sj Γ(I sj m ij +m sj ) Γ(I jd m vj )Γ(I sj m ij )s 1 !s 2 !(f 2 -s 1 )!Γ(m sj -s 2 ) B(I jd m vj +s 1 +s 2 +1, f 2 +I sj m ij +m sj -s 1 -s 2 -1), ϕ 20 (x) = x m sj -1 Ω sj m ij m sj Ω ij + x -I sj m ij -m sj +1 2 F 1 1, 1 -m sj ; I sj m ij + 1; - m ij Ω sj m sj Ω ij x , ϕ 21 (x) = (x+1) s 2 +1 x m sj +s 1 m sj Ω ij m sj Ω ij x+m ij Ω sj I sj m ij +m sj m jd Ω vj m jd Ω vj x+m vj Ω jd I jd m vj +s 1 +s 2 +1 2 F 1 (I sj m ij + m sj , I jd m vj + s 1 + s 2 + 1; f 2 + I sj m ij + m sj + I jd m vj ; 1 - m jd m sj Ω vj Ω ij x(x+1) (Ωjdmvjmij+mjdΩvjmijx)(Ωsjmvjmij+Ωijmsjxmvj)
).

Again, (39) is a very good closed-form approximation for the exact CDF when INR is large in the case of i.i.d fixed interferences. Using [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF], the outage probability is obtained. Simulations in Fig. 11 in Section V will show that this approximation has a very good match with the exact outage probability when INR is large or interferences are dominant.

E. When the signal experiences Rayleigh fading

In another special case when the signal experiences Rayleigh fading channel, (39) is further specialized to

F Γ j (x) = 1 -I sj m ij Ω I jd m vj jd Ω I sj m ij sj B(I jd m vj + 1, I sj m ij ) Ω vj x m vj + Ω jd -I jd m vj Ω ij x m ij + Ω sj -I sj m ij ϕ 22 (x) 2 F 1 I jd m vj + 1, I sj m ij + 1; I sj m ij + I jd m vj + 1; 1 -ϕ 22 (x) (40)
DRAFT June 12, 2018 where ϕ 22 (x) =

Ω ij Ω vj x(x+1) (m ij Ω sj +Ω ij x)(m vj Ω jd +Ω vj x
) . In the high SINR condition such that Ω sj → ∞ and Ω jd → ∞, using [14, (9.1)] as

lim ϕ→0 ϕ • 2 F 1 (a + 1, b + 1; a + b + 1; 1 -ϕ) = Γ(a + b + 1) Γ(a + 1)Γ(b + 1) (41)
in (40), one can get the high SINR approximation of (40) as

F Γ j (x) = 1 - Ω ij x Ω sj m ij + 1 -I sj m ij Ω vj x Ω jd m vj + 1 -I jd m vj . (42)
Using [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF], the outage probability is obtained. From ( 42) and ( 24), one has several insights as follows: [START_REF] Laneman | Energy-efficient antenna sharing and relaying for wireless networks[END_REF] with the increase of the number of interferers at the relay I sj or at the destination I jd , the outage probability increases; (2) with the increase of the average power of signal at the relay Ω sj or at the destination Ω jd , the outage probability decreases; (3) with the increase of the average power of interferers at the relay Ω ij or at the destination Ω vj , the outage probability increases.

V. NUMERICAL RESULTS

K ij P ij = K vj P vj = 1.
In the examples where the interferers have random number and locations, we let λ = λ I A I and assume that the distances l ij and l vj follow the uniform distribution as f l (l) = 2l L 2 , 0 < l < L, where L is the maximum radius of the disc. We assume λ = λ sj = λ jd and L = L sj = L jd in this case. In the examples where the interferers have fixed number and locations, the distances l ij , l vj and the number of interferers I sj , I jd are constant. Therefore, we assume l = l ij = l vj = 2, I = I sj = I jd = 10 and β = 3. In this case, there still exists path loss if the distances between nodes are large. The path loss is determined by l and β. This influence can be checked by examining Ω ij and Ω vj , as l and β can be absorbed by Ω ij and Ω vj as part of the average powers of the interference. In the calculations for both cases above, we assume the number of relays J = 2, b sj = c sj = c jd = 1 and b jd = 10. Also, let the June 12, 2018 DRAFT values of m sj be the same for any j; m jd be the same for any j; m ij be the same for any i, j; m vj be the same for any v, j. Note that our results are general enough to include other cases but these settings are used here as examples.

Figs. 123456show the outage probability vs. γ th in the case when the interferers have random number and locations. The GGA curve is obtained by using ( 15), ( 16) and ( 20) in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF] with numerical integration, the lower bound curve is obtained by using ( 22) in ( 24) while the asymptotic curve is obtained by using ( 23) in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF]. In general, one can see that the outage probability decreases when the value of γ th decreases or when the channel condition changes from Rayleigh Fig. 6. However, the outage performances for the Rayleigh case in these two cases above keep almost unchanged for the reasons explained below [START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF].

In all these cases above, the results based on GGA match very well with the simulation results, showing the accuracy of the approximation and the usefulness of our results. Also, from Figs. 23456, one can see that the lower bounds have considerable match with the simulation while the asymptotic curves match well with the simulation for small γ th . On the other hand, one can see that the gap between lower bound and simulation decreases when γSINR increases, as expected, when comparing Fig. 4 with other figures above. 29) in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF] with numerical integration, the lower bound curve is obtained by using (37) in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF] while the asymptotic curve is obtained by using (37) and ( 17) in [START_REF] Moschopoulos | The distribution of the sum of independent gamma random variables[END_REF]. In general, one sees that the outage probability decreases when the value of γ th decreases or when the Nakagami-m parameter increases. Also, our derived exact results match very well with the simulation results and our June 12, 2018 DRAFT derived lower bounds and asymptotic curves have considerable matches with the simulation, especially for small γ th in these figures, which verify the accuracy of our analysis. The influence of m ij and m vj is examined in Fig. 7 for γSINR = 15 dB, γINR = 0 dB, m sj = 4 and m jd = 5.

One can see that the curves with m ij = 1 and m vj = 1 has a slightly worse outage probability than the curves with m ij = 2, m vj = 3 and m ij = 6, m vj = 7 while the curves with m ij = 2, m vj = 3 are almost the same as the curves with m ij = 6, m vj = 7. Fig. 8 shows the result for γSINR = 15 dB, γINR = 0 dB while Fig. 9 shows the result for the same conditions except that γINR is increased from 0 dB in Fig. 8 to 20 dB in Fig. 9. One can see that the outage probability from simulation remains almost unchanged, as the SINR dominates the outage probability in the case of fixed interferers and the influence of changing INR can be ignored in this case. Fig. 10 shows the same conditions as Fig. 9 except that γSINR is increased from 15 dB in Fig. 9

E(Y c sj ) = a c sj Γ d sj +c p sj Γ d sj p sj (43)
where c represents the c-th order moment. Denote

η 1,ij = ∞ 0 f l (l ij )η(l ij )dl ij , η 2,ij = ∞ 0 f l (l ij )η 2 (l ij )dl ij and η 3,ij = ∞ 0 f l (l ij )η 3 (l ij )dl ij . One can derive the first-order moment of Y sj as E{Y sj } = ∞ I=0 e -λ I A I (λ I A I ) I I! I i=1 E{Ω ij }E{|h ij | 2 } = ∞ I=0 e -λ I A I (λ I A I ) I I! I i=1 K ij P ij η 1,ij , (44)
the second-order moment of Y sj as

E{Y 2 sj } = ∞ I=0 e -λ I A I (λ I A I ) I I! E I i 1 =1 Ω i 1 j |h i 1 j | 2 I i 2 =1 Ω i 2 j |h i 2 j | 2 = ∞ I=0 e -λ I A I (λ I A I ) I I! I i=1 K 2 ij P 2 ij η 2,ij m ij + 1 m ij + I i 1 =1 I i 2 =i 1 =1 K i 1 j P i 1 j K i 2 j P i 2 j η 2 1,ij (45)
DRAFT June 12, 2018 and the third-order moment of Y sj as

E{Y 3 sj } = ∞ I=0 e -λ I A I (λ I A I ) I I! E I i 1 =1 Ω i 1 j |h i 1 j | 2 I i 2 =1 Ω i 2 j |h i 2 j | 2 I i 1 =1 Ω i 3 j |h i 3 j | 2 = ∞ I=0 e -λ I A I (λ I A I ) I I! I i=1 K 3 ij P 3 ij η 3,ij (m ij + 1)(m ij + 2) m 2 ij + I i 1 =1 I i 2 =i 1 =1 K 2 i 1 j P 2 i 1 j K i 2 j P i 2 j η 2,ij η 1,ij m ij + 1 m ij + I i 1 =1 I i 2 =i 1 =1 I i 3 =i 2 =i 1 =1 K i 1 j P i 1 j K i 2 j P i 2 j K i 3 j P i 3 j η 3 1,ij . (46)
Note that (44), ( 45) and ( 46) require an infinite series. However, in reality, one does not need to include many terms in the calculation as e -λ I A I (λ Assume independent random variables u, x > 0 in the equations below. Using [START_REF] Suraweera | Performance analysis of fixed gain relay systems with a single interferer in Nakagami-m fading channels[END_REF] and [START_REF] Papoulis | Probability, Random Variables and Stochastic Processes[END_REF] in [START_REF] Gradshteyn | Table of Integrals, Series, and Products[END_REF] and after some manipulations, one has (50)

f Γ sj (u) = p sj m sj Ω sj
This integral can be solved by using [START_REF] Prudnikov | The algorithm for calculating integrals of hypergeometric type functions and its realization in reduce system[END_REF] as [START_REF] Zhong | Dual-hop systems with noisy relay and interference-limited destination[END_REF]. (51)

By interchanging the order of integration and solving the integration over t first using [ where γ(•, •) is the lower incomplete Gamma function [START_REF] Gradshteyn | Table of Integrals, Series, and Products[END_REF]. Then, using [14, (8.352)] to expand the lower incomplete Gamma function as a finite series, one can get (53) By using [14, (3.381)] and using the same method as that for (50) twice, one can get the CDF of Γ sj as [START_REF] Suraweera | Performance analysis of two hop amplify-and-forward systems with interference at the relay[END_REF]. . Therefore, (20) can be solved as (35) with the help of [14, (3.197)].

F

  ) and Ω vj = K vj P vj |h vj | 2 η(l vj ) are the average power of the Nakagami-m fading gain in the channel between DRAFT June 12, 2018

(- 1 ) n 2 +n 4 a n 2 +r 3 sj σ 2 (r 2 -r 3 +n 4 ) sj m sj Ω sj n 2 +r 2 +n 4 Γ d sj +n 2 +r 3 p sj n 2

 12242 !n 4 !r 3 !(r 2 -r 3 )!Γ d sj p sj . Proof : See Appendix D. June 12, 2018 DRAFT

  and µ 0,jd = k -k jd jd l l jd jd a -k jd p jd jd m jd Ω jd -l jd . Also, if one inserts the CDFs of high SINR approximations into (21) and with the help of DRAFT June 12, 2018

  ), (45) and (46) in Appendix A, where the c-th order moment of interference power E (Y sj ) or E (Y jd ) remains almost unchanged when m ij or m vj are large. For small values of m ij or m vj (≈ 1), they still have some influence on the outage probability, as

  AND DISCUSSION In this section, numerical examples are presented to show the effects of the number and locations of interferers by using the outage probability expressions derived in the previous sections. We assume γSINR = b sj γSINR sj = b jd γSINR jd , γINR = c sj γINR sj = c jd γINR jd , β = β ij = β vj and

  Fig. 3 to 20 dB in Fig. 4. One can see that the outage probability decreases with the increase of γSINR , as expected. Also, comparing Fig. 3 with Fig. 5, one can see that the outage probability for m sj = 2, m jd = 3, m sj = 2, m jd = 3 and m sj = 4, m jd = 5, m sj = 6, m jd = 7 increases when the value of β changes from 3 in Fig. 3 to 5 in Fig. 5. Comparing Fig. 3 with Fig. 6, one can see that the outage probability for m sj = 2, m jd = 3, m sj = 2, m jd = 3 and m sj = 4, m jd = 5, m sj = 6, m jd = 7 increases when the value of L increases from 10 in Fig. 3 to 20 in

Fig. 1 .Figs. 7 -

 17 Fig. 1. Outage probability vs. γ th for random interferers when γSINR = 15 dB, γINR = 0 dB, λ = 50, L = 10, β = 3, m sj = 4 and m jd = 5.

Fig. 2 .

 2 Fig. 2. Outage probability vs. γ th for random interferers when γSINR = 15 dB, γINR = 0 dB, λ = 50, L = 10 and β = 3.

Fig. 3 .Fig. 4 .

 34 Fig. 3. Outage probability vs. γ th for random interferers when γSINR = 15 dB, γINR = 20 dB, λ = 50, L = 10 and β = 3.

Fig. 5 .Fig. 6 .

 56 Fig. 5. Outage probability vs. γ th for random interferers when γSINR = 15 dB, γINR = 0 dB, λ = 50, L = 10 and β = 5.

to 20 dB in Fig. 10 . 2 Fig. 7 .Fig. 8 .

 2010278 Fig. 7. Outage probability vs. γ th for fixed interferers when γSINR = 15 dB, γINR = 0 dB, m sj = 4, m jd = 5.

Fig. 9 .Fig. 10 .

 910 Fig. 9. Outage probability vs. γ th for fixed interferers when γSINR = 15 dB, γINR = 20 dB.

2 Fig. 11 .

 211 Fig. 11. Outage probability vs. γ th for fixed interferers when γSINR = 10 dB, m sj = 2, m jd = 3, m ij = 2 and m vj = 3.

m sj u m sj - 1 a 2 sj 0 x d sj +r 1 0 x d sj +r 1 - 1 G 1

 12010111 d sj sj Γ(m sj )Γ d sj p sj ∞ σ (x -σ 2 sj ) d sj -1 x m sj exp -x -σ 2 sj a sj p sj m sj xu Ω sj dx. sj )Γ d sj p sj r 1 !(m sj -r 1 )! × ∞ integral in (49) can be transformed by replacing the exponential functions with the Meijer's G-function as[27, pp. 346]f Γ sj (u) = m sj r 1 =0 p sj m sj !a -d sj sj m sj Ω sj m sj σ 2(m sj -r 1 ) sj u m sj -1 e -m sj σ 2 sj u Ω sj Γ(m sj )Γ d sj p sj r 1 !(m sj -r 1 )! × ∞

0 f Γ 1 ( 2 sj t m sj - 1

 0121 CDF OF Γ sj FOR RANDOM INTERFERERS Using the definition of CDF and (48), one hasF Γ sj (u) = u (x -σ 2 sj ) d sj -1 x m sj exp -x -

  Γ sj (u) = p sj (m sj -1)!a

1 n 1 (- 1 ) n 1 m m sj +n 1 sj sj n 1 µ 4 ,(- 1 ) n 2 a n 2 +r 3 sj σ 2 (r 2 -r 3 ) sj m sj Ω sj n 2 +r 2 Γ d sj +n 2 +r 3 p sj n 2 ϕ 6 ,

 111114122226 APPENDIX DDERIVATION OF THE HIGH SINR APPROXIMATIONS FOR PDF AND CDF OF Γ sj FOR RANDOM INTERFERERS Using Taylor's series expansion of[START_REF] Benevides Da Costa | Interference-limited relaying transmissions in dual-hop cooperative networks over nakagami-m fading[END_REF] and[14, (3.381)] into (49) and (53), one can calculate the PDF and CDF of Γ sj as (54) and (55), respectively,f Γ sj (u) = N =0 m sj r 1 =0 µ 3,sj,r 1 ,n 1 e -m sj σ 2 sj u Ω sj u m sj +n 1 -1 + o (u/Ω sj ) N 1 ,(54)whereµ 3,sj,r 1 ,n 1 = !Γ(r 1 +1)Γ(m sj -r 1 +1)Γ sj,r 2 ,r 3 ,n 2 u n 2 +r 2 e -m sj σ 2 sj u Ω sj + o (u/Ω sj ) N 2 , (55)where µ 4,sj,r 2 ,r 3 ,n 2 = !r 3 !(r 2 -r 3 )!Γ d sj p sj . Furthermore, using (17) in (54) and (55) again, one can get (18) and[START_REF] Lin | Outage performance of relaying with directional antennas in the presence of co-channel interferences at relays[END_REF].June 12, 2018 DRAFT where ϕ 5,sj,ij,i,r =I sj i * =1 -Ω i * j m i * j -m i * j m sj Ω sj m sj (-1) r b ir Γ(m sj +r) Γ(m sj )(r-1)! .Also, the CDF of Γ sj in (29) becomesF Γ sj (u) = 1sj,ij,i,r,f u f ( m sj Ω sj u + m ij Ω ij ) f +r ,(60)where ϕ 6,sj,ij,i,r,f =I sj i * =1 -Ω i * j m i * j -m i * jm sj Ω sj f (-1) r b ir Γ(f +r) (r-1)!f !

  as it is not possible to have the same interferences in the signals from different relays to destination. Otherwise, the relays have to transmit their signals at the same time in the same frequency band and the destination will not be able to tell which signal is from which relay.Based on discussions above, we assume independentNakagami-m fading channels such that the fading powers |h sj | 2 , |h jd | 2 , |h ij | 2 , and |h vj | 2 are independent Gamma random variables with shape parameters m sj , m jd , m ij , m vj and scale parameters 1/m sj , 1/m jd , 1/m ij , 1/m vj , respectively, where the Nakagami m parameters are assumed to be integers. Also, assume E{|x| 2 } = 1, E{|x ij | 2 } = 1 and E{|x vj | 2 } = 1 such that the actual average signal power is absorbed by Ω sj , Ω ij and Ω vj , respectively. Denote σ 2 sj = E{|n sj | 2 } as the noise power between June 12, 2018 DRAFT S and R j and σ 2 jd = E{|n jd | 2 } as the noise power between R j and D. Let Y sj = I sj i=1 Ω ij |h ij | 2 and Y jd = I jd v=1 Ω vj |h vj | 2 . Define the SNR between S and R j , R j and D as γ SN R

				sj	=	E{Ω sj |h sj | 2 } sj σ 2	,
	γ SN R jd	=	E{Ω jd |h jd | 2 } σ 2 jd	, respectively, the interference-to-noise ratio (INR) between S and R j , R j
	and D as γ IN R sj

, • • • , I jd , at D are independent of each other, respectively. We also assume interferers change from broadcasting phase to relaying phase such that interferences at the destination are independent of those at the relays. Note also that I sj i , |h ij | are independent of I jd v , |h vj |, respectively, for different values of j,

  Exact closed-form expressions for the PDFs of Y sj and Y jd are not available and are difficult to obtain. As a result, only moment generating functions (MGFs) of Y sj and Y jd for independent and identically distributed (i.i.d.) channels are available in the literature[START_REF] Di Renzo | Average rate of downlink heterogeneous cellular networks over generalized fading channels: A stochastic geometry approach[END_REF], which use either two-or three-dimensional integrals and thus are very complicated and not convenient to use. Therefore, it is nearly impossible to get the closed-form expressions for the PDFs and CDFs of Y sj and is the Gamma function, d sj > 0, p sj > 0 are shape parameters and a sj > 0 is the scale parameter to be determined[START_REF] Chen | Novel approximations to the statistics of products of independent random variables and their applications in wireless communications[END_REF]. It is shown in Appendix A that one can calculate p sj and

	difficult to obtain, if not impossible. Then, closed-form expressions of the PDF and cumulative d sj in (10) by solving the two equations as
	distribution function (CDF) for Γ sj and Γ jd are derived. In the third subsection, outage probability
	is derived when the relay with the largest instantaneous end-to-end SINR is selected. Note that
	the analysis of outage probability in this section can also be considered as an exact result when
	each relay and the destination have only one fixed interferer following generalized Gamma
	distribution.
	A. PDF of Y sj and Y jd
	Y jd from these MGFs because the inverse Laplace transform is further needed. Thus, in the
	following, we will use GGA by matching the first-order, second-order and third-order moments
	of Y sj and Y jd to the first-order, second-order and third-order moments of a generalized Gamma
	random variable. To the best of the authors' knowledge, none of the works in the literature have
	considered using GGA to approximate the distribution of random interferers. Numerical results
	in Section V will show that the GGA approximation has a very good match with the simulation
	results. As Y sj and Y jd have the same distribution but with different parameters, we approximate
	the distribution of Y sj first. One can get the approximate PDF of Y sj using GGA as
	f Y sj (x; a sj , d sj , p sj ) =	p sj a -d sj sj x d sj -1 e -x a sj p sj Γ d sj	p sj	, x > 0	(10)
	DRAFT where Γ(•) June 12, 2018				June 12, 2018 DRAFT

  by simulation in Section V. Also, one can see from[START_REF] Al-Qahtani | Performance analysis of dual-hop af systems with interference in Nakagamim fading channels[END_REF] that, when the SINR increases (i.e. Ω sj or Ω jd increase, or σ 2 sj or σ 2 jd decrease, or E (Y sj ) or E (Y jd ) decrease), the outage probability decreases accordingly, which will be shown in Fig.3and Fig.4in Section V. Further, for fixed SINR, if one decreases INR (i.e.

decreases the interference power E (Y sj ) or E (Y jd )), the outage probability still decreases. This is because the order of the interference power E Y n 2 +r 3 sj or E Y n ′ 2 +r ′ 3 jd in (23) also affects the outage probability. The rate of decrease becomes small when the diversity order (determined June 12, 2018 DRAFT by m sj or m jd ) is small. Therefore, when the signal experiences Rayleigh fading (i.e. m sj = 1

  n 2 +r 2In another special case when the signal experiences Rayleigh fading channel, by setting m sj =

	sj	=	a n 2 +r 2 sj	Γ Γ	d sj +n 2 +r 2 p sj d sj p sj	, E Y	n ′ 2 +r ′ 2 jd	=	a n ′ 2 +r ′ 2 jd	Γ Γ	d jd +n ′ 2 +r ′ 2 p jd p jd d jd	.

DRAFT

June 12, 2018

  Therefore approximations of (44), (45) and (46) can be made by choosing finite series. Then, one can calculate the values of a sj , p sj and d sj in (10) by solving Furthermore, with the help of Beta function B(•, •)[14, (8.384)], one can simplify (47) as[START_REF] Di Renzo | Average rate of downlink heterogeneous cellular networks over generalized fading channels: A stochastic geometry approach[END_REF], that can be solved numerically by using popular mathematical software packages, such as MATLAB, MATHEMATICA and MAPLE.APPENDIX B DERIVATION OF THE PDF OF Γ sj FOR RANDOM INTERFERERS

	        	E{Y sj } =	a sj Γ Γ	d sj +1 p sj d sj p sj
	     	E{Y 2 sj } =	a 2 sj Γ Γ	d sj +2 p sj p sj d sj	(47)
	       	E{Y 3 sj } =	a 3 sj Γ Γ	d sj +3 p sj p sj d sj	.

I A I ) I I!

decreases quickly with I.

June 12, 2018 DRAFT

APPENDIX E DERIVATION OF THE PDF AND CDF OF Γ sj FOR FIXED INTERFERERS

Similarly, the PDF of Γ sj can be calculated using [START_REF] Papoulis | Probability, Random Variables and Stochastic Processes[END_REF] and [START_REF] Prudnikov | More Special Functions[END_REF] as

This integral can be solved by using [14, (3.351)] as [START_REF] Prudnikov | The algorithm for calculating integrals of hypergeometric type functions and its realization in reduce system[END_REF].

Using the definition of CDF and (56), one has

By interchanging the order of integration and solving the integration over t first using [14, (3.351)], one further has

Then, one can derive the CDF of Γ sj by using [14, (8.352)] and [14, (3.381)] in closed-form as