Precursor analysis of ionospheric GPS-TEC variations before the 2010 M7.2 Baja California earthquake

ABSTRACT Total electron content (TEC) data obtained from receivers of global navigation satellite systems (GNSSs) are increasingly being used to detect pre-earthquake ionospheric anomalies. In this study, the pre-earthquake ionospheric anomalies that occurred before the 2010 Baja California (32.128° N, 115.303° W, 4 April 2010, 22:40:42 UTC, Mw = 7.2) earthquake were investigated using TEC from Global Positioning System–Total Electron Content (GPS-TEC) data from five International GNSS Service (IGS) stations (SIO3, GOL2, QUIN, AMC2 and DRAO) near the epicentre and one IGS station (GRAZ) away from the epicentre. An analysis of the time of occurrence and duration yielded the following conclusions: both positive and negative anomalies are likely to occur and earthquake-related ionospheric anomalies occurred one to five days before the associated earthquake. The potential causes of these results were discussed. The observed anomalous variations in GPS-TEC may be attributed to the earthquake.


Introduction
Data obtained from dual-frequency global navigation satellite system (GNSS) receivers have recently been used both for ionospheric studies and for accurate positioning. Dual-frequency global positioning system (GPS) receivers are used to derive total electron content (TEC), which provide users with additional information on the upper atmosphere. The great Alaska earthquake in 1964 was the first earthquake whose seismo-ionospheric anomalies were observed (Davies & Baker 1965;Leonard & Barnes 1965). Earthquake precursors have been a research interest for a significant number of scientists (Silina et al. 2001;Liperovsky et al. 2000;Singh & Singh 2007;Zhou et al. 2009;Dogan et al. 2011;Jyh-Woei 2011;Pundhir et al. 2014). Some studies have focused on the relationship between seismic events and unforeseen variations in ionospheric activity (Parrot et al. 1993;Liperovsky et al. 2000;Afraimovich et al. 2001;Depueva 2012). Both satellite-and ground-based instruments have been used to examine seismo-ionospheric effects (Astafyeva & Heki 2011;Aggarwal 2015). Several hypotheses on the seismic and electromagnetic mechanism of geochemical and geophysical processes have also been explained by Pulinets and Boyarchuk (2004) and Namgaladze et al. (2009). Seismo-ionospheric anomalies generally occur in the D, E and F layers of the ionosphere. These anomalies may appear 1À10 days before the earthquake (Hayakawa & Molchanov 2002;Pulinets & Boyarchuk 2004;Priyadarshi et al. 2011 where f 1 and f 2 are the carrier frequencies of L1 and L2 signals of GPS satellites. A ð40:3 m 3 =s 2 ) and STEC m u indicate the slant total electron content in units of TECU ð1 TECU D 10 16 electrons=m 2 Þ on the slant signal path from the receiver u and satellite m. DCB m and DCB u are the differential code biases (DCBs) defined for pseudo-range measurements, index n denotes the time sample ranges from 1 to N (total number of time samples in a record) (Komjathy 1997;Liao 2000;Leick 2004;Dach et al. 2007;Jin et al. 2012).
VTEC can be obtained using a thin-shell approximation (Klobuchar 1986) of the single-layer ionosphere model for which the relation between STEC and VTEC is given as follows: (2) where the mapping function z m n ð Þ is defined as the satellite zenith angle at the receiver position. The mapping function M z ð Þ is defined as follows: where z 0 is the zenith angle at the ionospheric pierce point (IPP) where the line-of-sight between the satellite and ground receiver intersects the thin shell. R is the earth radius (6378.137 km), a D 0:9886 is a scaling factor from the modified single-layer mapping function (Schaer 1999) and H is the ionospheric shell height (350 km) (Mannucci et al. 1993;Langley et al. 2002;Rao et al. 2006;Spogli et al. 2013;Chakraborty et al. 2014). The calibrated STEC variations are obtained by removing the estimated DCBs from each satellite arc in equation (1). The VTEC values are determined using equation (2) for each continuous arc. In this study, the hourly VTEC values are estimated with the calibrated VTEC values at IPP over each GNSS station by fitting the approximation of second-order polynomial surfaces to IPP points for a 24-hour run (Durmaz & Karslioglu 2014): where ' IPP and s IPP are the spherical coordinates of the IPP in the sun-fixed reference frame; a 0 , a 1 , a 2 , a 3 , a 4 and a 5 are the polynomial surface coefficients. The topside hourly VTEC values at each station are calculated by substituting the sun-fixed spherical coordinates of each station alongside the estimated polynomial surface coefficients in equation (4). Insufficient data will reveal problems in estimating the hourly VTEC value. However, the missing VTEC values above the stations should also be completed by interpolating the grid TEC data obtained from IGS IONEX files (Schaer et al. 1998).

