Stochastic Geometry-Based Analysis of the Impact of Underlying Uncorrelated IoT Networks on LoRa Coverage

IoT networks are more and more present nowadays. Some IoT protocols share the same bandwidth leading to interference on neighboring networks and decrease of overall coverage. To contribute to this problem, an analytical study of the coverage of a LoRa network with underlying uncoordinated IoT networks for uplink transmissions is presented in this paper. Using stochastic geometry, closed form analytical expressions are proposed allowing to analyze the success and coverage probabilities for a LoRa network. An appropriate model of the path loss including real-life values is used to characterize the log-distance propagation parameters. The interference comes from both the LoRa network itself and the underlying IoT networks, modeled with an $\alpha $ -stable distribution based on recent measurements. It is shown that for an environment with a huge amount of surrounding uncoordinated IoT networks, the gateways deployment should be doubled to reach a decent coverage probability, compared to an environment where the underlying interfering networks are not considered.


I. INTRODUCTION
The importance of Internet of Things (IoT) in today's wireless communications systems is increasing day after day. Indeed, many applications related to these new types of Machineto-Machine (M2M) communications are implemented and integrated in our everyday life: temperature monitoring at given locations in a city [1], waste management [2], smart city usages [3] such as public surveillance, intelligent traffic management [4], [5], smart building monitoring [6], [7], or intelligent healthcare control [8]. With the environmental and health issues that have emerged over the past few years, IoT is quickly becoming a cornerstone of our well-being.
It is estimated by CISCO that more than 500 billion devices will be part of IoT networks by 2030 [9]. Several protocols that are classified according to their ranges have emerged in recent years. First of all, for short radii, so called WPAN for Wireless Personnal Area Networks, protocols like Zig-Bee [10] or z-Wave [11] have a range of a few dozens of meters. Other protocols classified in the WLAN (for Wireless Local Area Networks) category have longer areas from The associate editor coordinating the review of this manuscript and approving it for publication was Rupak Kharel .
one to a few kilometers, such as 802.11ah [12]. Finally, protocols like LoRa [13] or SigFox [14] have a range that reaches ten kilometers, and are classified in the Low Power Wide Area Networks category (LPWAN). All these protocols use unlicensed bands (or ISM bands), with a carrier frequency of 868 MHz in Europe and 915 MHz in North America.
Cohabitation between these different networks that share the same frequency bands creates interference on the signals received by the gateways, despite the fact that the modulations for each protocol are different. Therefore, in this paper, we propose to analyze the impact of underlying uncoordinated IoT networks on a LoRa network.

A. RELATED WORKS
Many works about Internet of Things [15], and more especially LoRa communications, emerged a few years ago. In the following of this section, we sort the related works dealing with the three main topics involved in our works: • LoRa networks modeling, • path-loss modeling and • interferential networks identification and modeling.

1) LoRa NETWORKS MODELING
Different means to model a network are possible. In this work, we propose to use stochastic geometry [16]. In [17], the authors use a Poisson rain (that is a space-time Poisson process of packets) to model a LoRa network, and analyze the successful reception probability versus the node density for all the spreading factors. The authors in [18] model a unique LoRa cell with end devices (EDs) that are uniformly distributed within a circular area. In [19], the authors model the LoRa network with Poisson Cluster Process (PCP). They calculate theoretically the intra-and inter-cluster interference, and analyze both the coverage probability and the energy efficiency for various values of LoRa nodes and cluster radii. The coverage probability is also calculated using stochastic geometry in [20]. A similar model where both the EDs and the gateways (GWs) are deployed uniformly with homogeneous Poisson-Point Processes (PPPs) is applied by the authors in [21] to optimize the coverage probability of a LoRa network with dual hops. The authors in [22] consider one LoRa cell and use the LoRaWAN capacity model given in [23] to calculate the Packet Delivery Ratio (PDR) for an inhomogeneous distribution of EDs around one gateway. In [24], the authors provide closed-form expressions of coverage probability and area spectral efficiency for a LoRa network with one gateway. Finally, in [25], the success probability and the coverage probability of a LoRa network are calculated litterally and analyzed with the help of stochastic geometry. The results are given versus the densities of GWs and EDs. This work is a good starting point of the present work, even if it takes into consideration a few approximations that are specified in the following sections.

2) PATH LOSS MODELS
Modeling the path-loss for LoRa communications is of crucial interest. Indeed, the path-loss conditions the RSSI of a signal, and therefore the spreading factor that is used. In this context, the authors in [26] evaluate the accuracy of RSSI for different LoRa devices under laboratory conditions, and during a measurement campaign in a semi-urban area. Their results lead to a first path-loss model. The works presented in [27] show the inaccuracy of the parameters of the path-loss models for a relatively dense environment (in Dortmund) for 433 and 868 MHz bands. Various measurements allow the characterization of a new path-loss model that is very close to the observed reality. The modeling is mainly based on the determination of the path-loss exponent (denoted η) and the intercept path-loss. In [28], the authors draw on the work cited above in order to propose path-loss models in the case of rural and urban settings in Lebanon. These models are very accurate for distances of up to 8 km in an urban environment, and 45 km in a rural environment. In [29], the authors propose path-loss models in urban, forest and coastal environments. Finally, the authors in [30] propose to perform similar methods in order to characterize the path-loss for indoor environments.

