Earthquake-induced precursory signals are always a topic discussed by scientists. Conventional seismic monitoring techniques along with knowledge of earthquake-associated signals had widely utilized to disclose an association between anomalous geophysical phenomena and earthquakes. However, due to the complexity of earthquake preparation process, identifying the real earthquake precursors from numerous and complicated geophysical signals is still a challenging problem.

Many reports have presented the possible associations between anomalous electromagnetic emissions in a broad frequency range (DC-VHF) and earthquakes^[1~13]. The generation of seismo-electromagnetic emissions has been verified by the experim ents^[14~21]. In the recent 40 years, the ionospheric disturbances possibly associated with seismic activity have also increasingly attracted attention of geophysicists. Since the Great Alaska Earthquake in 1964, the anomalies of the Vertical Total Electron Content (VTEC) or the peak electron density N_mF₂ a few days before strong earthquakes have extensively been reportedE^[22~36]. In these seismo-ionospheric anomalies reported, the VTEC or N_mF₂ around the epicenters showed anomalous enhancement or decrease several days before earthquakes. In identifying the disturbance source of ionospheric anomalies, some scholars gave valuable reference information, which provide the possibility for validating the seismo-ionospheric signals^{[31, 32, 37~41]}. Besides, many theoretical and experimental studies have been made, and some possible physical mechanisms of seismo-ionospheric effects have also been put forward. Among these possible physical mechanisms, the enhanced lithosphere-·atmosphere-ionosphere (LAI) coupling and acoustic gravity waves (AGW) generation in earthquake preparation zone are the most popular interpretations of seismo-ionospheric effects^{[32, 42~51]}.

In this paper, we will use GPS technique to examine the ionospheric variations during the M_s8. 0 Wenchuan earthquake (epicenter: 31. 0° N, 103.4° E) which occurred at 06: 28UT (14 : 28LT) on 12 May 2008. in Sichuan province, China. As one kind of ionospheric sounding technique with high efficiency and low cost. GPS has been widely used in space weather research. Owing to the utilization of two ultra high frequency waves. the VTEC derived from GPS measurements is a feasible parameter in investigating the seismo-ionospheric disturbances. Moreover, the increasing number of ground-based GPS receivers provides the possi-bility of simultaneously, continuously and precisely monitoring the variations of ionospheric VTEC^[52~54]. Therefore, if the pre-earthquake ionospheric anomalies do exist. GPS VTEC technique should have the ability to detect them.

In the following parts. we will introduce the procedure to derive GPS VTEC and examine the ionospheric VTEC variations during Wenchuan earthquake. Further, the features of the pre-earthquake ionospheric anomalies are analyzed and the sources probably related with the ionospheric anomalies are also discussed.

Total Electron Content (TEC) is an important ionospheric parameter for space weather research and defined as the line integral of electron density along the propagation path of the radio waves from satellite to ground receiver. Since sionosphere is a dispersive medium, we can derive the TEC from the geometry-free linear combination of phase-smoothed pseudo-range P₁/P₂^{[55, 56]}:

(1)

(2)

(3)

where ϕ_i(i=1, 2) is the carrier phase measurement after gross error is eliminated and cycle slip is repaired, and is the mean difference between all the accepted code and phase measurements in a continuous observation arc.

Seen from the above, the derivation of TEC from phase-smoothed pseudo range is rather simple. However, the TEC contains the differen-tial instrumental biases introduced by GPS satellite and receiver. Therefore, the estimation of the differential instrumental biases is vital in calcula-ting precisely the real TEC. In order to estimate the instrumental biases, a single-layer ionospheric model (VTEC model) in a sun-fixed coordinate system is usually introduced^[56].So, the TEC from formula (1) may be written as:

(4)

