Mehdi Anteur

Nathalie Thomas
email: nathalie.thomas@enseeiht.fr

Vincent Deslandes

André-Luc Beylot

On the performance of UNB for machine-to-machine low earth orbit (LEO) satellite communications

Keywords: Aloha, doppler effect, IoT, M2M, satellite, UNB

Ultra narrow band (UNB) is a widely used technology for machine-to-machine and low-power wide-area communications. Its properties, long range with small RF power, make it naturally attractive for satellite communications but also draw new challenges compared with terrestrial systems where this technology is already deployed. Indeed, the main advantage of UNB signals, their small bandwidth, makes them more sensitive to frequency drifts that are particularly present in the case of LEO satellite systems. It also implies the use of a random access method in which the carrier frequency is a parameter unknown by the receiver. In this paper, we propose a general semianalytical model to evaluate the performance of a terrestrial or satellite system using UNB technology, taking into account the multiuser interference and the frequency drift. This model is then used to assess the performance (packet loss ratio and throughput) on the return link medium access control (MAC) of a representative LEO satellite system. With our model, we analyze the effect of frequency drift on the system performance. This paper also proposes to investigate more deeply the multiuser interference modeling in order to estimate accurately the performances of a UNB system in terms of bit error rate (BER). We propose a semianalytical approach to study the interference in presence of arbitrary power imbalance that includes the effect of frequency offset and frequency drift and applicable for any linear modulation and any pulse-shaping filter. The expression of the multiuser interference is established in the general case. We then propose a methodology to compare this exact model to the Gaussian interference approximation (mainly used through the central limit theorem) in order to assess its accuracy.

INTRODUCTION

Machine-to-machine (M2M) [START_REF] Pursley | Direct-sequence spread-spectrum communications for multipath channels[END_REF][START_REF] Anteur | Ultra narrow band technique for low power wide area communications[END_REF] covers communications between devices that do not require direct human intervention. Although they have long existed, the services made possible by this type of communications are rapidly expanding and are one of the most dynamic markets in the information technology sector. The satellite field contributes to the growth of this promising market. Various M2M systems are operational or in preparation, offering the prospects of M2M systems with global coverage. Most of them, such as Orbcomm, Iridium, Globalstar, or Argos, use a low earth orbit (LEO) satellite constellation. [START_REF][END_REF] Among the variety of M2M applications, a new type of need is clearly emerging for very cheap service (< 10$ per year) on cheap terminal devices (< 10$) with very long autonomy (years with a single AA battery) and long range (> 10 km), for the transmission of very small-sized uplink messages (< 100 bytes). [START_REF][END_REF] A lot of examples could be given; the monitoring of devices by manufacturers to diagnose problems is a well-known one. Existing communication systems are not well suited to this new type of need. Sensor-oriented protocols (as ZigBee, Z-Wave, and Wibree), used in wireless personal area networks (WPAN), are not suitable for sending traffic at long range. In contrast, wireless wide area networks (WWAN) or cellular networks (as GPRS, UMTS, and LTE) offer a rich variety of long-range high data rate services but are not designed for the transmission of small messages for low consumption terminals (autonomy for a single battery in the order of days, not years). Thus, low-power wide-area (LPWA) communication systems have been introduced.

The air interface design of a LPWA transceiver essentially relies on two main waveforms: spread-spectrum or ultra narrow band (UNB).

Spread-spectrum 1 is a well-known technique and widely studied in the literature, which consists in spreading the signal into a wider bandwidth, allowing a good tolerance to interference. The principle of UNB consists in transmitting a signal with a very small bandwidth (< 1 kHz), which is particularly well suited for small uplink traffic and for LPWA communications. We will focus in this paper on this last solution, 2 used in Sigfox's network. Indeed, the long-range/low-power property of UNB [START_REF] Lassen | Long-Tange RF Communication: Why Narrowband is the De Facto Standard[END_REF] makes it very attractive for satellite communications. However, their small bandwidths make UNB signals more sensitive to frequency drifts. And the frequency drift is particularly present in the case of LEO satellite systems where, because of high satellite velocities (several kilometers per second), the Doppler rate is in the order of hundreds of hertz per second, which is in the same order of magnitude as a UNB signal bandwidth. An example of UNB signals, received from a satellite at an altitude of 720 km, is represented on Figure 1. The maximum Doppler rate is -100 Hz/s for a carrier frequency of 400 MHz or -700 Hz/s for a carrier frequency of 3 GHz (the values are negative because the frequency shift is decreasing as the satellite passes by).

Regarding access methods, classical frequency division multiplexing access (FDMA) techniques are not adapted for UNB transmissions. Indeed, it is not possible to reach the needed level of precision in synchronization (eg, few hertz for a 100-Hz channel) with low-cost electronic components embedded in end terminals. Thus, the use of UNB signals implies a random access (RA) method in which the carrier frequency is a parameter unknown by the receiver. There are many RA techniques proposed in the literature. [START_REF] Abramson | The ALOHA system-another alternative for computer communication[END_REF][START_REF] Roberts | ALOHA packet system with and without slots and capture[END_REF][START_REF] Choudhury | Diversity ALOHA-a random access scheme for satellite communications[END_REF][START_REF] Casini | Contention resolution diversity slotted Aloha (CRDSA): an enhanced random access scheme for satellite access packet networks[END_REF][START_REF] Pateros | direct sequence spread spectrum multiple access technique[END_REF][START_REF] Del Ro Herrero | High efficiency satellite multiple access scheme for machine-to-machine communications[END_REF] We will focus in our study on extended versions of Aloha for UNB signals, which we call time/frequency Aloha (TFA) for asynchronous transmissions and slotted TFA (STFA) when UNB message transmissions are synchronized in time: See Figure 2. In both cases, the frequency is chosen randomly within a specific bandwidth.