3) INTERFERENCE IDENTIFICATION AND MODELS
For this paper, the interference coming from neighboring uncoordinated non-LoRa networks must be characterized. In a few works, the authors analyze theoretically and experimentally the interference on a LoRa network, in the case of one gateway. In [31], the authors show the impact of a IEEE 802.15.4g interfering device on a LoRa link, and viceversa. It is shown that both links have some resilience to interference. In [32], the authors propose two optimization algorithms to determine the best LoRaWAN configurations given an underlying interfering network (here IEEE 802.15.4G). The model that is used in this paper takes into account only one gateway.
In our work, we aim to characterize the underlying IoT networks on a global network, in a multiple gateway scenario, whatever the underlying networks. Recently, this characterization has been the subject of several studies. The packets sent by IoT devices are very short compared to common cellular communications. Then, the noise induced by the underlying IoT networks (e.g. SigFox, 802.11ah, LoRa) can be seen as impulsive [33]. The power and the amplitude of such interference are considered as heavy-tailed random variables, and can be modeled by α-stable distributions [33]. The authors in [33] show that the modeling of interference in IoT networks by α-stable distributions is still valid by adapting the characteristic exponent of the distribution, even when we introduce a guard zone around the receiver and when we limit the power to a finite value. This validation is carried out by comparing the quantiles between the theoretical model and the empirical model. Finally, the works exposed in [34] show a very strong match between the theoretical interference model that is based on an α-stable heavy-tailed distribution and the measurements carried out in different zones (park, hospital, residential and industrial zones) in the city of Aalborg. In view of these related works, α-stable distribution is a good choice to model the interference from the underlying networks.

B. CONTRIBUTIONS AND ORGANIZATION
As seen in the previous sections, the related works that deal with LoRa networks modeling mainly focus on co-SF interference to calculate the coverage probability. These works do not take into account the uncoordinated neighboring networks, and do not model the path-loss with real parameters. Note that the premises of the studies presented in this article were introduced in [35]. The present paper proposes an improved and more complete system model than the one included in [35], and more thorough results and analysis.
In this paper, a LoRa network in an uplink configuration that considers various end devices densities and gateways densities is introduced. This system permits to model both dense and sparse networks. Then, the main contributions of this work are as follows: • The presence of all possible underlying IoT networks (e.g. SigFox, 802.11ah, LoRa) concurrently to the considered LoRa network is taken into account for a coverage probability calculation in a multiple gateway scenario. The interference induced by these uncoordinated networks are modeled by a heavy-tailed distribution. The parameters of the α-stable distribution are obtained by in-situ measurements and have been validated in [33], [34]. This model has been shown to fit the data.
• The path-loss used for all the communications is modeled by a regression curve that aggregates real-life pathloss exponent and intercept. Indeed, most stochastic geometry-based papers dealing with LoRa networks usually consider common path-loss model. The contribution of our approach lies in the use of such a more realistic path-loss model with experimental data proven and validated in [27] in order to get as close as possible to the experimental values of path loss for LoRa links at 868 MHz.
• The coverage probability for this LoRa network with the improved path-loss model and the α-stable distributed uncoordinated networks interference is given for various values of ED and GW densities. Note that the analytical results for one gateway have been validated by numerical simulations in [25], and the system model with multiple gateways based on stochastic geometry calculations was compared and validated with simulations in our previous works [16], [36].
• The values of the SF-based ED densities are not obtained with analytical approximations (as made in related works) but with Monte-Carlo simulations that converge to results that are much closer to the true values.
In Section II, we introduce the system model. All the analyzed metrics and the theoretical approach on the interference models are presented in Section III. In Section IV, we depict and discuss the analytical and numerical calculations for SNR, SIR and coverage probability. Finally, conclusions are given in Section V.

II. SYSTEM MODEL
In this section, we introduce the system model through several assumptions. The system includes a LoRa network with another underlying IoT network. Note that for this paper, we focus on the uplink. This system is schematized on Fig. 1. All the notations introduced in this model are given in Table 1. The system models an infinite LoRa network, made up of gateways (GWs) and end devices (EDs) in an infinite R 2 space [37]. The GWs are uniformly and randomly located in the R 2 space, and are modeled by a homogeneous Poisson-Point Process (PPP) G with intensity λ G . The EDs are modeled by an independently marked PPP [36] denoted as where {e i }, P, d ij , and {SF i } denote the sets of ED locations, the ED transmit power, the length of the LoRa radio links (i.e. the distance between the i-th ED and its receiving j-th GW), and the assigned Spreading Factor (SF). The locations e i are determined according to an unmarked PPP E ∈ R 2 with intensity λ E (with λ G λ E ). We consider that all EDs transmit with a constant power P = 19 dBm [38].