where is a mapp-ing function with the elevation angle e(satellite seen from GPS station), R is the mean radius of the Earth and H is the effective ionospheric height, while ((b_r+b_s) is the instrumental biases of satellite and receiver pair at L1/L2 (1ns of instrumental biases=2. 852TECu). The detailed aspects about the estimation of the instrumental biases can be found in the related papers^[57~59]. Generally, the instrumental biases are very stable over time on a scale of days to months so that they can be treated as constants during ionosphere modeling^{[55, 60]}. In this paper, adopting the vertical TEC model of spherical harmonic expansion, we Used the least square fitting to estimate the instrumental biases. After the real TEC is inversed, with the single-layer ionospheric model which makes an assumption that all the electrons are concentrated in a thin spherical layer at a height of 350 km above the Earth's surface. we can convert the real TEC into VTEC at the Ionospheric Pierce Point (IPP) which is the intersection of the receiver-satellite line-of-sight with the single-layer spherical shell. After all IPP VTECs were computed, the VTEC over any given point in the observation area can be obtained through the inverse distance squared weighting method^[61]. For reducing the effects of large ionospheric horizontal gradients, GPS signals at elevation angles less than 50 degrees will not be considered. After an observation session, the time series of VTECs over any given point in the observation area will be obtained.

Due to the comparability of daily variations in the ionosphere, if there are not sudden changes in the solar-terrestrial environment, the VTEC over a given point could not generate large day-to-day variations in a short period. So, under the assumption of a norma]distribution for the VTEC, we computed the mean value of the previous 10-day VTECs and the associated standard deviation σ to construct the upper bound u+2σ and the lower bound u-2σ. During the process of identifying ionospheric anomalies, the upper/lower bounds form the normal reference threshold. If an observed VTEC in the following day goes beyond either the associated lower or upper bound, we declare that a negative or positive ionospheric anomaly is detected.

In this part, using continuous GPS stations the observations 0f 25 in China (Fig. 1), we examined the variations of ionospheric VTEC during W enchuan earthquake. The sampling rate of each GPS receiver is 30 seconds. The unit of VTEC is TECu (1TECu=10 ¹⁶electrons/m²) in the present paper.

图 1

Fig. 1

Fig. 1 The distribution of the GPS stations (black hollow triangles). The black pentacle denotes the epicenter of Wenchuan earthquake

In view of the strong influence of solar and geomagnetic activities on the ionosphere, the solar-terrestrial environment should be taken into account in identifying the disturbance source of ionospheric anomalies. Fig. 2 presents the variations of solar F10. 7 flux and geomagnetic activity indices (K_p, Dst) during the earthquake. One can see that the solar activity was rather quiet within 18 days before the earthquake, while a series of geomagnetic disturbances and three weak geomagnetic storms (K_p≥4) were recorded from April 24 to May 6, 2008. After May 7, the geomagnetic activity was low and rather quiet (K_p≤2, Dst≥-20 nT).

图 2

Fig. 2

Fig. 2 From April 24 to May 17, 2008 (UT) : (a) The solar F10. 7 flux; (b) The geomagnetic Kp index; (c) The geomagnetic equatorial Dst index. The vertical dashed line indicates the tme of Wenchuan earthquake.

Using the method mentioned above, we computed the time series of ionospheric VTECs over LUZH GPS station which is the nearest one to the epicenter of Wenchuan earthquake. Fig. 3 demonstrates the statistical analysis results of VTEC variations during Wenchuan earthquake. Seen from Fig. 3, the VTEC showed the negative or positive anomalies 13, 9, 6 and 3 days before Wenchuan earthquake. Especially on May 9, namely 3 days before Wenchuan earthquake, the positive anomalies are the most remarkable. Moreover, two positive anomalies were detected on May 3, 2008.

图 3

Fig. 3

Fig. 3 The tme series of VTEC over LUZH GPS station from April 25 to May 14, 2008 (UT). The red/blue and black curves represent the upper / lower bounds and observed VTEC, respectively. The black arrows point to the anomalies. The little black circle dot points out the time of Wenchuan earthquake

In order to display these anomalies of VTEC more roundly and clearly, taking the corresponding upper or lower bound as a reference, we created the di ferential two-dimensional VTEC (ΔVTEC) maps around the epicenter in the periods of ionospheric VTEC anomalies. Fig. 4, 5, 6, 7 and 8 describe the distributions of ΔVTEC 13, 9, 6 and 3 days before Wenchuan earthquake, respectively. Here, ΔVTEC stands for the absolute magnitude with which the observed VTEC goes beyond the upper or lower bound. ΔVTEC < 0 stands for the negative anomalies. ΔVTEC=0 implies that the observed VTEC just lies between the upper bound and lower bound. ΔVTEC > 0 stands for the positive anomalies.

图 4

Fig. 4