Currently, no model can be used to assess system performance in a UNB system with frequency drift. del Ro Herrero and De Gaudenzi [START_REF] Del Ro Herrero | Generalized analytical framework for the performance assessment of slotted random access protocols[END_REF] and all the reference therein propose a general framework to study the performance of random access schemes for system using forward error correction (FEC) and taking into account power imbalance. This model is more general than any previous model, but it relies on a Gaussian interference approximation, which does not hold well for asynchronous UNB signals (see Section 4.3). An attempt has been made in Do et al [START_REF] Do | On the benefits of random FDMA schemes in ultra narrow band networks[END_REF] to model multiuser interference for UNB signals, but the model is restricted to iso-power BPSK signals without FEC. In our contribution, we propose a semianalytical model that applies for many Aloha-based random access schemes, and which supports arbitrary power distribution, frequency drift, and FEC. This model relies on the definition of a vulnerability area in which the message, or packet, of interest can interfere with other ones (see Figure 3 representing a zoom on Figure 1). This allows to compute the probability of collision and then estimate the packet loss ratio (PLR)

and the throughput for the UNB transmission. To accurately investigate the bit error rate (BER) performance, the effect of packets' carrier frequency offset on the system capacity has to be taken into account but is not so trivial, because, on the one hand, it increases the number of collision, but on the other hand, the collisions are less destructive since two interfering messages are rarely completely overlapping. Even if the interference can be generated by simulating a full system, it could be tedious to implement multiple systems and compare them, hence the need of simple, easy to use, and accurate interference modeling. The main difficulty is to evaluate the probability density function (PDF) of the multiuser interference, and even though it is possible to derive it in some specific scenarios (for instance, the PDF of BPSK interferer with rectangular pulse shaping is given in previous studies [START_REF] Chiani | Analytical distribution of linearly modulated cochannel interferers[END_REF][START_REF] Yu | A review of estimating the shape parameter of generalized gaussian distribution[END_REF][START_REF] Kharrat-Kammoun | Performance analysis of ir-uwb in a multi-user environment[END_REF][START_REF] Ilow | Analytic alpha-stable noise modeling in a poisson field of interferers or scatterers[END_REF] ) in the general case, it is not easy to derive it in a tractable form. Another approach is to approximate it with a known distribution and fitted parameters.

In this paper, we propose a semianalytical approach to study the interference in presence of arbitrary power imbalance that includes the effect of frequency offset and frequency drift and applicable for any linear modulation and any pulse shaping filter. The expression of the multiuser interference is derived in the most general case. We then propose a methodology to compare this exact model to the Gaussian interference approximation (GIA), mainly adopted in virtue of the central limit theorem, in order to assess its accuracy.

The rest of the paper is organized as follows. In Section 2, we describe the system under consideration and how it is modeled. In Section 3, we describe the semianalytical model used to evaluate its random access performance. In Section 4, an exact model of the multiuser interference is derived and used to compare with Gaussian interference approximation in terms of obtained bit error rates. Finally, some conclusions are drawn in Section 5.

SYSTEM MODELING

The system under consideration is a LEO mobile satellite system operating in L band, able to provide a global coverage and an uplink traffic to millions of low-power and low-cost M2M terminals. The satellite has a single antenna beam on which it transmits a beacon to indicate its presence to the terminals. Terminals in visibility may then transmit their message to the satellite. On the ground, terminals are considered uniformly distributed and transmitting small-sized packets, low data rate, and low duty cycle traffic packets using UNB signals. The aggregate pure Aloha (PA, the users send their messages at the same frequency and at any random times), slotted Aloha (SA, the time is chosen within a discrete set of values; the frequency is the same for all users), time/frequency Aloha (TFA, time is random, and the frequency is chosen randomly within a specific bandwidth), or slotted TFA (STFA, the time is chosen within a discrete set of values; the frequency is chosen randomly within a specific bandwidth). Figure 2 summarizes these different random access schemes with and without frequency drifts.

Let us call s n (t) the complex envelop associated to a transmitted message (terminal number n). It can be written as

s n (t) = N s -1 ∑ k=0 a n,k h(t -kT s ), (1)
where a n,k represents the kth complex symbol of the nth user and is taken in the modulation M-ary alphabet. h is the pulse shaping filter impulse response,

T s = 1 R s
is the symbols period, and N s is the number of symbols in each considered burst.

Each of the s n signals goes through a different channel that consists of an amplitude factor A n (including transmit power, free path loss and fading), a delay n (including asynchronism between users and the propagation delay), a phase offset n (asynchronism between the users and the receiver and channel), a frequency shift f n relative to the carrier frequency (caused by the access method and the Doppler effect), and a frequency drift d n (assumed to be caused by the Doppler rate only). All of the parameters are independent between users, but for a specific user, the parameters might not be independent. For instance, in a system without power control when the delay is high (long distance), the amplitude is low (high free space loss). Without loss of generality, we can assume that the receiver is synchronized on the useful signal; let us say s 0 (t), which means that we can assume 0 = 0, 0 = 0, f 0 = 0, and d 0 = 0. The complex envelope of the received signal r(t) can be written as

r(t) = L ∑ n=0 A n r n (t) + n(t) (2)
with n(t) being a white complex Gaussian noise of double-sided power spectral density 2N 0 , L the number of users, and r n (t) the complex envelope of the signal received from the nth user defined as

r n (t) = N s -1 ∑ k=0 a n,k h(t -kT s -n )e j n +2 f n (t-n )+2 d n (t-n ) 2 2 . ( 3
)
Because of the asynchronism between the terminals and the satellite, n follows a uniform distribution. Since the traffic is assumed to be Poisson, the instant when the signal is received t n (including n ) is also uniformly distributed (in the case of slotted random access, t n is simply discretized). The frequency shift f n will be assumed to be uniformly distributed and independent from the other parameters since it is caused by the random access method. We also assume that the modification of the distribution caused by the Doppler shift is compensated on the transmitter side. However, the power and the frequency drift cannot be modeled with a simple distribution since these physical parameters are not independent. Therefore, to generate these two physical parameters, we consider the point of view of a single satellite whose parameters are summarized in Table 1. We first randomly pick a point under the satellite coverage, and then we compute the physical parameters at this point.