B. ASSUMPTION 2 (EUCLIDEAN DISTANCES BETWEEN END DEVICES AND GATEWAYS)
Each GW is characterized by its index j (G j ) and its position g j on the plane. We consider that G 0 is located at the origin of the coordinate system, i.e. g 0 = (0, 0). Likewise, each ED is characterized by its index i (E i ) and by its position e i . Each ED is connected to the closest gateway (as explained by the authors in [37], the distance r to the closest gateway is distributed according to f 1 (r) = 2π λ G re λ G πr 2 ). Assuming that E i is connected to G j , the Euclidean distance (in km) between E i and G j is d ij = e i − g j .

C. ASSUMPTION 3 (BASEBAND RECEIVED SIGNAL BY A GATEWAY)
As we focus on the uplink, the GW receives the useful signal, the interference coming from the other LoRa EDs, the aggregate interference coming from the underlying uncoordinated IoT networks, and the Additive White Gaussian Noise [25], [32]. Then, the baseband received signal by G j can be given by (2), as shown at the bottom of the next page, where U j [n] denotes the useful signal, I j [n] the interference created by the other LoRa transmitters from the studied LoRa network, p(d ij ) the path loss attenuation, g t,i the i-th transmitter ED antenna gain, g r,j the j-th receiver antenna gain, h ij the fast fading coefficient between the i-th transmitter and the j-th receiver, S i [n] the unit-variance signal sent by the i-th ED,  In this paper, we assume that the signals sent by the EDs experience a path loss attenuation that is a function of the distance between each ED and the receiving GW. The path loss in dB varies with distance according to the following equation [28], [29]: where η, d ij and PL 0 denote the path loss exponent, the distance between E i and G j in km and the path loss at a reference distance d 0 = 1 km, respectively. Then, for a given ED i, the path loss attenuation between this ED and the j-th GW is expressed as follows: The path loss parameters directly depend on the environment. These parameters are derived from fitting curves over measured data. In our work, we propose to take into account the path loss parameters measured and modeled in the city of Dortmund for LoRa links at 868 MHz, given in [27], i.e. η = 2.65 and PL 0 = 132.25 dB. Moreover, we assume that both the receiving and the transmitting antennas are isotropic. Consequently, both the receiving and transmitting antennas have a unit gain, i.e. g t,i = g r,j = 1, ∀(i, j) [25].

E. ASSUMPTION 5 (SMALL-SCALE FADING AND AWGN NOISE)
Small-scale fading describes the rapid fluctuations in the amplitude and phase of a signal over a short period of time or a short distance [39]. In this model, the envelope of the received signal is statistically described by a Rayleigh distribution. The small scale fading induced between E i and G j denoted as h ij follows a Rayleigh distribution. Then, in terms of power, h ij 2 ∼ Exp (1). We assume Additive White Gaussian Noise (AWGN), denoted as A[n] in (2), with zero-mean and power P A = −174 + NF + 10 log BW dBm, where NF = 6 dB is the receiver noise figure and BW = 125 kHz is the bandwidth of the signal [25].

F. ASSUMPTION 6 (SPREADING FACTOR ALLOCATION)
In the case of LoRa communications, an SF is allocated to each ED according to its distance from the GW to which it is connected [25]. This allocation is typically done by the NetServer which sends feedback in response to short test frames transmitted by the ED after connecting to the network. The SF are between 7 and 12. Tab. 2 and Fig. 2 list the allocations of SF (SF i ) as a function of the distance d ij between an ED (E i ) and its associated GW G j . For example, if the distance d ij between E i and its associated gateway G j is comprised between l 3 and l 4 (i.e. between 3 and 4 km), the SF associated with this ED will be SF i = 10. This allocation of SF according to the distances EDs-GWs therefore leads to a spatial formation of rings centered on the positions of the GWs, as shown in Fig. 3.

G. ASSUMPTION 7 (NON-LoRa INTERFERING NETWORKS)
We consider that the studied LoRa network is located in a space where other underlying IoT networks using the same European ISM frequency bands (868 MHz) such as Sigfox, z-Wave, 802.11ah or Zigbee coexist. Even though they do not use the same technology, these non-LoRa networks interfere with LoRa communications. As studied in [34], this (2) VOLUME 10, 2022  interference can be regarded as noise. However, the modeling of this interference by an AWGN is not a good choice. Therefore, we propose to model the interference Z j by an α-stable noise, the measurement data of [34] indeed suggesting that the distribution of the power of the interference is heavy tailed [40]. It has been shown in [41] and in [34] that network interference are well-modeled by impulsive noise, and more particularly by α-stable noise. The measurement in [34] details a LoRa-like configuration with frequency bands of 125 kHz. The power is shown to be heavy tailed. Considering the α-stable model from [41], conditioned on the set of active interfering devices α , the interference is Gaussian and its mean power depends on this set α . Finally, the power, unconditioning on α , is a totally skewed α-stable random variable, as it has been proven and validated analytically and experimentally in [33], [34]. Note that even if the powers of the LoRa emitters is fixed, the powers of the non-LoRa emitters are random. Then, the use of α-stable distribution to model the interference coming from the underlying uncorrelated IoT networks is valid according to [41].
Thus, in (2), the term Z j is considered as an interferential noise for the j-th GW, and the power P Z of this interferential noise follows an α-stable distribution: where α, γ , β and δ are the characteristic exponent (with 0 < α < 2), the scale parameter (with γ ∈ R + ), the skew parameter (with β ∈ [−1, 1]) and the shift parameter (with δ ∈ R), respectively. This model is not exactly supported by the data from [34] because the α should be given by the inverse of the channel attenuation (so less than one) but this is generally not the case in the measurement with estimated values ranging from 0.94 in the Hospital complex to 1.77 in the industrial area. However this is consistent with [40] that shows that a more realistic device location gives an interference adequately modeled by an α-stable distribution but with a modified α. In the following we will choose a value of alpha of 1.69 which is obtained for the residential area. Measured dispersion (γ ) values are found between 10 −11 and 10 −10 . They are expected to increase with the increase in device density as it is shown in the theoretical results from [34]. Consequently we chose a value of 10 −10 and 2.10 −10 for a moderately interfering scenario and a denser scenario.