Fig. 4 The ΔVTEC distribution around the epicenter of Wenchuan earthquake at: (a)06:00UT; (b)07:00UT; (c)08:00UT; (d)09:00UT; (e)10:00UT on April29, 2008.The black pentacle represents the epicenter

图 5

Fig. 5

Fig. 5 The ΔVTEC distribution around the epicenter of Wenchuan earthquake at: (a) 03 : 00UT; (b) 04 : 00UT; (c) 05 : 00UT; (d) 06 : 00UT; (e) 07 : 00UTon May 3, 2008. The black pentacle represents the epicenter

图 6

Fig. 6

Fig. 6 The ΔVTEC distribution around the epicenter of Wenchuan earthquake at: (a) 13 : 00UT; (b) 14 : 00UT; (c) 15 : 00UT; (d) 16 : 00UT; (e) 17 : 00UT; (f) 18 : 00UT on May 3, 2008. The black pentacle represents the epicenter

图 7

Fig. 7

Fig. 7 The ΔVTEC distribution around the epicenter of Wenchuan earthquake at: (a)06:00UT; (b)07:00UT; (c)08:00UT; (d)09:00UT; (e)10:00UTon May6, 2008. The black pentacle represents the epicenter

图 8

Fig. 8

Fig. 8 The ΔVTEC distribution around the epicenter of Wenchuan earthquake at: (a)08:00UT; (b)09:00UT; (c)10:00UT; (d)11:00UT; (e)12:00UT on May9, 2008.The black pentacle represents the epicenter

From Fig. 4, 5, 6, 7 and 8, we see clearly the amplitude-size and spatial coverage of the VTECs below the lower bound or beyond the upper bound. One obvious feature is that the VTEC anomalies drifted from east to west with the rotation of the Earth in all the 4 anomalous days. In other words, the area of the VTEC anomalies in any one of the 4 days almost located over the same place in the sun-fixed reference frame. Moreover, relative to the non-perturbed state of the ionosphere, at the tive peak times of ionospheric anomalies, the maximum variability amplitude of the ionospheric VTEC around the epicenter in the 4 anomalous days exceeded 40%, , 40%, , 80%, 40%, and 70%, (Fig. 9) respectively. Besides, all the ionospheric anomalies lasted beyond 4 hours in the 4 days. However, due to the limited observation area, we can't determine exactly how long the onospheric VTEC anomalies lasted.

图 9

Fig. 9

Fig. 9 In comparison with the mean value of the previous 10-day VTECs, the distributions of the relative variability of the VTEC around the epicenter of Wenchuan earthquake at five peak times: (a) 08 : 00UT, April 29, 2008; (b) 05 : 00UT, May 3, 2008; (c) 16 : 00UT, May 3, 2008; (d) 08 : 00UT, May 6, 2008; (e) 10 : 00UT, May 9, 2008. The black pentacle represents the epicenter

On the basis of the above analysis, we may declare that the ionospheric VTEC displayed the anomalous variations before Wenchuan earthquake. However, are all these ionospheric anomalies related with Wenchuan earthquake? To search the source of the ionospheric disturbances, we reviewed the level of solar and geomagnetic activity. Seen from Fig. 2, the solar-terrestrial environment was quiet after May 7, so the ionospheric anomalies on May 9 couldn't come from the sudden variations of the solar or geomagnetic activities. Moreover, the ionospheric anomalies have a local effect at the given time point during the anomalous period. Generally speaking, the low and relatively quiet activities of the sun and geomagnetism could not cause the large ionospheric disturbances. Therefore, the positive anomalies on May 9 could be related to Wenchuan earthquake. As for the ionospheric anomalies on the other 3 days, they just appeared in the period when a series of geomagnetic disturbances were happening. From Fig. 4, 5, 6, and 7, all these ionospheric anomalies look like local effects at the given time point. In a few days before April 29, although strong geomagnetic disturbances (3 < K_p < 4) happened, it still did not reach a level of weak storm. Moreover, the geomagnetic conditions were quiet (Dst≥-20 nT, K_p≤3 and the daily sum of K_p is lower than 20) as a whole on April 29, so it may be difficult to infer that the ionospheric anomalies on April 29 were a result of geomagnetic-effects. However, the anomalies appeared 13 days before Wenchuan earthquake, which is far away from the time scale (0~5 days) of seismo-ionospheric signal summarized by Pulinets et al. ^[32] and Liu et al.^{[28, 62]}. Therefore, it can not be determined whether these ionospheric anomalies were associated with the earthquake. On May 3, the ionospheric anomalies (Fig. 5) in the afternoon is probably related to the Wenchuan earthquake because the geomagnetic conditions were almost quiet (K_p≤3 and Dst≥-20 nT) during the ionospheric anomalies and a few hours before, while the ionospheric anomalies (Fig. 6) in the night had a long-span spatial effects, especially along the east-west direction, so it could be caused by the continuous geomagnetic disturbances (2 < K_p≤3 and Dst < -20 nT) on May 3. On May 6, the geomagnetic activity reached a level of weak storm (K_p≥4 and Dst < -20 nT) and the ionospheric anomaly had a large spatial distribution (Fig. 7), therefore, the ionospheric anomalies on that day may be a result of geomagnetic effects.