The Doppler rate and the signal to noise ratio per symbol,

E s N 0
, under the satellite coverage are shown in Figure 4. We can see for instance that high Doppler rate is highly correlated to high

E s N 0
. The dashed line represents the radio visibility, and the area in which the user are picked. Indeed, since terminals only transmit when the satellite beacon is detected, we consider that all the messages received by the satellite have a Relying on the definition of a vulnerability area in which the message, or packet, of interest can interfere with other ones, this section proposes to estimate the packet loss ratio (PLR) and the throughput for a UNB transmission. Our semianalytical model applies for many Aloha-based random access schemes and supports arbitrary power-distribution, frequency drift, and FEC.

Performance model

The two random access performance metrics that we are interested in are the throughput or capacity (in bit/s/Hz) and the packet loss ratio (PLR).

In this section, we describe the semianalytical model used in order to estimate the PLR. The throughput is linked to the PLR through the relation:

Throughput(λ) = λ(1 -PLR(λ)), ( 4
)
where λ is the average MAC load. To derive the PLR, it is assumed that the traffic is composed of packets that have the same length T and the same two-sided frequency bandwidth B. The number of transmissions within a unit time of T and frequency bandwidth of B follows a Poisson

PLR(λ) = ∞ ∑ k=0 P(k, λ)PLR k , ( 5
)
where P(k, λ) is the probability for a packet to interfere with k other packets and PLR k is the probability to lose a packet when it is interfering with k messages.

The probability P(k, λ) for a message s n (t) (packet of interest or PoI) to interfere with k other ones ({s m (t)} , m = 1, • • •, k with m ≠ n) for a given average MAC load λ can be expressed as the probability that the k other messages arrive within s n (t) vulnerability area. The vulnerability area is a geometric zone in the time/frequency space in which, if a packet begins, it will overlap with the PoI. The size of the vulnerability area depends on the packet length T, the packet bandwidth B, the frequency drift d n of the PoI, the frequency drift of the messages in collision d m , m ≠ n, and on the used access method. Note that the vulnerability area may also depend on modulation, coding, and filtering. When the frequency drift is distributed according to a continuous pdf p d (see Section 2), P(k, λ) can be expressed as

P(k, λ) = ∫ +∞ -∞ p d (d n ) (g(d n )λ) k k! exp (-λg(d n )) dd n , (6)
where

g(d n ) = + T B ∫ +∞ -∞ p d (d m )|d n -d m |dd m (7)
with = 1 for SA, = 2 for PA and STFA, and = 4 for TFA. g(d n ) representing the average number of messages interfering with a packet of frequency drift d n for a load of 1 Erlang. It is worth noting that when the frequency drift can be neglected (ie, when p d is the delta Dirac function), this formula can be simplified into the well-known expression [START_REF] Abramson | The ALOHA system-another alternative for computer communication[END_REF][START_REF] Roberts | ALOHA packet system with and without slots and capture[END_REF] :

P(k, λ) = ( λ) k k! e -λ . ( 8
)
Since P(k, λ) decreases rapidly as k increases, the PLR can be approximated with a truncated sum:

PLR(λ) ≈ K ∑ k=0 P(k, λ)PLR k , ( 9
)
where K is chosen such as

∑ K k=0 P(k, λ) > 1 -10 -3
, meaning that we consider 99.9% of the possible number of colliding packets. In the studied scenarios and for the loads under consideration, this criterion leads to K < 80. PLR k is computed numerically using a Monte-Carlo simulation.

First, the signal of a useful packet; let us say s 0 (t) is generated with a set of physical parameters

E s N 0 0 , t 0 , 0 , f 0 , d 0 .
Then the signals of k interfering packets are generated with the parameters

E s N 0 m , t m , m , f m , d m for m ∈ {1, .
., k} in a way that leads to a collision with the useful packet. A packet is counted as lost if there is a least one erroneous bit at the output of the FEC decoder.

It is assumed that the receiver performs a perfect detection of the packets and a perfect estimation of the physical parameters of the useful signal. In a practical system, these two aspects can become challenging, and the performance results obtained in the following are therefore representing an upper bound. Detailed signal detection processing considerations are beyond the scope of this paper and will be treated in future work, but the detection could for instance be done using a multidimensional search based on preamble symbols and the frequency drift estimated with a specific algorithm (for instance, see Abatzoglou [START_REF] Abatzoglou | Fast maximnurm likelihood joint estimation of frequency and frequency rate[END_REF] and references therein). These two estimators imply an increased receiver complexity, and some performance losses should be expected, especially in high interference conditions, nevertheless for the loads under consideration their impact is assumed to be negligible. To generate interfering messages, the frequency drift d m is first drawn using an arbitrary distribution. Then the instant t m is drawn, using a uniform distribution, such as a collision occurs. For two messages to interfere, they must share the same time/frequency space. The time condition for two messages of length T to interfere is straightforward:

Δt = t m -t 0 ∈ [-T, T]. ( 10
∀x ∈ R, H(x) = ⎧ ⎪ ⎨ ⎪ ⎩ 0 x < 0 1 2 x = 0 1 x > 0 , ( 11
)
it is possible to write the bounds of the frequency shift in a single simple expression:

f min (d 0 , d m , Δt) = d 0 -B + (d m -d 0 )H((d m -d 0 )Δt)Δt -(d m -d 0 )H((d m -d 0 ))T + d m -d 0 Δt. (12)
With the symmetry in frequency, it can be demonstrated that

f max (d 0 , d m , Δt) = 2f 0 -f min (-d 0 , -d m , Δt). ( 13
)
Without loss of generality, we can choose t 0 = 0 and f 0 = 0, but no specific value can be assumed for d 0 since its distribution depends on the system.

Performance analysis

The equations in the previous section enable to compute the throughput and PLR by computing two terms: PLR k and P(k, λ). In the following paragraphs, those two terms are evaluated for the system described in Section 2 and the effect of the frequency drift is analyzed.

Figure 8 shows PLR k for k taking values from 1 up to 50 for both TFA and STFA and for different symbol rates. It illustrates two phenomena.

First and as expected, we can observe that for the same number of interfering messages, TFA has a lower PLR than STFA. Indeed, in the slotted case, messages are aligned in time so interferences occur potentially during the full length of the packet while in the nonslotted case, only part of the burst is affected by the interfering packet. Secondly, we can also observe that the PLR decreases with the symbol rate. This is also expected since, as the symbol rate decreases, the signal bandwidth becomes narrower and the signal length increases. Both these aspects make the effect of the frequency drift more significant since the area where messages overlap shrinks and at the same time the number of symbols involved in the collision is reduced (see Figure 9). It is worth noting that in order to do a fair comparison, the link budget is kept the same when the symbol rate is changed (this means for instance that the terminal RF transmission power decreases with the symbol rate to keep the same E s ∕N 0 range).

However, even though the PLR difference for various symbol rates is present, it is not significant and the reduction of the PLR offered by the frequency drift can be neglected with our system parameters for R s ≥ 400 Bd (and more generally, when the signal bandwidth becomes large compared with the worst case frequency drift). Furthermore, as the number of interferers increases, the difference becomes less significant.

Combining the probability of collision with the probability of losing a packet with (5), we obtain the overall PLR (Figure 10). The figure on top (resp bottom) shows the throughput (resp the PLR) as the function of the input load. The load λ expressed in bit/s/Hz is obtained with the following formula:

λ bit∕s∕Hz = λ Erlang N b TB , ( 14
)
where N b is the payload size in bits and B is the frequency bandwidth of a packet taking into account the filter roll-off.

In the commonly admitted domain of practical use (typically PLR ≤ 10 -1 ), the PLR performances of TFA are consistent with the performance obtained for the PLR k . However for STFA, the performances are reversed: The PLR is smaller for larger symbol rates. This reflects that in STFA, the increase of collision caused by the frequency drift is not compensated by the reduction of the severity of multiple access interference. The PLR of PA and SA, given for reference, are computed assuming a perfect synchronization and a perfect correction of the frequency drift. The throughput is given for illustration purpose and for comparison with results in the literature only. For a PLR of 10 -2 , TFA can achieve a spectral efficiency of 0.013 bit/s/Hz (with R s = 100 Bd), which is +67% compared with STFA, more than two times of what PA can achieve in the same conditions and 2.4 times the spectral efficiency of SA. These figures, even though they are an improvement over those of classical Aloha schemes, are several order of magnitude away from state of the art random access scheme such as E-SSA [START_REF] Del Ro Herrero | High efficiency satellite multiple access scheme for machine-to-machine communications[END_REF] or CRDSA. [START_REF] Casini | Contention resolution diversity slotted Aloha (CRDSA): an enhanced random access scheme for satellite access packet networks[END_REF] Nevertheless, the results presented in this paper are interesting because they demonstrate that it is not necessary for a UNB transmitter using TFA to compensate for the frequency drift since the PLR is improved by it. Furthermore, these results are not limited to systems impaired by Doppler rate. They could be extended to systems in which the frequency drift is not naturally present but voluntarily introduced in order to improve performances. It would then be possible to tailor specifically the frequency drift distribution in order to maximize the gain in performance. Our study focuses on simple random access schemes, but it would be interesting to measure the benefit of such measure on state of the art random access schemes that can achieve much higher capacity. For the study of such cases the model would show its limitation since the effect of the packet detection and the frequency synchronization would no longer be negligible. The model can also be improved to take into account different packet lengths, different symbol rates or different code rates.

MULTIUSER INTERFERENCE MODELING FOR UNB TRANSMISSIONS

Introduction

This section investigates more deeply into multiuser interference modeling in order to estimate accurately the performance of a UNB system in terms of bit error rate (BER). We propose a semianalytical approach to study the interference in presence of arbitrary power imbalance that includes the effect of frequency offset and frequency drift and applicable for any linear modulation and any pulse shaping filter. The expression of the multiuser interference is established in the general case. We then propose a methodology to compare this exact model to the GIA in order to assess its accuracy.

Exact interference model

The complex envelope of the received signal r(t) can be written as in Equation 2:

r(t) = L ∑ n=0 A n r n (t) + n(t) ( 15
)
with n(t) being a white complex Gaussian noise of power spectral density 2N 0 , L the number of users, and r n (t) the complex envelope of the signal received from the nth user given in Equation 3. For the following, r n (t) will be rewritten as

r n (t) = N s -1 ∑ k=0 a n,k h(t -kT s -n )Ψ (t -n , n , f n , d n ) (16)
with the function Ψ:

Ψ (t, , f, d) = e j +2 ft+2 d t 2 2 . ( 17
)
At the output of the matched filter, the signal is given by

y(t) = r(t) ⊗ h * (-t) = A 0 y 0 (t) + L ∑ n=1 A n y n (t) + w(t), ( 18
)
where ⊗ represents the convolution product, h * is the complex conjugate of h, w(t) is the filtered noise, and y n (t) the complex envelope of the nth user signal after the matched filter. The first term is considered to be the useful signal, the second term is the multiuser interference, and the last term is the filtered thermal noise. As stated in the introduction, the main difficulty is to find the distribution (pdf) of this multiuser interference term. In some specific case, it is possible to derive the exact pdf but in our general case, the task becomes very strenuous. Our approach is thus to establish a simplified expression of the interference to make its implementation friendly and then use it to run Monte-Carlo simulations. This semianalytic approach is halfway between a full implementation of the physical layer, which is more complex and a full analytical approach that is out of reach.

The nth interference signal at the output of the matched filter, y n (t), can be expressed as

y n (t) = r n (t) ⊗ h * (-t) = ∫ ∞ -∞ r n (x)h * (x -t)dx = N s -1 ∑ k=0 a n,k ∫ ∞ -∞ h * (x -kT s -n )h * (x -t)Ψ (x -n , n , f n , d n ) dx. (19)
It can be demonstrated with simple algebra that the Ψ function verifies the following relation:

Ψ (x -n , n , f n , d n ) = Ψ (x -n -kT s + kT s , n , f n , d n ) = Ψ x -n -kT s , n + 2 f n kT s + 1 2 d n (kT s ) 2 , f n + d n kT s , d n . ( 20
)
Using (20) and the change of variable u = x -kT sn , we can rewrite (19) into its final form:

y n (t) = N s -1 ∑ k=0 a n,k g (t -kT s -n , f n + d n kT s , d n ) × e j n +2 f n kT s + 1 2 d n (kT s ) 2 (21)
with

g(t, f, d) = ∫ ∞ -∞ h(u)h * (u -t)e fu+d u 2 2 du. ( 22
)
The function g represents one transmission channel (among L) end-to-end impulse response parameterized by the input frequency shift and drift for the considered channel. With this expression, the implementation of ( 21) is convenient because g can be computed using a discrete convolution, tabulated on a grid, and interpolated outside the grid. In practice, the integral is bounded because implementable filters have finite impulse responses.

Equation 21 enables the computation of a signal exact contribution to the interference for any filter, linear modulation, and parameter distribution. It will be used in the following to estimate the exact distribution of the interference and compare it to a Gaussian distribution.

Comparison with the Gaussian approximation

Accuracy assessment methodology

In order to compare GIA approximation with the exact interference provided by our model, we need to define a criterion. Since it is the performance of the system that we ultimately want to estimate, in our study, the accuracy is assessed by evaluating the disparity between the bit error rate (denoted P e ) computed with the exact interference provided by our model and the BER computed with the Gaussian approximated distribution (denoted Pe ). The metric chosen to measure the discrepancy between the two BERs will be the log-ratio between both dB = 10log 10 ( Pe ∕P e ). The approximation can be considered accurate when dB ≈ 0 dB. It is difficult to set an objective threshold for the accuracy since it belongs to the user of the model to define an error tolerance depending on the need, but for the sake of argumentation, we choose to consider it accurate when | dB | < 1 dB (which correspond roughly to less than 25% of relative error). Using this metric, we can also quickly assess if the distribution is overestimating/underestimating the BER by looking at the sign of dB . This criterion alone is however not sufficient since the amplitude of the interference is bounded while the probability distribution considered are not. Therefore, for asymptotic BER curves (interference to noise ratio:

INR → ∞)
, there is a signal to interference ratio (SIR) value over which P e = 0, whereas for considered distributions Pe > 0. This means that such criterion would lead to the conclusion that approximations are never accurate. For this reason, we also defined a BER threshold P thld e under which we consider that the BER is small enough to be considered as null. Said differently, we restrain the accuracy assessment in practical BER/SINR ranges. In the rest of this article, a model is considered accurate when | dB | < 1 dB for SINR values such as P e < 10 -7 . Note that, given the chosen metric, we will conduct the calculations considering an uncoded modulation.

BER computation for the exact interference model

From Equations 18 and 21, it is possible to write the lth sample of the received signal (sampling at R s = 1∕T s ) as

y(lT s ) ≜ y[l] = A 0 a 0,l + A 0 2N t ∑ k≠l k=-2N T a 0,k g ((l -k)T s , 0, 0) + L ∑ k=1 A k y k [l] + w[l], ( 23
)
where w[l] is the lth sample of noise after the matched filtering that we assume to be complex white Gaussian of variance 2 w . The impairment term is thus

z[l] = A 0 2N t ∑ k≠l k=-2N T a 0,k g ((l -k)T s , 0, 0) + L ∑ k=1 A k y k [l] + w[l]. ( 24
)
We chose to restrain the study to independent identically distributed interference with unit variance:

E A 2 k = 1, ∀k ∈ [1, L], with E[x]
being the expectation of x. The average power of the useful signal is set to E A 2 0 = P. In the following, we will write

A 2 0 = PA ′ 2 0 where A ′ 2
0 is a random variable with unit variance that has the same distribution as A k , ∀k ∈ [1, L]. Furthermore, we define the different signal/interference/noise ratios as

SIR = E A 2 0 ∑ L k=1 E A 2 k = P L (25) SINR = E A 2 0 2 w + ∑ L k=1 E A 2 k = P L + 2 w (26) INR = ∑ L k=1 E A 2 k 2 w = L 2 w . ( 27
)
The average amplitude of the useful signal P will thus serves to change the SIR while the noise level 2 w will serve to change the INR. Under these assumptions, the BER for an uncoded BPSK modulation is

P e = E A 0 Φ(A 0 ) = E A 0 Φ A ′ 0 LSINR 1 + 1 INR , ( 28
)
where Φ(x) = Pr z[l] > x is the complementary cumulative distribution function (CCDF) of the interference. The expectation is computed numerically using a Riemann sum:

E A 0 Φ A ′ 0 x = lim K→+∞ 1 ∑ i p( i ) ∑ i=0 Kp( i )Φ( i x) ( 29
) i = 0 + i k -0 K , ∀i ∈ [0, K] (30)
with K = 10 6 , 0 = F -1 p (10 -12 ), K = F -1 p (1 -10 -12 ), and where p (resp F -1 p ) is the probability density function (resp the inverse cumulative distribution function) of the random variable A ′ 0 .

BER computation for the Gaussian interference approximation

In the GIA, the impairment term is considered to have a complex Gaussian distribution with zero mean and a variance 2 z = E |z[l]| 2 . Since E y k [l] = 0 and A k , y k [l] and w[l] are independent random variables, [START_REF] Anteur | Ultra narrow band technique for low power wide area communications[END_REF] z can be written as

2 z = E |z[l]| 2 = 2 ISI E A 2 0 + y L ∑ k=0 E A 2 k + 2 w (31) with 2 y = E |y 1 [l]| 2 = ... = E |y L [l]| 2 (32) 2 ISI = 2N t ∑ k=-2N T ,k≠l |g ((l -k)T s , 0, 0)| 2 ( 33
)
2 y is the average interference power at the output of the matched filter for a unit amplitude input signal and 2 ISI is the variance of the inter-symbol interference. [START_REF] Anteur | Ultra narrow band technique for low power wide area communications[END_REF] y can rarely be written in a simple closed form expression, and it depends on the filter h and on the pdf of each of the physical parameters. In addition 2 y 0 = E |y 0 [l]| 2 = 1 in our model because we assumed perfect time, phase, and frequency synchronization for the useful signal. The variance caused by ISI, [START_REF] Anteur | Ultra narrow band technique for low power wide area communications[END_REF] ISI , is null for ideal Nyquist filters and negligible for practical filters at regular SINR. For instance,

) 34
For the GIA, the BER is given using the Gaussian CCDF:

Pe = E A 0 ⎡ ⎢ ⎢ ⎣ Q ⎛ ⎜ ⎜ ⎝ √ A 2 0 2 z ∕2 ⎞ ⎟ ⎟ ⎠ ⎤ ⎥ ⎥ ⎦ = E A ′ 0 ⎡ ⎢ ⎢ ⎣ Q ⎛ ⎜ ⎜ ⎝ √ √ √ √ 2PA ′ 2 0 2 w + L 2 y ⎞ ⎟ ⎟ ⎠ ⎤ ⎥ ⎥ ⎦ (35) Pe = E A 0 ⎡ ⎢ ⎢ ⎣ Q ⎛ ⎜ ⎜ ⎝ √ A 2 0 2 z ∕2 ⎞ ⎟ ⎟ ⎠ ⎤ ⎥ ⎥ ⎦ = E A ′ 0 ⎡ ⎢ ⎢ ⎣ Q ⎛ ⎜ ⎜ ⎝ √ √ √ √ 2PA ′ 2 0 2 w + L 2 y ⎞ ⎟ ⎟ ⎠ ⎤ ⎥ ⎥ ⎦ = E A ′ 0 Q √ 2A ′ 2 0 SINR 1 + INR 1 + 2 y INR . ( 36
)

Comparison with perfect frequency synchronization