III. STOCHASTIC GEOMETRY METRICS AND INTERFERENCE CHARACTERIZATION
A. PRELIMINARY DEFINITIONS 1) SUCCESS PROBABILITY In many works dealing with stochastic geometry applied to cellular networks such as [36], the performances are measured from the Signal-to-Interference-plus-Noise-Ratio (SINR). Nevertheless, it is not possible to calculate the coverage probability of a LoRa network with this SINR-based methodology [23].
Thus, in a LoRa network, two conditions must be met to specify the success of a transmission from an ED to a GW: i) the Signal-to-Noise Ratio (SNR) should be greater than a certain threshold q SF i depending on the used SF and ii) the Signal-to-Interference Ratio (SIR) should be greater than a certain threshold w. Then, the success of the transmission between the i-th ED and the j-th GW can be written as follows [25]: where SNR ij , SIR ij , q SF i and w denote the SNR and the SIR of the LoRa link between the i-th ED and the j-th GW, the SNR threshold for the SF used for the communication between E i and G j , and the SIR collision threshold between different SF, respectively. Obviously, the success of a transmission depends on the distance between the transmitter ED and the receiver GW. We introduce the success probability metric that is depicted as follows [25]: which means that both conditions on the SNR and the SIR for the i-th ED are met with at least one GW. Thus, the success probability corresponds to the probability that an ED can be connected to at least one GW.

2) COVERAGE PROBABILITY
As seen in the previous section, the success probability of a transmission is linked to the distance between the transmitting ED and the receiving GW. Recall that the coverage probability is the probability that a random element of the transmitters is not in outage, at any chosen time. Thus, the coverage probability can be expressed as follows [25]: where f d (x) = 2π λ G xe −λ G π x 2 is the distance distribution of a random ED to its nearest GW.

B. SIGNAL-TO-NOISE RATIO (SNR)
As we have seen in the previous section, the calculation of the coverage probability is possible only if we retrieve the SNR and SIR. Thus, it is necessary to characterize both the interference and the noise implied in each LoRa communication, i.e. the powers of the elements I j [n] and N j [n] in (2). As shown previously and on Fig. 1, the interference come from both the LoRa EDs and the non-LoRa EDs.

1) TOTAL NOISE
In our work, we propose to incorporate this interferential noise Z j with an α-stable distributed power P Z in the total noise N j with a power P N . As seen in the previous sections, the total noise also includes the AWGN A, with a power P A . We consider that AWGN and interferential noise are strictly independent. The total noise power can be expressed as: 2) SNR SUCCESS PROBABILITY As explained previously, the SNR can be expressed as follows: The success probability relative to the SNR directly depends on a SF-based threshold q SF which are given in Tab. 2.
The success probability relative to the SNR corresponds to the probability that the SNR is greater than or equal to the threshold (for a given SF). Then, The first equality comes from the characterization of the SNR in Section III-B2, the second equality comes from the fact that h 2 ij ∼ Exp (1), the third equality comes from the independance between the AWGN noise P A and the interferential noise P Z . In [34], it is proven that the β parameter of the α-stable distribution equals 1. This assumption leads us to the fourth equality that comes from the definition of the Laplace transform of a random variable that follows an α-stable distribution with β = 1 [42].

