Granger causality analysis of deviation in total electron content during geomagnetic storms in the equatorial region

The changes in the solar wind and interplanetary medium’s physical conditions due to the solar activity result in several space weather phenomena such as geomagnetic storms and substorms which causes large magnetic field perturbations and disturbances in the near-Earth environment [1]. The technologies, such as the Global Positioning System (GPS), which play an important role in navigation, are severely affected by these disturbances. Therefore, it is important to mitigate the damages and errors caused due to these phenomena. Also, it is necessary to get a deeper knowledge of the physical processes responsible for generating such disturbances in the near-Earth environment and model and forecast their complex behavior. This is attempted by investigating dependencies between the parameters defining the geomagnetic storms namely disturbed storm index (Dst) and auroral electrojet index (AE) and the TEC which defines the dynamics of the ionosphere which impacts the positional accuracy. The causality between the variables is evaluated in this study. The details about the ionosphere TEC and geomagnetic indices are explained in the next section.

Ionosphere

The ionosphere is a region in the upper atmosphere which extends from around 50 to 1000 km height and is characterized by partially ionized plasma [2]. The ionosphere is described by the TEC in the layer. The TEC is the total number of electrons present along a path between a radio transmitter and receiver. It is measured in electrons per square meter. By convention, 1 TEC Unit or 1 TECU = 10¹⁶ electrons/m².

The TEC is estimated from Global Navigation Satellite System (GNSS) observables and is an important tool in studying the space weather impacts. The TEC in the ionosphere depends on the solar activity like solar flares, coronal mass ejections, high-speed solar wind, solar cycle, solar maxima, and minima. As these solar activities vary with time and have different impacts at different locations, therefore, the TEC also varies due to local time, latitude, longitude, season, geomagnetic conditions, and solar cycle and exhibits temporal and spatial variation. The TEC is found to be maximum near the equator and tapering at poles. The seasonal effects are also observed in TEC due to the movement of the Earth around the Sun. Furthermore, the electron density is linked to the 11-year solar cycle, and during this cycle, it goes through a maximum when the ionosphere is more likely to be disturbed, and the electron density much higher and unpredictable as compared to a quiet day [3]. The daily distribution of TEC also frequently gets affected by geomagnetic storms, during the high solar activity period.

The ionosphere attributes to one of the largest errors in GPS positioning. Apart from positional error, the ionosphere also causes Faraday rotation and bending of radio waves of GPS signal. The irregularity in the ionosphere also leads to rapid fluctuations in signal amplitude/phase or scintillations. The dispersive nature of the ionosphere adds to the complexity and makes the positional error dependent on the frequency of the incoming signal [4]. This dependence is described in Eq. 1 and is obtained from the Appleton and Hartree general equation for the ionosphere’s refractive index [5].

where n_iono is the refractive index, f_plasma the plasma frequency, N is the electron density, and f_signal is the incoming signal frequency. This Eq. 1 is modified and used to evaluate the group delay of a ray path crossing the ionosphere, which is given by Eq. 2

For a single-frequency receiver on L1 frequency, positional error due to group delay is minimized using a correction code. This code emulates the spatial and temporal variations and is broadcasted to the receivers. Several models have been proposed and the Klobuchar model is being currently used in GPS receivers. However, in a dual-frequency receiver, the TEC is computed at two different frequencies and error is eliminated [6]. The TEC estimation for a dual-frequency receiver at L1 and L2 frequencies is shown in Eq. 3 where L₁ is given as 1575.42 MHz and L₂ 1227.6 MHz and P1 and P2 are group path lengths.

Thus, TEC is an important parameter to understand the dynamics of the ionosphere. The TEC has a linear relationship with the positional error and 1 TECU of electron content produces a range error of 0.16 m at L1 frequency [7].

The ionosphere is one of the largest obstacles for the Global Positioning System (GPS) to become the primary navigational aid for critical applications and can cause positioning errors, which may be more than 50 m. As seen above, these errors can be eliminated in dual-frequency receivers. However, for single-frequency receivers, these errors can be only reduced by applying fixed corrections based on GNSS observables. For equatorial regions, due to equatorial anomaly and complex spatio-temporal variations in TEC, furthermore, the space weather phenomena like geomagnetic storms cause unpredictable irregularities in the ionosphere causing deviation in TEC pattern. Although there are several geomagnetic indices available to explain the strength of the geomagnetic storms, they are of little use in describing the deviation in the TEC pattern directly. Hence, there is a need to devise a method that can explain the impact of geomagnetic storms on TEC. This paper investigates the causality method to study the impact of a geomagnetic storm on TEC.

Geomagnetic storms and geomagnetic indices