In this paragraph, we assume that transmitters and receiver are perfectly synchronized in frequency, meaning f n = 0 and d n = 0. Perfect frequency synchronization scenario has already been studied in Chiani [START_REF] Chiani | Analytical distribution of linearly modulated cochannel interferers[END_REF] for a rectangular filter in the absence of noise. The exact PDF has been derived, and it is shown that the GIA is only accurate when the signal to interference ratios are low or for large value of L when the Central Limit Theorem applies. The accuracy of GIA for the baseline scenario is depicted on Figure 11 in terms of BER, and the corresponding dB metric is also depicted. The BER curve of the GIA is represented in black dashed line, and the BER curves for the exact interference are in straight lines for different number L of interfering signals. Using our criteria, we conclude that the GIA is inaccurate because on the SIR domain for which the Pe < 10 -7 , | dB | > 1 dB. However, as a consequence to the central limit theorem, as L increases, the SIR domain for which | dB | ≥ 1 increases.

Thus, for L large enough, the GIA can always be considered accurate.

Comparison with frequency shift and drift

In high speed mobility or in low orbit satellite communications, the frequency varies over times because of the Doppler effect, and this affects the interference distribution. The accuracy achieved with frequency drift and uniformly distributed frequency shift are depicted in Figure 12 for different satellite altitude. The results show that the frequency drift reduces the accuracy of the GIA. Indeed, as the satellite orbit altitude decreases, the satellite speed increases and the Doppler rate becomes more significant and | dB | increases. Even for large values of L, the GIA is only accurate in small SINR region. These results show that in power controlled systems, in the presence of frequency error or not, the interference cannot be considered as a white Gaussian noise except in high SIR region (typically SIR ≥ 1 dB for GEO and SIR ≥ -2.5 dB for LEO) and/or for large number of interfering signals. Furthermore, it appears that the presence of frequency errors degrades the accuracy of GIA and increase the weight of the tail.

CONCLUSION

In this paper, we first investigate the performance of UNB systems in the frame of satellite and Aloha based random access systems. We proposed a semianalytical model to derive the PLR and throughput of random access protocols even when the signal is subject to frequency drift. This model was then used to produce absolute performance on a realistic scenario. It shows that even though the frequency drift increases the number of collisions in a contention based system, the collisions are less destructive because they are localized in a smaller fraction of the message. In asynchronous random access, the reduction of multiple access interference compensates the increased number of interfering packets, and in our reference, scenario TFA shows a twofold increase of capacity over PA. These encouraging results show that the Doppler effect experienced with LEO satellite communications actually increases the performance of TFA with UNB signals, and the advantages offered by the UNB technique could thus also profit satellite networks. In future work, it could be interesting to investigate more deeply how UNB signal with frequency drift behaves with high throughput random access scheme and how advanced signal processing technique such as successive interference cancelation can improve the performance. Regarding interference, the second part of the paper investigates the legitimacy of replacing multiuser interference by a white Gaussian noise in BER analysis. A metric measuring the discrepancy of the BER predicted by the GIA and the exact BER is defined to gauge the accuracy of the GIA. By evaluating this metric taking into account the frequency shift and the frequency drift, it appears that multiuser interference can only be considered Gaussian distributed in restricted cases: low SINR, low INR, and large number of interference. It could be interesting in the future to test more scenarios, with different pulse shaping filters, power imbalance, fading, and check when multiuser interference can be considered Gaussian distributed or not. Even though the numerical results of this study are focused on Gaussian approximation, we proposed general exact model that served as a reference for BER comparisons. Furthermore, the methodology employed to assess the accuracy of Gaussian approximation can also be used to study other distributions and, in particular, generalized Gaussian distribution, -stable distribution, or distributions with heavier tail.