C. SIGNAL-TO-INTERFERENCE RATIO (SIR) 1) LoRa INTERFERENCE CHARACTERIZATION
As depicted on Fig. 1 and in (2), the interference taken into account for the calculation of the SIR come from the LoRa EDs that transmit their data in the same time interval as the considered link. In this section, we propose to characterize these interference I j [n], and their relative power I j .
LoRa technologies enable the reception of signals from numerous EDs thanks to the orthogonality of the transmission subbands and the quasi-orthogonality of SFs. Despite the fact that some works such as [43] and [44] show that the imperfect orthogonality between each concurrent SFs results to a further coverage loss, we propose to take into account the co-SF interference (i.e. the interference coming from the EDs that use the same SF). Indeed, the practical rejection gain between quasi-orthogonal SFs ranges from 16 to 36 dB [25], which is quite high for our purpose.
We consider E the set of the EDs in R. This set is formed thanks to a PPP with density λ E . In order to characterize the interference on a typical link, it is necessary to identify the set of the interferers: • the set of the EDs that use the same typical link SF ( E k with density λ E k ), • the set of the EDs that both use the same SF and transmit at the same time slot (¨ E k with densityλ E k ). The duty cycle policy enforced by ETSI for LoRaWAN communications is set to ξ = 1% [13]. Then, the density of the EDs that transmit at the same time slot isλ E k = ξ λ E k .
As the density of the effective interferers directly depends on the intensity of the EDs that use the same SF, we propose to characterize the density λ E k .

a: CHARACTERIZATION OF λ E k
The PPP of the EDs E can be divided into six PPPs E k , with k ∈ [7,12] that correspond to the PPP of all the EDs that use the same SF. For instance, if E i belongs to E 11 , it means that E i transmits its data with SF 11. According to Tab. 2, it also means that it is located between l 4 and l 5 km away from its related GW.
Nevertheless, the calculation of the densities λ E k is not so simple. Indeed, as shown in Fig. 4, these densities depend on the density of the GW deployment. For sparse GW deployment (λ G = 0.005 GWs per km 2 ), the number of EDs that use SF 12 is prominent (i.e. λ E 12 ), as shown on the left plot of Fig. 4. This is due to the fact that the distances between the GWs are large, and thus, most of the EDs are also very far from their related GW. In the contrary case, for a dense GW deployment (λ G = 0.05 GWs per km 2 ), the GWs are very close. As the EDs are linked to their nearest GW, this leads to the fact that most of the EDs use low SFs as shown on the right plot of Fig. 4.
In [25], the authors propose a theoretical approximation of the densities λ E k in both sparse and dense GW deployment. For sparse GW deployment (i.e. when λ G 1), the densities of the EDs for each SF are given by: and for dense GW deployment, these densities are approximated as follows: where v k are the fitting parameters for SF k defined in [25].
In order to characterize more deeply the interfering devices sets, we propose to i) identify practically the densities of the different SF-based interfering sets, ii) evaluate the relative error between the simulated densities and the approximated ones taken from [25], and iii) choose between both methods (simulated or approximated) for the theoretical calculations.

b: SF-BASED EDs DENSITIES IDENTIFICATION
We identify the densities for the different SF-based EDs sets using Monte-Carlo simulations. For this purpose, λ E is set to 5 EDs per km 2 . The values of the GW densities vary between λ G = 0.001 GWs per km 2 (i.e. sparse GW deployment) and λ G = 0.1 GWs per km 2 (i.e. dense GW deployment) and the densities λ E k are calculated in a 40.000 km 2 space. The convergence of the densities values is reached after 3.500 Monte-Carlo rounds. The results of the various SF-based densities λ E k are given in Tab. 3. Fig. 5 shows the relative error between the simulated ED densities for all the SFs depicted in Tab. 3 and the approximated values taken from (12) and (13). This figure clearly shows the limits of the approximations. Indeed, for SFs 8, 9 and 10, the error reaches a maximum error rate of 12% for all the GW deployment densities. Nevertheless, the error rate for SF7 is higher (21%) in the case of a very sparse GW deployment (λ G = 0.001), which is mainly due to the huge density of the SF12 EDs (as seen in the left plot of Fig. 4). On the other hand, for dense GW deployment, the error rate is quite high for SF11 and SF12, which is linked to the fact that the density of low SFs is very high in these cases, as shown in the right plot of Fig. 4.
In the following of this paper, we propose to use the values of λ E k given in Tab. 3 that are better than the ones given in [25] to calculate the related densitiesλ E k .

2) SIR SUCCESS PROBABILITY
Thanks to the previous results on the interference characterization, the SIR can be expressed as follows: The success probability relative to the SIR directly depends on a SF-based collisions threshold. In our works, as we only take into account the co-SF interference, we assume that the required aggregate SIR for a successful transmission is set to w = 1 dB, so w = 1.259 [43]. Proposition 1: The success probability relative to the SIR is: where L I j w Pp(d ij ) is the Laplace transform of the random variable I j at w Pp(d ij ) conditioned on the locations of the ED at e i and the GW at g j , and where GHF (x) is the Gauss Hypergeometric function defined as GHF (x)

D. SUCCESS AND COVERAGE PROBABILITIES
Proposition 2: For a given E i that is effectively connected to G 0 , the success probability of the transmission is given in (18), as shown at the bottom of the next page.
Proof: See Appendix B. Finally, the coverage probability is calculated with (8). Note that the coverage probability has no closed form, and is calculated numerically.

IV. NUMERICAL ANALYSIS
In this section, the stochastic geometry metrics are presented through analytical results obtained from the equations given in Section III. For more clarity, we assume that E i is connected to G 0 that is located at the origin of the plane.