In addition, due to the thermospheric dynamics caused by various types of atmospheric disturbances from lower atmosphere, the VTEC may generate signiticant day-to-day variations^[63]. However, compared with a non-perturbed state of the ionosphere, the relative variability associated with thermospheric dynamics usually does not go beyond 30%^{[32, 64~66]}. However, seen from Fig. 9, the maximal relative variability of the VTEC obviously went beyond the limk of 30%, during the ionospheric anomalies in all the 4 days, in comparison with the mean value of the previous 10-day VTECs. Moreover, the spatial sizes of the VTEC decrease and enhancement by more than 30%, extended larger than 1500 km in latitude and 4000 km in longitude at the five peak times of the ionospheric anomalies. According to the large relative variability of VTEC, we may exclude the influence on the ionosphere from meteorological effects.

Summarizing the above analysis, we can draw a conclusion that the positive and negative VTEC anomalies had appeared before Wenchuan earthquake, while the positive anomalies in the afternoon of May 3 and on May 9 are possibly the seismo-ionospheric signals.

Analyzing the ground-based GPS VTEC data, we found the positive and negative ionospheric VTEC anomalies before Wenchuan earthquake. These ionospheric anomalies appeared 13, 9, 6 and 3 days before Wenchuan earthquake. In the 4 days of ionospheric anomalies, the spatial sizes of the VTEC decrease or enhancement by more than 30% were larger than 1500 km in latitude and 4000 km in longitude. Moreover, the ionospheric anomalies lasted beyond 4 hours and drifted from east to west in any one of the 4 days. In view of the solar-terrestrial environment and the characters of the ionospheric disturbances, we think that the VTEC anomalies in the afternoon of May 3 and on May 9 around the epicenter are most likely the ionospheric precursors of Wenchuan earthquake.

As for the mechanism of the seismo-ionospheric effects, two popular theories were put forward. The first one is that the acoustic gravity waves (AGW) from the lithospheric gas discharge in the vicinity of an active tectonic fault penetrate into the ionosphere to cause ionospheric oscillation^[43]. However, the theoretic calculations^[67] and experi-mental result^[24] of the M=6. 7 Northridge earthquake of 17 January 1994 showed that the VTEC variations associated with AGW were 2~2. 5 order of magnitude lower than the background one. Based on this point, Pulinets and Boyarchuk ^[32] thought the penetration of AGW into the ionosphere is quite improbable and its precursory abilities may be neglected. The second one is that one extra electric field on the Earth's surface penetrates into the ionosphere and then generate the E×B drift. In earthquake preparation zone, the increasing radon emanation from active faults and cracks before earthquakes can ionize the near ground atmosphere producing a large vertical electric field on the Earth's surface^[32]. It was also experimentally validated by Freund^[68] that the mobile positive holes can be activated in the crust by microfractures during the dilatancy stage of earthquake preparation, and the diffusion and outflow of these holes can generate high electric fields on the Earth's surface. Through a quasi-electrostatic model of atmosphere-thermosphere-ionosphere coupling, Kim and Hegai^[69], Pulinets et al. ^[70] and Sorokin et al. ^[71]showed that a strong vertical electric field on the Earth surface could penetrate into the ionosphere and modify its dynamics and electron density distribution prior to the earthquake onset.