FIGURE 1 FIGURE 2

 12 FIGURE 1 Example of multiple UNB signals received by a LEO satellite [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 3

 3 FIGURE 3 Time/frequency detection in white, time/frequency/frequency drift detection in black on two UNB signals with frequency drifts of -400 and -200 Hz/s [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 4

 4 FIGURE 4 Physical parameters of the signal received as a function of the positions of the users under the satellite coverage. The dashed line represents the satellite radio visibility (E s ∕N 0 > 0.25 dB) and the whole surface represents the geometrical visibility [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 5

 5 FIGURE 5 Marginal distribution of the Doppler rate (left) and of E s ∕N 0 (right) [Colour figure can be viewed at wileyonlinelibrary.com]

  demodulated correctly. The obtained marginal distributions of the frequency drift (denoted p d in the following) and of the Introduction

FIGURE 6

 6 FIGURE 6 Examples of vulnerability areas for a TFA access method, a PoI whose d n = 0, and three different interfering packets whose frequency drifts d m are set to d 1 , d 2 , and d 3 [Colour figure can be viewed at wileyonlinelibrary.com]

  Figure 6 represents examples of vulnerability areas for a TFA access method, a PoI whose d n = 0, and three different interfering packets whose frequency drifts are set to d 1 , d 2 , and d 3 . It shows that the vulnerability area can be decomposed into two parts: a part formed by shapes of the same size as the packets (in plain color on the figure) and a second part resulting from the presence of a difference of frequency rate (middle part in hatched grey). With simple algebra, it can be demonstrated that the size of the vulnerability area is equal to 4BT + T 2 |d n -d m | (see Section 6). It can be shown, in the same way, that the vulnerability area would be equal to 2BT + T 2 |d n -d m | in the case of a STA access method and BT + T 2 |d n -d m | in the case of a SA access method.

FIGURE 7

 7 FIGURE 7 Frequency shift interval [f min , f max ] that leads the interfering packet (in red) to interfere with the packet of interest (in blue) [Colour figure can be viewed at wileyonlinelibrary.com]

)

  Then the minimum and maximum values of the frequency shift (f min , f max ) are computed for a message with starting instant t m and frequency drift d m to interfere with the useful packet. Finally, a frequency shift value f m is drawn, using a uniform distribution, within the interval f min , f max . Such bounds are represented on Figure7for all possible sign combinations of the time difference Δt and drift difference Ed = d m -d 0 . This figure shows a schematic time-frequency representation of packets in which the packet of interest is represented in blue while interfering packets at the collision limits are represented in red. We can note that packets are represented by parallelogram instead of rectangles because of the frequency drift. Points of collision are represented green triangle (resp circle) for the minimum (resp maximum) frequency shift value. Using the Heaviside function H(x)

FIGURE 8

 8 FIGURE 8 PLR k as a function of the number of interfering message k [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 9

 9 FIGURE 9 Schematic representation of interfering messages for different symbol rate (for fixed payload size) [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 10

 10 FIGURE 10 Comparison of the throughput and PLR of PA∕SA∕TFA∕STFA with different symbol rates in presence of frequency drift [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 11

 11 FIGURE 11 Comparison between the GIA and the exact interference model for L interfering signals assuming perfect frequency synchronization [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 12

 12 FIGURE 12 Accuracy of the GIA for L interfering signals in presence of frequency drift and uniformly distributed frequency shift (INR = 10 dB) [Colour figure can be viewed at wileyonlinelibrary.com]

4BT + dT 2 .

 2 Let us take the example of TFA access. First we can compute the vulnerability area for a PoI whose frequency drift d n = 0 using geometrical considerations. Figure13represents in blue the PoI and in red an interfering packet with a frequency drift denoted d. This allows to compute area A as a function of d and T, T being the packet duration. A is then used to compute the vulnerability area of the PoI (see Figure14): For a PoI with a frequency drift d n and an interfering packet with a frequency drift denoted d m , the vulnerable area becomes 4BT + |d n -d m | T 2 .

FIGURE 13 FIGURE 14

 1314 FIGURE 13 PoI in blue, interfering packet in red for a frequency drift denoted d [Colour figure can be viewed at wileyonlinelibrary.com]

TABLE 1

 1 Main system parameters

	Parameter	Value
	Satellite orbit	700 km
	Carrier frequency	1.6 GHz
	Satellite noise temp	450 K
	Satellite antenna pattern	Gaussian (G 0 = 4 dBi, 3dB = 110 • )
	Terminal antenna pattern 0 dBi (isotropic)
	Terminal RF power	25 mW
	Symbol rate (R Pulse shaping filter (h)	SRRCF, Roll-Off 0.5
	Coding	3GPP Turbo Code, Rate 1/3
	Demodulation threshold	0.25 dB
	Payload size (N b )	200 bits
	Total overhead	138 symbols ( 46% ): address (32 bits),
		CRC (16 bits), signaling, and pilot symbols

s = 1∕T s ) 100, 200, 400 Bd Modulation QPSK (M = 4, a n,k ∈ {±1±j})


Léa J Becker

Sarah H Journée

Pierre-Eric Lutz

Dr Ipek Yalcin
email: yalcin@inci-cnrs.unistra.fr

Comorbidity of chronic pain and anxiodepressive disorders: Deciphering underlying brain circuits

published or not. The documents may come

challenging topic by focusing on brain structures implicated in their bidirectional relationship. By reviewing lesion, pharmacologic, optogenetic or chemogenetic studies, they nicely recapitulate the role of individual brain structures (such as the prefrontal cortex -PFC, the amygdala, the nucleus accumbens -NAc, the locus coeruleus -LC, and the hippocampus) in the anxiodepressive disorders-pain dyad.

Even though studies focusing on individual brain regions improve our understanding of pathophysiology, it is also essential to reach circuitry level of understanding. Indeed, most central nervous system disorders are linked to alterations in brain circuits rather than to changes occurring in single structures. While this issue was difficult to explore until recently, the rapid development of optogenetic and chemogenetic tools now allows the targeting of anatomically defined neuronal pathways (Fig 1A). Regarding chronic pain-induced emotional dysregulation (CPED), a few preclinical studies have started to unravel critical brain circuits (Fig. 1B), which link together some of the brain regions highlighted by Bravo et al. Lee and collaborators (2015) first showed that optogenetic activation of PFC neurons terminals within the NAc blocks both nociceptive and negative affective behaviors in neuropathic animals (Fig. 1A). To further identify the neurochemical identity of this pathway, they combined optogenetic and pharmacological manipulations, and infused an AMPA receptor antagonist into the NAc, which blocked the effects of the optogenetic modulation.

The amygdala is another major hub for CPED. Llorca-Torralba and colleagues (2019) demonstrated that chemogenetic silencing of LC noradrenergic neurons projecting to the basolateral amygdala (Fig. 1A) relieved pain-induced anxiety-like behaviors in neuropathic rats, without affecting mechanical hypersensitivity.

To decipher potentially divergent contributions from subpopulations of neurons located in the same brain structure but sending projections to distinct target sites, Hirschberg and colleagues highlighted another approach (2017). The authors used retrograde viruses injected in the PFC or spinal cord (SC) to target ascending and descending projections from noradrenergic neurons of the LC, respectively (Fig. 1A).

They showed that chemogenetic activation of the LC-SC pathway had an anti-nociceptive effect, while the activation of LC-PFC projections enhanced aversive and anxiety-like behaviors. This study along with that of [START_REF] Llorca-Torralba | Chemogenetic silencing of the locus coeruleus-basolateral amygdala pathway abolishes pain-induced anxiety and enhanced aversive learning in rats[END_REF] point out that refining the analysis of specific neuronal subpopulations based on their efferences can reveal distinct and even functionally opposite effects in the modulation of the sensory and affective components of pain.

While the different studies discussed so far assessed pathways involving two brain structures, circuits underlying pain perception and emotional regulation are likely to be more complex. Recently, Zhou and collaborators (2019) dissected an ascending pathway involved in pain-induced depressive-like behavior that is composed of the dorsal raphe nucleus (DRN), the central nucleus of the amygdala (CeA), and the lateral habenula (LHb), using both chemo-and optogenetics. Their results demonstrated that inhibition of serotonergic DRN neurons terminals in the CeA (Fig. 1A) resulted in depressive-like behavior in naive mice, while activation of this pathway in neuropathic animals produced antidepressant-like and antinociceptive effects. DRN serotonergic neurons preferentially target somatostatin neurons in the CeA, which in turn project to the LHb. Surprisingly, electrophysiological recordings showed that activation of CeA somatostatin fibers elicited excitatory responses in the LHb. Indeed, even though most somatostatin neurons are GABAergic, a small population of them is in fact excitatory.

Silencing those fibers in the LHb (Fig. 1A) rescued chronic pain-induced depressive-like behavior, while activation of this pathway in naive animals conversely induced depressive-like behavior. Finally, by combining optogenetics with neuronal tracing and electrophysiology, they assessed the whole DRN-CeA-LHb circuit and showed that CPED in their model stems from the loss of serotonergic inhibition of LHb-projecting neurons of the CeA.

Taken altogether, these results deepen our understanding of the interconnected pathways that differentially modulate affective and sensory components of pain.

Nevertheless, more research will be needed to better understand the big picture of CPED. While existing studies have proven the strength of opto-and chemogenetics in combination with pharmacological or electrophysiological approaches, future studies should consider additional levels of analysis. Within this line, cell-type, pathway-specific, or single-cell RNA-sequencing all have the potential to reveal underlying molecular mechanisms. On the other hand, techniques such as fiber photometry or calcium imaging, which are now more and more used in vivo and even in freely-moving rodents, give more precise information at cellular and cell-specific levels. Finally, the use of functional magnetic resonance imaging (fMRI) or diffusion tensor imaging could give a more global insight in brain pathway alterations in CPED. Such non-invasive approaches allow for longitudinal studies and offer promising perspectives when considering that temporality is crucial when it comes chronic pain and its affective consequences.

We believe that the leap forward in novel technics and their combination with more classical approaches will lead to a better understanding of the circuits underlying the anxiodepressive disorders-pain dyad. [START_REF] Lee | Activation of corticostriatal circuitry relieves chronic neuropathic pain[END_REF]. Upper right panel: Cell-type specificity by promoter-dependent expression: injection in the PFC or SC of a retrograde virus expressing DREADD under the synthetic promoter PRS targeting noradrenergic (NA) neurons. Cell bodies of NA neurons expressing DREADD and projecting either the PFC or SC are stimulated by CNO injection in the LC through a cannula [START_REF] Hirschberg | Functional dichotomy in spinal-vs prefrontal-projecting locus coeruleus modules splits descending noradrenergic analgesia from ascending aversion and anxiety in rats[END_REF].

Figure 1 .

 1 Figure 1. Tools and approaches used to dissect neuronal circuitry underlying chronic pain-induced emotional dysregulation. A. Optogenetic and chemogenetic technics used to dissect brain circuits. Upper left panel: mPFC-NAc specificity reached by injection of anterogradely transported opsin in the PFC and by an optic fiber

  Abbreviations: 5-HT: serotonin; AMY: amygdala; BLA: basolateral amygdala; CeA: central nucleus of the amygdala; CNO: clozapine-N-oxide; DREADD: designer drug exclusively activated by a designer drug; DRN: dorsal raphe nucleus; LC: locus coeruleus; LHb: lateral habenula; PFC: prefrontal cortex; PSAM: pharmaco-selective actuator module; NA: noradrenaline; NAc: nucleus accumbens; SC: spinal cord; SOM: somatostatin; TH: tyrosine hydroxylase.

Acknowledgments

This work was supported by the Fondation de France (FdF, N° Engt: 00081244), the Fondation pour la Recherche Médicale (FRM, "Équipe FRM", EQU201903007809) and the French National Research Agency (ANR) through the Programme d'Investissement d'Avenir (under the contract "ANR-17-EURE-0022").

A

Promoter-dependent expression [START_REF] Hirschberg | Functional dichotomy in spinal-vs prefrontal-projecting locus coeruleus modules splits descending noradrenergic analgesia from ascending aversion and anxiety in rats[END_REF] Figure