A. SNR SUCCESS PROBABILITY ANALYSIS
We can see in (11) that the SNR success probability (denoted as P SNR ij ≥ q SF i ) depends on the distance between the EDs and their GW (i.e. the SF that they use for the transmission), and the interferential noise power. Fig. 6 shows the Signal-to-Noise Ratio success probability for AWGN-only noise (blue curve) and for AWGN plus interferential α-stable noise (red and green curves) vs. distance d i0 between the i-th ED and its relative 0-th GW. The α-stable noise is modeled with the values given in Tab. 1. The red curve describes the behavior of the SNR success probability in a sub-urban

FIGURE 6.
Signal-to-Noise Ratio success probability for AWGN-only noise (blue curve) and for AWGN plus interferential α-stable noise vs. distance d i 0 in km. environment (i.e. for γ = γ 1 = 1 · 10 −10 ) while the green curve stands for an urban environment (i.e. for γ = γ 2 = 2 · 10 −10 ). First, in all cases, we can see that the curves decrease with the distance d i0 . The curves are discontinuous at all l k step values. These two observations are relatively obvious in view of (11) which takes into account the allocations of SF as a function of the distance from E i to G 0 . It is obvious that the addition of interference decreases the SNR success probability. Indeed, this probability reaches a level of 0.8 for d i0 = 1.7 km (with SF 8) and d i0 = 2.2 km (with SF 9) for an environment without α-stable noise, for d i0 = 1.6 km with αstable noise with γ = γ 1 and d i0 = 1.4 km for γ = γ 2 . Even more explicitly, the probability is never lower than 0.2 for an environment without α-stable noise, but reaches 0.2 from d i0 = 7.85 km with an α-stable noise with γ = γ 1 and from d i0 = 6.4 km for γ = γ 2 (i.e. with SF 12). Despite the fact that the α-stable signal strength with γ 2 is only twice as strong as that with γ 1 , the SNR success probability is drastically reduced. It can therefore be concluded that the presence of many underlying IoT networks can lead to a very strong decrease in SNR, and therefore a decrease in the SNR success probability, even for very short ED-GW distances. In an very dense urban environment (in terms of IoT networks), this drastic decrease in the SNR success probability will inevitably decrease the overall coverage probability for the studied LoRa network. Fig. 7 shows the Signal-to-Interference Ratio success probability vs. distance d = d i0 between the i-th ED and its relative 0-th GW. The curves are obtained thanks to (16), for λ E = 5 EDs per km 2 and λ G = 0.005, 0.01, 0.05 (orange curve, green curve and purple curves, respectively) GWs per km 2 .

B. SIR ANALYSIS
Like the SNR, the SIR is decreasing as a function of d i0 for all the values of λ G , in a discontinuous way. This discontinuity is due to the different SIR thresholds for all the SFs. We clearly see that for a sparse network (i.e. λ G = 0.005), the SIR coverage probability is almost continuously decreasing xdx . (18) 8798 VOLUME 10, 2022 until a value of 0.4 for d i0 = 5 km, and then reaches a value close to 0. Both this continuity and the value of 0 are mainly due to the densities of EDs for each SF. Indeed, in a sparse network, the densities λ E k are very small for k ∈ [7 ; 11] (which leads to a small amount of co-SF interference for SF 7 to SF 11) and very big (i.e. close to λ E ) for SF 12 (which leads to a huge amount of interference). On the contrary, for a dense network, the density of EDs that are using SF 12 (i.e. λ E 12 ) is much lower than in a sparse network. Thus, the number of co-SF interferers with SF 12 (i.e. for a distance between the ED and the GW greater than 5 km) is smaller, and then the SIR coverage probability is better than for a sparse network. Nevertheless, as the GWs are closer one to each other compared with a sparse network, the number of EDs that share the same SF (from SF 7 to SF 11) is higher, which leads to a worse SIR coverage probability for d i0 ∈ [0, 5] km.