A geomagnetic storm is one of the major space weather activities which affects the TEC and causes deviation in the TEC. A geomagnetic storm is a disturbance in the magnetosphere that may cause a sudden change in electron density. The Earth’s magnetosphere, thermosphere, and ionosphere are driven by the energy emitted from the Sun. The solar wind transfers its wind energy to the Earth’s magnetosphere through magnetic reconnection which leads to geomagnetic storms. These slow and fast solar winds from the coronal region also lead to powerful solar events like coronal mass ejections (CMEs) from the Sun [8] and the corotating interaction regions (CIRs). The CMEs are the result of plasma outbursts from the Sun’s active region [9]. The CMEs interact with the solar wind and interplanetary magnetic field of the Earth. The southward-directed solar magnetic field interacts strongly with the oppositely oriented magnetic field of the Earth and results in geomagnetic storms. The severe geomagnetic storms lead to anomalous changes in the ionospheric TEC, resulting in frequent amplitude and phase fluctuations. They may also cause cycle slips, amplitude, and phase scintillations, or even loss of lock. Such events not only affect the determination of the position of the receiver, but also the velocity and time of GPS receivers.

A geomagnetic storm may lead to an increase or decrease in the electron density as compared to quiet days when solar and geomagnetic activities are low. Thus, a geomagnetic disturbance may cause a positive ionospheric storm or a negative ionospheric storm. The impact of the geomagnetic storm on TEC depends on the phase and origin of the storm. A positive ionospheric storm is seen during the main phase, while in the recovery phase, negative storms are pronounced at all latitudes [10]. Positive storm effects with enhanced TEC are observed at geomagnetically low and mid-latitudes in the daytime, and negative storm effects are observed near the geomagnetic equator [11].

The TEC in the equatorial region is also impacted by the equatorial anomaly, which causes TEC accumulation at certain latitudes due to the formation of crests. This is primarily due to the equatorial electrojet (EEJ), which is caused due to vertical EXB drift leading to the fountain effect. The entire phenomenon is dependent on the EEJ and found to be more pronounced during the high solar activity period or equinox months. Hence, the geomagnetic storm effects are far more pronounced in the equatorial regions.

The strength and impact of geomagnetic disturbances are estimated using geomagnetic indices like Kp, Dst, and AE, to name a few [12]. In this study, two indices AE and Dst are selected. Both the indices are available at 1-h interval while Kp index is a 3-hourly index. Furthermore, there is a good correlation between Dst and AE; hence, AE and Dst are selected for the study. The magnitude of these indices is determined using the horizontal H component of the geomagnetic field. These indices have a pattern characteristic pattern during quiet and disturbed conditions.

The AE index characterizes the intensity of the auroral zone currents or auroral electrojet. It is the difference between the largest negative and positive H component variations, the AL and AU indices. The AE index uses magnetograms of the H component. This is collected from twelve observatories located over the longitude in the northern hemisphere at auroral or subauroral latitudes [13]. In a quiet time, this index’s value is tens of nT, and during storms and substorms, it increases to several hundred and more than a thousand nano-Tesla (nT).

The Dst index is the globally averaged value of the horizontal component of the Earth’s magnetic field at the magnetic equator from a few magnetometer stations [14]. The Dst is computed once per hour and reported in near-real-time. During quiet times, the Dst value is between + 20 and − 20 nT. Based on a geomagnetic storm’s strength, it can be classified as a moderate storm for Dst between − 50 and − 100 nT, intense for Dst between − 100 and − 200 nT, and severe or super-storm for Dst less than − 200 nT [15].

In the proposed work, an attempt is made to see if TEC can be used to study and understand the impact of space weather phenomena. This study of the dependency of TEC on AE and Dst indices can be helpful to understand the impact of space weather phenomena on the satellite-based system. The advantage of using TEC is its high temporal resolution as compared to other indices used for measuring geomagnetic storms like Dst and AE, which are available at 1-h intervals, or Kp, which is available at 3-h intervals. Furthermore, the equatorial ionosphere is characterized by large ionospheric gradients (even within 5°X 5° latitude and longitude). The deviations and perturbations in the TEC at different latitudes due to geomagnetic storms are also different. Thus, investigating the causality between the geomagnetic storm and TEC at the regional level can be useful in improving the existing methods used for correcting positional errors. This can be achieved with the high spatial resolution regional TEC data available from the GNSS receivers which have a wide global coverage. As causal inferences can result in the selection of physical quantities which are more informative, hence, the proposed study can be further combined with data-driven models for improved estimation of positional forecasting errors in the propagating signal.

Granger causality test

Causality refers to the dependency between variables and is different from correlation. Although there is a well-known correlation between variation in TEC during the occurrence of geomagnetic storms and substorms, there is no clear, direct cause and effect relation between them which can be modeled to forecast TEC.