Heretofore, a few scholars had also studied the seismo-ionospheric effects of Wenchuan earthquake. Zhao et al. ^[72] analyzed the N_mF₂ from ionosondes and the VTEC from GPS measurements. and found the ionospheric anomalies over the epicenter of Wenchuan earthquake on May 2, 3, 6, 7 and 9. Examining the solar-terrestrial environ-ment and characters of ionospheric anomalies, they thought the enhancement of electron density or VTEC on May 2 and 3 could be caused by the geomagnetic disturbances through a global-scale disturbance map (RTEC, subtract the 12-day median value) which is formed by singling out the maximum disturbed TEC at every grid point during 10 hours, and the decrease of electron density on May 6 and 7 is not clear whether they are related to the storms or the Wenchuan earthquake, while the anomalous enhancement accompanied by its coun-terpart in the southern hemisphere on May 9 is probably associated with seismo-ionospheric coup-ling process. Later, according to the global ionosphere map (GIM, ftp://cddisa.gsfc.nasa.gov) and the ionospheric occultation observation of FORMOSAT 3/COSMIC satellites, Liu et al. ^[73] also reported that the ionospheric VTEC significantly reduced in the afternoon of 4~6 days and the late evening of 3 days before Wenchuan earthquake, while extremely enhanced in the south and south-east of the epicenter along with their geomagnetic conjugate points in the afternoon of May 9. Further, considering the low Dst index during Wenchuan earthquake, they thought that all these ionospheric anomalies are the indications of seismo-generated electrom agnetic emissions. However, they didn't analyze the decrease of ionospheric VTEC on April 29, although the anomalous decrease is extremely striking, which was showed in their paper. Recently, Zhang et al. ^[74] found the dramatic decrease of the o⁺ density recorded by DEMETER satellite around the epicenter of W enchuan earthquake during the period of 10:00~13:00 LT on May 9 and thought that this abnormal phenomenon may be related to the earthquake. All these scholars mentioned above thought that the E × B drift was probably the explanation for seismo-ionospheric effects, while the E was one extra electric field which possibly originated from the increasing radon emanation, the photoelectrons provided by seismo-generated thermal radiations or the mobile positive holes activated in the crust bv microfractures during the dilatancy stage of earthquake preparation. Besides, some anomalies have been found in VLF transmi-tter signals and other parameters^{[75, 76]}.

Based on the above research results, it won widespread recognition that the ionospheric anoma-lies on May 9 were probably the precursors of Wenchuan earthquake. As for the ionospheric anomalies in the other periods, due to the differences of observations or methods for identifying anomalies, scholars have different views. In the present paper, relative to the future epicenter, the ionospheric anomalies probably related to the earthquake showed a long-distance drift from east to west (Fig. 5 and 8), which is hardly explained only using one extra electric field. Therefore, the possible coupling process of the earthquake-ionosphere needs a further investigation.

Up to now, earthquake forecast, especially short-term earthquake forecast, has not make a substantive breakthrough. Identifying different precursor signals is still a mainstream method for earthquake forecast. Compared with other known precursor signals of earthquake, seismo-ionospheric precursors have an advantage in determining the approximate time of earthquake. Recently, the statistical studies have proved that the anomalous variations in the ionosphere appeared within the time interval 0~5 days prior to the main shock^{[28, 32, 62]}. So. the ionospheric precursors are very valuable in the short-term earthquake forecast.

Out of question, the seismo-ionospheric effects increase our confidence in fulfilling earthquake forecast. However, only from the great earthquake in specific earthquake preparation zone (e. g., Wenchuan earthquake), we can not draw a universal conclusion that the pre-earthquake ionospheric VTEC anomalies related with strong earthquake must appear. Before the real physical mechanism of the seismo-ionospheric effects has been founded, we still need to make statistical analysis of more earthquakes to sum up the characters of ionospheric disturbances related to earthquakes.

This work has been supported by the Director Foundation of the Institute of Seismology, Chinese Earthquake Administration (grant No. IS200916012), Nation Key Technology R & D Program (grant No. 2OO8BAC35BO2) and National High Technology Research and Development Program of China (grant No. 2007AA12Z169). In addition, we thank the Crust M ovement Observation Network of China for providing the GPS observation data. the World Data Center for Geomagnetism, Kyoto for providing geomagnetic equatorial Dst index, and the NOAA National Geophysical Data Center for providing the data of solar F10. 7 flux and geomagnetic K_pindex.