C. SUCCESS PROBABILITY
Figs. 8a and 8b show the success probability H (e i ) vs. the distance between the ED transmitter and its GW receiver with AWGN+α-stable noise with γ = γ 1 and γ = γ 2 , respectively. The curves are taken from (18) for λ E = 5 EDs per km 2 and λ G = 0.005, 0.01, 0.05 GWs per km 2 (in blue, orange, green, yellow,purple and red, respectively). The solid and the dotted lines represent the success probability with AWGN noise only, and the success probability for AWGN noise + interferential α-stable noise, respectively. As shown in Figs. 8a and 8b, the success probability is null for d i0 > 5 km in sparse and moderately dense networks. We have seen in the previous sections that for these GW deployments, the number of EDs using SF 12 is smaller than for other SFs. As the number of SF 12 EDs is very small, the integral taken into account in (18) is almost null.
For all the noises, the success probability follows curves with steps that correspond to the different SFs distances. This is due to both the SNR and SIR success probabilities explained in (11) and (16). Nevertheless, for an AWGN noise, the probability in a dense network is almost constant for low values of SFs, which is due to the fact that the number of EDs with low SFs is quite small compared to the ones in high SFs (and thus the SIR is high, as shown in Fig. 7).
With AWGN noise, the difference in success probability is becoming very high between dense and (moderately) sparse networks from SF 8. Indeed, the value of 0.8 is reached at d i0 = 1.6 km (SF 8) in (moderately) sparse networks, and at d i0 = 6.2 km for dense networks (SF 12). Such a difference denotes inevitably the need to increase the GW deployment density to reach a better success probability. Moreover, the maximum difference between dense and (moderately) sparse networks is reached at d i0 = 5.1 km, for a value of 0.93.
In Fig. 8a, we can see the effect of underlying IoT networks. Indeed, the addition of α-stable noise decreases the success probability in all the cases. This decrease is directly linked with the SNR success probability depicted in Fig. 6. For sparse and moderatly dense networks, the difference between LoRa-only network (i.e. with only AWGN) and LoRa+underlying IoT network is very thin for SFs 7, 8 and 12. However, for SFs 9 to 11, this difference is getting bigger, with a maximum difference of 0.06 for d i0 = 3.95 km in a sparse network, and 0.04 for d i0 = 3.98 km in a moderately dense network. Thus, we can say that the impact of underlying IoT networks on a sparse or moderately dense LoRa network is relatively weak in the case of a sub-urban environment (i.e. with γ = γ 1 ). On the contrary, in the case of an urban environment, i.e. with γ = γ 2 , the impact of underlying networks is way bigger, as seen in Fig. 8b. Indeed, the success probability is reaching a value of 0.19 for d i0 = 1.97 km (SF 8) (whereas it reaches only 0.06 for the same d i0 ). The maximum difference between AWGN and AWGN+α-stable noise is 0.24 for d i0 = 2.9 km (SF 9) in a sparse network, and 0.2 for the same value of d i0 in a moderately dense network. We can conclude that the addition of underlying IoT networks in an urban environment has a big impact on the success probability of a LoRa communication for a sparse or moderately dense GW deployment.
For a dense GW deployment, the main differences between AWGN and AWGN+α-stable noise are quite similar. In a sub-urban environment, the difference for SFs from 7 to 11 are quite thin. This difference reaches a value of 0.09 for d i0 = 4.96 km, which is acceptable. For SF 12 (recall that in a dense network, the number of EDs using SF 12 is very small), this difference reaches a value of 0.2 for d i0 = 7.8 km. Then, the addition of underlying IoT networks has a relatively big impact on the LoRa network for high SFs in this configuration. Nevertheless, as the number of EDs using high SFs is very negligible, we can say that alpha-stable noise impact is moderate. On the contrary, as shown in Fig. 8b, the impact of underlying IoT networks in an urban environment is big. Indeed, even for low SFs, the success probability decreases drastically (for SF8, at d i0 = 1.9 km, H(e i ) = 0.98 with AWGN noise, and H(e i ) = 0.86 with AWGN+αstable noise). The difference beween the two values is even bigger for higher SFs, reaching a difference of 0.32 for d i0 = 4.95 km, i.e. for SF 11. For SF 12, the maximum difference is reached for d i0 = 7 km with a value of 0.5. However, as explained previously, the number of EDs using SF 12 is very small, and thus, even if this difference is high, the impact on the overall LoRa link success probability is quite moderate. We can conclude that in an urban environment, i.e. with γ = γ 2 , the impact of underlying IoT networks on a LoRa link is severe in the case of a dense GW deployment. This impact is high for both high SFs (that are however seldom in use), and for low SFs.
Thus, we can conclude that underlying uncoordinated IoT networks have a big impact on the success probability of a LoRa communication for a sparse or a dense GW deployment and will undoubtedly have an effect on the overall LoRa coverage probability. Fig. 9 shows the coverage probability C vs. density of GWs over density of EDs. The curves are taken from (8) for λ E = 5 VOLUME 10, 2022  . Coverage probability C vs. density of GWs (in GWs per km 2 ) for λ E = 5 EDs per km 2 . The first curve (AWGN) considers both LoRaWAN internal interference and AWGN noise. The other two curves (AWGN + α-stable noise γ 1 and γ 2 ) consider internal LoRaWAN network interference, AWGN noise and external interference (modeled by α-stable distributions) in the case of a moderate interference scenario (γ 1 ), and a dense interference scenario (γ 2 ). Note that external interference modeled by alpha-stable distributions does not affect SF allocation.
We can see that for all the cases (AWGN and AWGN+α-stable noise), the coverage probability increases along the density of the gateways λ G . The difference between the AWGN and AWGN+α-stable with γ 1 is quite thin for all the gateways deployment densities, with a maximum difference of 0.04. Note that the α-stable noise with γ 1 is quite small. Thus, we consider a small number of underlying IoT networks. Nevertheless, with a bigger α-stable noise (i.e. with γ 2 ), the difference is more significant. Indeed, a difference of 0.2 in the coverage probability is reached between the AWGN and the AWGN+α-stable noise. Moreover, the cov-erage probability reaches the value of 1 when λ G = 0.1. Note that for this density, the number of EDs that use high SFs is very small. Thus, the integral in (8) takes into account almost only the very small values of d i0 . Indeed, when λ G = 0.1, the density of GW is very high, as seen in Fig. 4. The calculation of the probability of coverage for a GW given by this equation is an integral which theoretically spreads over R+. However, the function f d (x) takes λ G into consideration. This function involves a real integral that occupies only the GW's coverage radius. Thus, the d i0 integration variable is very small. These results show clearly that in an urban environment, i.e. in an area where the number of underlying IoT networks is high, the coverage probability is largely decreased. Fig. 9 also shows that, for the considered conditions, the coverage probability reaches the value C = 0.95 when λ G = 0.048 for an environment without any underlying IoT network (i.e. with AWGN only). In an environment with a little amount of underlying IoT network (γ 1 ), this value of the coverage probability is reached for λ G = 0.071. Finally, for a dense deployment of uncoordinated IoT networks (i.e. with AWGN and α-stable noise with γ 2 ), C = 0.95 when λ G = 0.093, which is almost twice as big as for the first studied case. In this specific scenario, we can conclude that in order to reach a decent coverage probability (i.e. C ≥ 0.95) in an urban environment with a big amount of underlying IoT networks, the density of the gateways must be almost doubled compared to an environment where there is no uncoordinated network.