Several attempts to forecast TEC using geomagnetic disturbances (in terms of both geomagnetic indices AE and Dst measurements) during magnetic storms and substorms have been developed using Artificial Neural Networks and linear or nonlinear regression models [16, 17]. However, most of these models are based on using a large historical dataset of these physical quantities. Many feature selection methods have also been combined with these models to identify the most relevant physical quantity. However, there is little work done in the area to identify the most informative physical quantities. Most of the studies are based on the correlation between the physical quantities which may not be very indicative due to the nonlinear and abrupt nature of these variations [18].

In a stochastic system, Granger causality between the variables can be established if it is possible for variable X_t to cause Y_{t + 1} or for Y_t to cause X_{t + 1} where t is the time variable [19]. This paper investigates the causality between the variables—deviation in TEC, Dst, and AE. The Granger causality test or G test method proposed by the Nobel Economics Prize recipient Clive W. J. Granger is used to analyze the causality between the variables.

For a time series, the Granger causality is said to exist between two variables, X and Y, if variable X can help explain Y’s future values, considering both time series are stationary or steady. Therefore, before conducting the Granger causality test, it is necessary to conduct a unit root test of the time series’ stationarity, which ensures the stationarity of the time series. The Augmented Dickey-Fuller test (ADF test) is generally used to conduct this unit root test of stationarity of the series. The Granger causality is sensitive to the lag period, and under different lag periods, completely different test results can be obtained if a precondition of stationarity is not satisfied. Thus, a series of pretests must be performed on the data before the G test.

In the present study, the deviation in TEC denoted by DTEC is taken as Y or the dependent variable, and Dst or AE are explanatory variables X₁ and X₂. The causality test is performed to check if AE/Dst can cause deviation in TEC. Hence, if X₁/X₂ does not help predict variable Y, which is DTEC, then X₁/X₂ is not the cause for the deviation in TEC. On the contrary, if the Dst or AE is the cause for DTEC, then AE and Dst should be able to predict the variable DTEC. A statistical hypothesis is tested to establish the causality. This can be explained with a mathematical formulation of the test based on vector autoregression (VAR) modeling of stochastic processes based on the past value of two variables Y and X [20]. The regression equation for two variables can be expressed as shown below:

where p is the maximum number of lagged observations; the coefficients of the model are the contributions of each lagged observation to the predicted values of X₁ (t), X₂ (t), and Y (t); and E₁, E₂, and E₃ are residuals (prediction errors) for each time series. If the variance of E₁ (or E₂/E₃) is reduced by the inclusion of the Y (or X) terms in the equation, then it is said that Y (or X) Granger-(G)-causes X (or Y). In other words, Y G-causes X if the coefficients in Aij are significantly different from zero. This is tested by performing a t-test or chi-squared test of the null hypothesis that Aij = 0, given assumptions of covariance stationarity on X and Y.

Cointegration

The cointegration test is done to establish the presence of a statistically significant connection between two or more time series. It is seen that if two variables are cointegrated, then there exists causality between variables in at least one direction [21]. Thus, a cointegration test can be viewed as an indirect test of long-run dependence. It occurs when two or more non-stationary time series have a long-run equilibrium and move together so that their linear combination of variables results in a stationary time series. There is a linear combination of the variables with an order of integration less than that of the individual series. In this context, cointegration can help understand if there is a long-run equilibrium between deviation in TEC during the disturbed condition and Dst/AE. The cointegration test establishes a stationary linear combination of time series that are not themselves stationary.

Thus, the cointegration test indicates a long-run equilibrium relationship between variables, while the Granger causality test indicates a unidirectional causality. The results of cointegration determine the type of regression model to be implemented for the causality test. The regression results with non-stationary variables can be spurious if the variables are non-stationary and cointegrated. Furthermore, the regression with the first differenced variables is for short-run relationship; hence, it cannot capture the long-run information. In such cases, the causality is investigated through vector error correction model (VECM). It is an extension of the VAR model to include cointegrated variables that balance the short-term dynamics of a process with the long-term dependencies. The VECM expresses the long-run dynamics of the process including error correction terms that measure the deviation from the stationary mean at (t−1) time. Thus, linear Granger causality on VAR can be applied only to time series that are stationary. If data are not stationary and not cointegrated, then the VAR can fit to the differenced time series. For a cointegrated non-stationary time series, with a long-term equilibrium relationship, the time series have to be fitted with the VECM model to evaluate the short-run properties of the cointegrated time series.

In this paper, three variables namely deviation in TEC, Dst, and AE are investigated under different storm conditions and for two different locations in the equatorial region. The primary aim is to identify the extent of causality and identify causal variables that can cause a state transition. As per the Granger causality principles [22], forecasting is related to identifying causal variables responsible for state transitions. Therefore, Granger causality inferences between variables can be combined with forecasting and can improve forecasting.