Analysis of VTEC anomalies
Investigation on the pre-earthquake ionospheric VTEC anomalies can be performed by using 15day moving median (MM) methods that employ the quartile-based statistical analysis method (Liu et al. 2009). The MM values of GPS-TEC were calculated in advance by using the quartile-based statistical analysis method. The lower quartile (LQ) and upper quartile (UQ) were also calculated.
Assuming that the GPS-TEC values are in the normal distribution with mean (m) and standard deviation (s), the expected values of MM and LQ or UQ are m and 1.34s, respectively (Klotz & Johnson 1983). The lower bound (LB) and upper bound (UB) are calculated as LB D MM ¡ 1.5 (MM ¡ LQ) and UB D MM C 1.5(UQ ¡ MM), respectively. Anomalous variations can be detected in observed GPS-TEC greater than UB or lesser than the associated LB (Liu et al. 2009). For example, the VTEC values for the first 15 days were used to generate MM, UB and LB for the 16th day. Similarly, 15 days of VTEC data between the 2nd and 16th days are used to generate bounds for the 17th day. If more than one-third of the data (e.g. eight hours are anomalous in a day) are greater or lesser than these UBs and LBs in a day, this day is taken as anomalous (Liu et al. 2009).
4. Pre-earthquake ionospheric TEC anomalies of Baja California M w D 7.2 earthquake occurred on 4 April 2010 The California region has experienced a considerable number of earthquakes over the years (http:// earthquake.usgs.gov/earthquakes/eqinthenews/). The 2010 Baja California earthquake (32.128 N, 115.303 W, 4 April 2010, 22:40:41 UTC), with a magnitude of 7.2, was the latest one to hit the Baja California area since 1902 many strong earthquakes have been reported (http://earthquake.usgs. gov/earthquakes/eqinthenews/). The Baja California earthquake was caused by the northwest-trending strike-slip (transform) faults at shallow depth along the principal plate boundary between the North American and Pacific plates at about 4.6 mm per year. The faults are distinct from but parallel to the strands of the San Andreas fault system (http://earthquake.usgs.gov/earthquakes/eqinthenews/).
In this study, GPS-TEC data were used to investigate ionospheric abnormal behaviours prior to the 2010 Baja California earthquake. All maps, graphics, calculations and some statistical procedures were performed by MATLAB® scripts. For this earthquake, the radius of the earthquake preparation zone (EPZ) was estimated using the Dobrovolsky formula: r D 10 0:43ÃM , where r is the radius of the EPZ (km) and M is the magnitude of the earthquake on the moment magnitude scale (Dobrovolsky et al. 1979). The radius of the EPZ for the Baja California earthquake was calculated to be 1247.38 km from the epicentre of the earthquake. IGS stations in the vicinity of the earthquake zone are SIO3 (32.859 N, 117.249 W), GOL2 (35.425 N,116.889 W), QUIN (39.974 N,120.944 W) and AMC2 (38.803 N,104.524 W). According to Oyama et al. 2008, the disturbed region extends across a wide range of latitudes and longitudes such as 60 in the eastÀwest direction and 40 in the northÀsouth direction for large earthquakes of magnitude (M > 6). Hence, in this study we also selected one IGS station outside of the EPZ (DRAO (49.322 N, 119.625 W)). The earthquake (EQ) is a localized phenomenon and space weather conditions are a global phenomenon; therefore, we have selected one IGS station (GRAZ (47.066 N,15.493 E)) away from the earthquake epicentre to show a clear picture (Figure 1). The IGS RINEX files, SP3 (precise satellite orbits) files and IGS IONEX (ionospheric TEC maps and satellite DCBs) files were obtained from the GNSS data and products archive of Crustal Dynamics Data and Information System to investigate the anomalous variations on ionospheric TEC during the earthquake. These data are available from the website ftp://cddis.gsfc.nasa.gov. Hourly GPS-TEC data for the days between 20 March and 18 April 2010 were processed by following the methodology explained in Section 2. TEC data were analysed using the 15-day MM method analysed using the methodology explained in Section 3.
In order to check the reliability of DCB estimations, the estimated DCBs of AMC2 station were compared with the estimates in the IGS IONEX files. AMC2 station was used because common DCB values of the same station are included in the IGS final IONEX file. The DCBs of other three stations exist in the IGS final IONEX file discontinuously. The differences between the estimated daily receiver DCBs and IGS IONEX DCBs of AMC2 station are shown in Figure 2. In addition to this comparison, the reliability of the estimated GPS-TEC values of all IGS receivers were checked by comparing their values with GIM-TEC values obtained from IGS IONEX files. The daily mean values and RMS of differences between the GPS-TEC and IGS IONEX file values for IGS stations were then calculated and only AMC2 IGS station result was shown in Figure 2.
A comparison between the DCBs estimated from the AMC2 station and the DCBs from the IGS IONEX files revealed a difference of less than 1 ns. The calculated RMS (root mean square) of differences of the estimated AMC2 IGS receiver DCBs and value of the IGS IONEX file is §0.29 ns. The results indicate that the estimated AMC2 site receiver DCBs may be considered reliable. The accuracy of IGS IONEX TEC values was defined in the range of 2À8 TECU (https://igscb.jpl.nasa.gov/ components/prods.html). The results in Figure 2 show that the calculated mean and RMS of GPS-TEC values are within the acceptable range and hence can be considered reliable.
The ionospheric parameters are mainly affected by solar geophysical conditions and geomagnetic storms, especially in the equatorial and polar regions. When discussing the relationship between ionospheric anomaly and earthquake, the solar terrestrial environment must be considered to exclude anomalies that may have been caused by solar or geomagnetic activities (Pulinets et al. 2003;Ulukavak & Yalcinkaya 2014c, 2016b. Geomagnetic and solar indices (i.e. Dst, Kp, F10.7 and EUV) are obtained for space weather conditions. We accessed the geomagnetic storm indices (Dst) from the archive of the Data Analysis Center for Geomagnetism and Space Magnetism Graduate School of Science of Kyoto University via the link http://wdc.kugi.kyoto-u.ac.jp/dstae/index.html. The Kp indices data from the National Oceanic and Atmospheric Administration are archived in ftp://ftp.ngdc.noaa.gov/STP/GEOMAGNETIC_DATA/INDICES/KP_AP/. The solar radio flux (F10.7 cm) data from the Canadian Space Weather Forecast Center in Ottawa from http://www. spaceweather.gc.ca/solarflux/sx-en.php and EUV Flux were obtained from the University of Southern California Space Sciences Center archived from http://www.usc.edu/dept/space_science/semdata folder/semdownload.htm. In this study, we investigated the geomagnetic and solar indices (Dst, F10.7, EUV and Kp) to distinguish seismic anomalies from other anomalies related to space weather conditions during the period from 20 March to 18 April 2010 (Figure 3). Kamide et al. (1998), Rozhnoi et al. (2004) and Contadakis et al. (2012) established three classes of geomagnetic storms according to the magnitude of the storms, namely weak (Dst min > ¡50 nT), moderate (¡50 nT > Dst min > ¡100 nT), and intense (Dst min < ¡100 nT). In Figure 3(a), positively peaked Dst index values before the earthquake was C20 nT on 25 March, and negatively peaked Dst index values after the earthquake were ¡55 nT (moderate storm occurrence) on 5 April, ¡81 nT (moderate storm occurrence) on 6 April, ¡54 nT (moderate storm occurrence) on 7 April, and ¡67 nT (moderate storm occurrence) on 12 April 2010, respectively. In this regard, it indicates that ionosphere was quiet before the earthquake (26 MarchÀ4 April) and after the earthquake even ionosphere was disturbed due to the geomagnetic activities (5À7, 12 April). Vitinsky et al. (1986), Bruevich et al. (2014) and Coley et al. (2014) indicated the possibility of classifying solar activity in terms of the F10.7 index. In Figure 3 The Kp indices range from '0' (very quiet) to '9' (extremely disturbed). Kp indices between '0' and '3' values were classified as very low and quiet geomagnetic activity, respectively. Kp indices between '4' and '5' were classified as active and moderate storms, whereas Kp indices between '5' and '9' were classified as strong storms, respectively (http://www.spaceweatherlive.com/en/help/the-kp-index). In Figure 3(c), Kp index values were very high after the earthquake, reaching 7.7 (severe storm), 5.3 (moderate storm) and 5.7 (moderate storm) on 5, 6 and 12 April 2010, respectively. According to Figure 3(c), space weather conditions were very low between 20 March and 4 April 2010, and the earthquake occurred during geomagnetic quiet conditions. After the earthquake, ionosphere was disturbed due to the geomagnetic storms (5, 6 and 12 April).
All index (Dst, EUV, F10.7 and Kp) values are evaluated together, the space weather conditions have been identified as quiet between 26 March and 4 April before the earthquake.
The LB, UB and MM values were calculated using daily variations in the one-hour resolution of GPS-TEC data calculated from selected five IGS sites close to the earthquake epicentre and one IGS station away from the epicentre. IGS stations arranged with the increasing distance from the epicentre in Figure 4.
The earth's ionosphere is affected by various facts. Geomagnetic activity, solar activity, meteorological events and human-induced effects influence the earth's ionosphere from space to ground (Aggarwal 2015). The ionosphere also presents normal day-to-day, seasonal and diurnal variations making it difficult to identify possible pre-earthquake ionospheric anomalies as it changes the ionospheric parameters, such as GPS-TEC, on a regional scale .
The GPS-TEC time-series data shown in Figure 4 clearly indicate that the obvious first positive anomaly occurred »10 days (on 25 March 5) before the earthquake at all IGS stations near the earthquake epicentre. The second positive anomaly observed only away from the epicentre of the earthquake at GRAZ IGS station on 5 April. The third positive anomaly observed only at QUIN IGS station on 7 April. The fourth positive anomaly observed at all IGS stations on 12 April in the EPZ. GPS-TEC values were 2À8 TECU higher than the UB on certain days (25 March, 5, 7 and 12 April).  and Astafyeva and Heki (2011) indicated that the intensity of the solar UV radiation, which may cause increase in daily GPS-TEC values after the peak of the EUV flux was observed. Higher EUV flux raises the photoionization effect and has an important role in the production of the ionospheric anomalies. Photoionization can produce more electrons and therefore enhances the electron density in the ionosphere (Wu et al. 2004). In this study, two peak values observed in the variations of EUV flux values (Figure 3(b)). The first peak occurred on 21 March and the second smaller one observed on 25 March, which may cause increase in daily TEC values 4À5 days later . In this regard, we can expect the causes of EUV flux increase in GPS-TEC between 21 March and 25 March. Indeed the reason of the first positive anomaly could be the intensity of the solar UV radiation, which is seen on 25 March. Significant enhancements of GPS-TEC may be correspond to the mentioned positive Dst event increased the GPS-TEC instantaneously (Mendillo 2006). In this regard, another possible reason of the GPS-TEC enhancement on 25 March may be the positive Dst event occurrence (Figure 3(a)). The reason of the second, third and fourth positive anomalies could be affected by continuity of the geomagnetic storms which is explained by Mendillo (2006). Mendillo (2006) investigated the TEC disturbance patterns during various geomagnetic storms to describe seasonal and solar cycle effects. In his study, most space weather conditions have found that positive ionospheric anomalies in the electron density caused by geomagnetic storm effect occurred within the first 24 hours after onset while negative ionospheric anomalies in the electron density may persist for up to several days explained by Mendillo (2006). As can be seen from our results, the second, third and fourth positive anomalies occurred on 5, 7 and 12 April caused by geomagnetic storm effect occurred on 5 April 2010.
The negative anomalies started to occur observed at SIO3, AMC2, GOL2, QUIN and DRAO IGS stations one to five days before the earthquake under the quiet space weather conditions. Figure 4 clearly shows that VTEC values were 2À8 TECU lower than the LB on certain days (30 March to 3 April) before the earthquake. The first negative anomalies can be clearly observed from the QUIN IGS station on 30 March, GOL2 and AMC2 IGS stations on 31 March, DRAO IGS station on 1 April and SIO3 IGS station on 2 April, respectively. Negative anomaly could not be observed before the earthquake over GRAZ IGS station which is away from the epicentre. Figure 3(a) and 3(c) clearly shows that geomagnetic activity and geomagnetic storm started to occur the day after earthquake on 5 April 2010. Hence, negative anomalies observed after 4 April may be attributed to geomagnetic storm as shown in Figure 4. The reason of the positive and negative anomalies after the earthquake (5, 6, 7 and 12 April) could be affected by continuity of the geomagnetic storms which is explained by Mendillo (2006).
According to Liu et al. 2009, the significant anomalous day was determined from the percentage of number of anomalous hours in one-day GPS-TEC anomalies. The percentage limit of anomalous GPS-TEC was determined as »33% (eight hours are anomalous in a day) ( Figure 5). In Figure 5, the percentage of significant positive and negative anomaly days were marked separately. Figure 5 shows the percentage of first positive anomalous GPS-TEC variations with 33.33% over AMC2 and with 66.67% over DRAO, on 25 March 2010. The abnormality percentage is high only over DRAO IGS station than the other IGS stations could be the reason of the latitude of this station, because the geomagnetic activity effects on ionosphere increase from the low to high latitudes on the earth (https://www.spaceweatherlive.com/en/help/the-low-middle-and-high-latitude). The first negative anomalous GPS-TEC variations over IGS stations arranged with the increasing distance from the epicentre with 33.33% over SIO3 on April 2, with 33.33% over QUIN on 30 March, with 58.33% over GOL2 on 31 March, with 33.33% over AMC2 and with 50.00% over DRAO on 1 April, respectively, and the last positive anomalous GPS-TEC variation with 33.33% over the GOL2 and AMC2 stations on 12 April 2010. All percentages of significant positive and negative anomaly days were given in Table 1. Table 1 and Figure 5 show that significant negative anomalous days over IGS stations were determined in quiet space weather conditions before the earthquake. Negative anomalous days were determined on SIO3 IGS station two days, on GOL2 IGS station four days, on QUIN IGS station five days, on AMC2 and DRAO IGS stations three days before the earthquake with the increasing distance from the epicentre, respectively. The percentages calculated between 26 March and 4 April 2010 (in quiet space weather conditions) may be more related to the 2010 Baja California earthquake than space weather conditions.

Conclusion
The earth's ionosphere is affected by various facts, from space to ground due to the variation in geomagnetic activity, solar activity, meteorological events and human-induced effects (Aggarwal 2015). The ionosphere also presents normal day-to-day, seasonal and diurnal variations making it difficult to identify possible pre-earthquake ionospheric anomalies as it changes the ionospheric parameters, such as GPS-TEC, on a regional scale . In this work, the effect of GPS-TEC day-to-day changes before strong earthquake is analysed. We studied the response of the ionosphere to the large earthquake that occurred in Baja California, on 4 April 2010 of magnitude 7.2, using GPS-TEC data obtained from IGS network stations. The results indicated that anomalous variations of the GPS-TEC were detected during the one to five days before the 2010 Baja California earthquake. The two main well-known different physical quantities of the ionosphere are the VTEC and the critical frequency of the F2-layer (foF2 or NmF2). The NmF2 is the maximum electron density at the F2 layer of the ionosphere, approximately 300À500 km altitude, while VTEC is the integration from a ground-based GNSS receiver to the GPS satellites orbit approximately 20,000 km altitude (Liu et al. 2004a). Liu et al. (2001) emphasized that NmF2 and GPS-TEC had similar susceptibility. They are highly correlated and GPS-TEC can also be used to obtain seismo-ionospheric precursors if they exist. In this study, we investigated the ionospheric VTEC anomalies by analysing the possible causes of these ionospheric anomalies based on space weather conditions. Solar activity (EUV and F10.7), geomagnetic activity (Dst) and geomagnetic storm (Kp) indices were examined. The indices showed that the status of the space weather and magnetic field conditions were very quiet during 26 March and 4 April 2010. Hence, we can conclude that the 2010 Baja California earthquake occurred under quiet space weather conditions. The GPS-TEC values between 30 March and 4 April 2010 were approximately 2À8 TECU lesser than the LB values calculated for the same days. Similar negative anomaly conclusions were drawn in several studies (Liu et al. 2004a(Liu et al. , 2004b(Liu et al. , 2009) who have shown the significance of negative TEC anomalies occurred in Taiwan and China. Liu et al. 2004a observed negative TEC disturbances (TEC decreases) for Taiwan earthquake. In their study, the anomalous decreases in the GPS-TEC were found one to five days before the earthquakes. The equatorward displacement of the equatorial northern crest and the reduction of the electron concentration in the crest were probably the reasons for such effects. In our results presented above, we have emphasized the negative anomalies (depletions in TEC before 4 April 2010) as a result of large 2010 Baja California earthquake. In this regard, the abnormal change in GPS-TEC variations one to five days before the earthquake may be considered as the ionospheric precursor of the 2010 Baja California earthquake. Anomalous variations after the 2010 Baja California earthquake were related to occurrence of the geomagnetic storms (on 25 March, 5, 6, 7, and 12 April 2010) in the agitated space weather conditions. The development of GNSS technology and assimilation techniques allows us to improve our understanding of the ionosphere. In the near future, GPS-TEC will be increasingly utilized to discover the mechanism of seismo-ionospheric coupling.