V. CONCLUSION AND FUTURE WORKS
In this paper, we have introduced interfering underlying IoT networks in a LoRa network. The interfering networks have been modeled by α-stable distributions, and the coverage probability for various gateways densities has been calculated thanks to stochastic geometry-based analytic calculations. The theoretical calculations of the SNR, SIR, success probability and coverage probability for a LoRa link with underlying unccordinated IoT networks is the main contribution of this work. Thanks to numerical analysis based on the theoretical results and on the implementation of a real path loss model that has been taken from in-situ measurements, we have proven that the presence of underlying networks surrounding a LoRa network could induce drastic decreases in the coverage probability of the overall studied network. The coverage probability decrease is particularly high for sparse and moderately dense LoRa networks (in terms of gateway deployment), with a maximum difference of 0.2 compared to a network with no neighboring network. All the observations in this work have led us to the conclusion that for an urban environment, i.e. in a space with a huge amount of surrounding uncoordinated IoT networks, the GW deployment should be doubled to reach a decent coverage probability. The challenge of communications of EDs with multiple GWs could increase the coverage probability. In future works, it could be interesting to study other spreading factor allocation strategies that consider the channel interference [45], that choose the least solicited spreading factors [46], or that use tools from matching theory to maximize the minimal short-term average user rates [47]. These strategies could possibly increase the coverage probability of the overall LoRa network. Another perspective could be linked with the density of end devices. Indeed, in this work, we have dealt with the coverage probability for a moderately dense end devices density. It could be interesting to analyze the coverage probability for sparse and very dense networks.

APPENDIX A PROOF OF PROPOSITION 1
The success probability relative to the SIR corresponds to the probability that the ratio between the useful signal and the aggregated interference is greater than the value w. We consider that the studied i-th ED and j-th GW communicate via the SF k. Then, in (19), as shown at the bottom of the page, (a) comes from the definition of the SIR relative to the aggregated interference I j at G j . (b) comes from the fact that VOLUME 10, 2022 h ij ∼ Rayleigh 1 2 and then h 2 ij ∼ Exp (1). (c) is derived from the characterization of the interference. Note that mk is a unit indicator, i.e. if the m-th ED uses the SF k, it equals 1, otherwise it is null. In other words, mk indicates if the mth ED is an interferer ( mk = 1) to the useful signal, or not ( mk = 0). It can also be noted that we take into account all the EDs in the whole infinite R 2 space ( m>0 , as we consider that the useful link goes from the 0-th ED to the jth GW). (d) comes from the fact that the channel fadings are assumed independent. We also take into account the fact that

APPENDIX B PROOF OF PROPOSITION 2
The general success probability has been given in (7), and corresponds to the probability that a given E i can be connected to any GW G j : Then, the success probability of the transmission for a given E i that is connected to G j is equal to the complement probability of this transmission not decoded by the other GWs: We assume that we take into account only the connections to the target GW G 0 , then the success probability can be derived as given in (20), where (a) comes from the fact that we assume an independence between the SNR and the SIR for all the GWs and EDs, which leads to the lower bound inequality. In (b), we have decomposed the product in (a) with j = 0 and j ≥ 1. The success probability given in (18) can be finally found after inserting the analytical values of the probabilities linked to the SNRs and the SIRs given in (11) and (16). During his scholarship, he specialized in acoustics and advanced applied mathematics at the Technical University of Denmark in Kongens-Lyngby (Denmark). After spending three years in the music industry, he worked during four years as a Research Coordinator in the domains of music information retrieval (MIR) and cyber security. He integrated the Institut d'Électronique et des Technologies du numéRique (IETR), University of Nantes, in 2015 as a Ph.D. student. His research domains are mainly energy and spectral efficiencies in device-to-device communications and heterogeneous networks, and stochastic geometry.