## 1. Introduction

Ionospheric scintillations are fluctuations in amplitude and phase of the radio wave signal passing through the ionosphere on its path to the receiver, caused by small scale irregularities in the electron density structure [1]. Occurrence of scintillation depends on various factors, such as: solar activity, geomagnetic conditions, geo-location of the receiver, time of day, angle of signal arrival, ionospheric structure geometry and many others [1-6]. Causing signals disturbances and distortion, scintillation can significantly affect the GNSS accuracy and cause severe problems to commercial navigation systems. Scintillating signals can be classified by the intensity of fluctuations into categories of weak and strong scintillations. High and low magnetic latitudes represent the most affected regions by ionospheric irregular structures in F and E ionospheric layers. While phase scintillations are more pronounced in the sub-polar and polar regions, amplitude scintillations are significantly stronger and more pronounced in near equatorial regions [7-10].

In last five decades, statistical studies of the ionospheric scintillating radio signals brought different approaches and solutions. Some of these solutions are widely adopted and used in scintillation modelling and forecasting ionospheric dynamics, as joint Gaussian distribution of complex radio wave signal. This solution gives applicable results for weak scintillating signals [3], but in case of strong scintillations there is not an easy way to derive satisfactory results, leading to a need for further investigations [8]. While phase scintillations mainly follow the Gaussian distribution and do not represent a problem in ionospheric modelling, more complicated case is with amplitude scintillations. Revolution in ionospheric modelling and developing of the ionospheric scintillation theory has been made during '70s and ‘80s, covered mostly in papers by Rino and Fremouw [11-16]. During this period various researches had been performed providing different information on a probability distribution function (PDF) used as amplitude scintillation descriptor [11-12, 17]. Last two decades brought few intriguing studies [18-21] performed on statistics of scintillating signal leading to a need for more detailed and precise description of radio scintillation signal’s amplitude and phase PDFs for strong fluctuations.

The focus of this paper is on the statistical analysis of ionospheric scintillation in high and low latitude during strong and moderate geomagnetic activity. Gaussian and Nakagami-m PDFs of scintillating signals have been examined using real measured data and compared with theoretically derived PDFs. Further testing had been done on different data intervals, necessary for correct higher order statistical analysis, avoiding errors influence. The analysis results are presented with higher order moments, dependent on various parameters (scintillation indices, geo-location and solar/magnetic activity). Implementation of higher order moments, skewness and kurtosis, could give additional information about the ionospheric irregularities influence on the propagating signal and relation to the time delay of the signal.

## 2. Analysed data and methods

The data used in analysis have been measured by NovAtel GPS Ionospheric Scintillation/TEC Monitor (GISTM), model GSV4004B, one of three monitoring receivers installed at the Polish Polar Station in the region of the Hornsund Bay in southern Spitsbergen (approximately 77° N, 15.55° E). The second data source was derived using a Septentrio's PolaRxS Ionospheric Scintillation Monitoring (ISM) receiver located at Presidente Prudente (approximately-21.99° N, 308.59° E), Brazil, set during Concept for Ionospheric Scintillation Mitigation for Professional GNSS in Latin America (CIGALA) project. The geographical locations of both stations are indicated in the Fig. 1. The locations of the receivers made possible study of the ionospheric scintillations in the most affected regions characterized by very intense and complex behaviour of plasma - polar and auroral latitudes, due to interaction of the solar wind and interplanetary magnetic field with the Earth's magnetic field, and low latitudes, due to equatorial anomaly and electrojets. Parameters of interest in performed data analysis are signal amplitude and phase during quiet and active geomagnetic periods, recorded during measurement campaigns.

Both receivers collect raw amplitude and phase output data with a 50 Hz sampling rate, which are post-processed by removing trend and outliers, most probably caused by instrumental error. Mean and standard deviation method was used to identify and remove data outliers, defined as all data points taking absolute value greater than three standard deviations about the mean value. Following the widely used procedure proposed by Van Dierendonck [13] for detrending signal amplitude and phase, the raw data were detrended with the high pass sixth-order Butterworth filter with the cutoff frequency at 0.1 Hz. As additional conditions, study contained only measurements made at elevation angles greater than 15°, to avoid errors induced by most intense multipath, and with continuous satellite contact or a time of lock greater than 180 seconds. Computed ionospheric indices used in tracking fluctuations of amplitude and phase are S_{4}, defined as the standard deviation of the received signal amplitude during 60 seconds of measurement, normalized to the average signal amplitude, and σ_{φ}, defined as the standard deviation of the detrended carrier phase evaluated over 60 seconds.

Two case studies on statistical analysis of the scintillating signals were undertaken during strongest geomagnetic storm in 2010, on 5-6^{th} April, at Hornsund Bay, Spitsbergen, and moderate geomagnetic condition during 1^{st} November 2011, at Presidente Prudente, Brazil. Geomagnetic disturbance index Kp is used as indicator of magnetospheric and ionospheric influence to the H component of geomagnetic field, which is in close relation with the generation of ionospheric irregularities producing scintillations. Contributing parameters, as provisional 1 hour disturbance storm time (Dst) and 1 minute auroral electrojet (AE) indices are used in further data selection and in gaining detail information on all latitude magnetic activity. Data for all geomagnetic condition parameters were retrieved from the web archive of World Data Centre for Geomagnetism, Kyoto (wdc.kugi.kyoto-u.ac.jp). From Fig. 2 it is possible to notice that in case of 5-6^{th} April 2010 level of Kp index reaching a maximum value of 8, while Dst and AE reached values of-81 and 2291 nT, respectively. While in case of low latitude data Kp index achieved not so high level with maximum value of 5, whereas Dst and AE come to-72 and 1261 nT, respectively.

Analysis covered processing of 1, 3 and 12 hours data intervals for each of the tracked satellite PRNs, in order to test higher order moments on optimal amount of data mandatory for correct results. First analysis tests, not included in the paper, took data intervals of 3 and 12 hours, but due to averaging of large dataset, information from data points with pronounced scintillations had been lost and in all cases resulted to a close Gaussian PDF. Data analyses included Gaussian, and Nakagami-m PDF of scintillating signals and its dependences on higher order moments. Nakagami-m distribution is known as mostly used in the statistical analysis for description of behaviour of strong ionospheric scintillation influence on the signal amplitude. Skewness is calculated as follows (third moment about the mean divided by sample variance powered by 1.5):

and kurtosis as fourth moment about the mean divided by squared sample variance and normalized by 3:

*n* is number of values in a sample, *x*_{i} represents *i*-th value in the sample and *X* is sample mean value.

Choice of the PDFs in analyse is mainly relied on scientific and historical aspect [8, 14-20]. Nakagami-m distribution is known as mostly accepted solution in modelling of signal intensity, while Gaussian distribution is well known and widely used in statistical descriptive tool in science. In case of examination of signal intensity PDF, Gaussian distribution is used as reference distribution to measured data and Nakagami-m distribution. Gaussian distribution with zero mean, showed as acceptable, simple and effective solution for representation of the signal phase PDF for weak and moderate scintillations, is defined as

where *I* is signal amplitude, *σ* is a standard deviation and *µ*is mean value. Gaussian distribution have found use across the various scientific fields, due to simplicity of use, good approximation of variety of natural phenomena and well clarified theory behind.

Nakagami-m distribution, widely used in characterisation of ionospheric scintillation of signal intensity, radio links and wireless fading channels, is given by formula [9]:

*I* represent single value of the signal amplitude calculated per 60 seconds, function shape parameter *m* is equal to *1/S*_{4}^{2}, *Γ(.)* is Gamma function and *<. >* represents time average value of the signal amplitude. Spread parameter in Nakagami function, marked as *<I>*, is equal to signal amplitude mean value, which is set to 1. Nakagami distribution represents good approximation in describing multipath scattering with different groups of reflected wave [23].

Theoretical distribution of higher order moments, used in comparison with distribution of measured data, is calculated from the formula:

where *P(I)* is Nakagami PDF, *I* is signal amplitude and *n*=[1..4] is the order of the calculated statistical moment. From 3.381(4) formula of Gradshteyn and Ryzhik [23], integral was simplified to analytical problem.

## 3. Results

In order to obtain the information about probability distribution of scintillating signal, statistical analysis was performed by using higher order moments of signals phase and amplitude. Analysis was performed in two steps, first one included testing calculated higher order moments for received phase and amplitude for several time intervals (1, 3 and 12 hours), and the second step included comparison of the PDF of measured data with PDFs obtained from theory. Initial test cases included data set lengths of 3 and 12 hours of measurements, fulfilling the condition that amount of data points should be large enough for correct calculations of higher order moments. In these cases analysis showed good agreement with Gaussian distribution function, but only due to averaging over whole data set and losing valuable information about strong amplitude and phase fluctuations. This fact was confirmed with analysis of data samples of 1 hour length, which are shown in this paper. Additional condition is introduced, after detrending and removing outlier from measured data, only samples with more than 80% of usable data was taken in further processing.

Results of signal phase analysis for measurements in all three cases are shown in Fig. 3, with randomly chosen satellites PRN11, 4 and 2 (from left to right in figure) from group of cases with good alignment with Nakagami PDF in case of amplitude distribution, and PRN24, 29 and 32 from group showing deviant characteristic. Most of analyzed samples showed good agreement with Gaussian distribution, but as it could be seen from the last row of plots in some cases (large skewness and/or kurtosis), phase distribution deviated from Gaussian PDF.

Results of analysis for amplitude scintillation are displayed in Fig. 4. Presented are PDFs of measured amplitude data with fitted theoretical Gaussian and Nakagami distributions. As in Fig. 3 first two lines illustrates good agreement between experimental PDF and Nakagami fitted distribution. The agreement is more pronounced between 10:00 and 19:00 of local time, while the second case, showing deviations from Nakagami distribution, appears more in evening and early morning hours, especially in case of low latitudes. More detail analysis on larger data set is required for any further conclusions made on statistics of PDF deviation appearance and possible links to geo-physical parameters controlling these phenomena.

Fig. 5 depicts hourly values of higher order moments dependence on mean values of σ_{φ} and S4. All the σ_{φ} and S4 values are calculated every minute from randomly chosen PRN with full data set, which means that in specific hour PRN was visible all the time and measurements were performed without loose of lock. In graphs showing dependence from σ_{φ}, skewness (blue y-axis and marks) and kurtosis (red y-axis and marks) are mainly taking values around zero, proving good agreement with Gaussian distribution. Only in case of 5^{th} April 2010, it is possible to notice significant deviation from kurtosis, which occurs most possibly due to instrumental error. Bottom graphs displays higher order moments of amplitude, where skewness of measured data (circular markers) follows fitted distribution of moments (solid red line-skewness, and blue line-kurtosis) for smaller values of S4, while higher S4 skewness characteristic of measured data start rising much faster than theoretical one. Kurtosis controverts theoretical distribution, because of very sparse and chaotic behaviour, which could mean that Nakagami PDF is not ideal solution for representing amplitude distribution.

## 4. Conclusions

Paper presents results of a GPS signal measurements statistical analysis at high and low latitudes. Focus was on the probability distribution of phase and amplitude under disturbed geomagnetic conditions for moderate and strong scintillations. Results show, as expected, a good agreement between the measured scintillating signal phase distribution and Gaussian distribution. More interesting case is for the signal amplitude, where a systematic increase of skewness and kurtosis with the S_{4} values has been observed. In the case of weak scintillation (small S_{4}) the probability distribution of amplitude is close to the Gaussian, while for strong scintillations (large S_{4}) the skewness and kurtosis indicate considerable departure from the Nakagami distribution. This might be indicative of the non-Gaussian distribution of ionospheric electron density fluctuations [19]. Dispersion in case of kurtosis could be due to the natural spread, but further experiments are required.

Future research should consider more probability function models and χ^{2} goodness-of-fit tests for checking the precision of the distribution functions fit to the measured data. Also, more data sources, especially in low latitudes, and comparisons with in-situ measurements and other ground instruments measurements would be valuable reference in development of final model of scintillation signal distribution.

## Acknowledgments

This work was undertaken in the scope of the TRANSMIT ITN (www.transmit-ionosphere.net), funded by the Research Executive Agency within the 7th Framework Program of the European Commission, People Program, Initial Training Network, Marie Curie Actions--GA No. 264476. The author would like to thank for GPS ground data from Presidente Prudente, Brazil, obtained through the CIGALA and CALIBRA project facility, and geomagnetic indices data provided by World Data Centre for Geomagnetism, Kyoto. The author is grateful to Prof. Andrzej W. Wernik from Space Research Centre for helpful discussions and comments.