A statistical study on the effect of earthquakes on the ionosphere, based on the subionospheric LF propagation data in Japan

A superimposed epoch analysis has been undertaken, in order to find the correlation of the ionospheric perturbations with seismic activity. We take the wave path from the Japanese LF transmitter (frequency=40 kHz) and an observing station of Kochi (wave path length of 770 km), and a much longer period (of five years) than before, is considered. This subionospheric LF propagation can be called “an integrated measurement” in the sense that any earthquakes in the LF sensitive area just around the great-circle path can influence the observed LF signals, so that we define the “effective magnitude” (Meff) by integrating the total energy from different earthquakes in the sensitive area on a current day and by converting it back into magnitude. A superimposed epoch analysis for the effective magnitude greater than 6.0 has yielded that the ionosphere is definitely disturbed in terms of both amplitude and dispersion, and that these perturbations tend to take place prior to an earthquake. The statistical z-test has also been performed, which has indicated that the amplitude is definitely depleted 2–6 days before the earthquake day and also that the dispersion is very much enhanced during the same period. This statistical study has given strong support to the existence of seismo-ionospheric perturbations for high seismic activity.


Introduction
It is recently recognized that electromagnetic phenomena would take place prior to an earthquake (e.g.Hayakawa and Fujinawa, 1994;Hayakawa, 1999;Hayakawa and Molchanov, 2002), and a further seismic effect appears not only in the lithosphere, but also in the atmosphere and iono-Correspondence to: M. Hayakawa (hayakawa@whistler.ee.uec.ac.jp) sphere (Hayakawa, 2004).Subionospheric VLF/LF propagation is widely used in recent years to study such seismoionospheric perturbations.Two possible methods have been proposed for such a study of subionospheric VLF/LF propagation data.The first is based on the analysis of nighttime amplitude and/or phase anomalies (Gokhberg et al., 1989;Gufeld et al., 1992), and further developed as the "fluctuation method" (Shvets et al., 2002(Shvets et al., , 2004a, b;, b;Hayakawa et al., 2002Hayakawa et al., , 2004a, b;, b;Horie et al., 2006).The second method is called the "terminator time" method, which is based on the determination of the characteristic times of minimums in the amplitude/phase diurnal variations during sunrise and sunset (Hayakawa et al., 1996;Maekawa and Hayakawa, 2006).The first report on the use of subionospheric VLF propagation was published by Gokhberg et al. (1989), who studied the bay-like phase anomalies of VLF Omega signals for long paths and suggested the relationship of such VLF propagation anomalies to the earthquakes.Hayakawa et al. (1996) then presented the first convincing evidence on the ionospheric perturbations for the famous Kobe earthquake (on 17 January 1995) by means of the VLF data of reception at Inubo of the Japanese Omega signals transmitted from Tsushima, Kyushu, and they found anomalous shifts in the terminator time from a few days before the earthquake until the earthquake date.Hayakawa et al. (1996) and Molchanov et al. (1998) interpreted this terminator time shift in terms of the lowering of the ionosphere by a few kilometers, and this idea was further supported by the following full-wave computation by Soloviev and Hayakawa (2002).Using the same terminator time method, Molchanov and Hayakawa (1998) further investigated the data during thirteen years for the same propagation path between Tsushima and Inubo, and found that the ionospheric perturbations appear for a majority (∼80%) of large earthquakes whose magnitude is greater than 6.0, whose epicenter is very close to the great-circle paths and also whose depth is shallow.Later many papers have been published on the event studies for recent large earthquakes including Tokachi-oki (August 2003) (Shvets et al., 2004b(Shvets et al., ), 2004 Mid Niigata Prefecture (Hayakawa et al., 2006), 2004Indonesia Sumatra (Horie et al., 2006) earthquakes.These kinds of event studies would be of essential importance for the study of seismoionospheric perturbations.However, in addition to these event studies, it is highly required to undertake any statistical study on the correlation between ionospheric disturbances and earthquakes based on abundant data sources.There have been very few reports on the statistical correlation between the ionospherc perturbations and earthquakes (Shvets et al., 2002(Shvets et al., , 2004a;;Rozhnoi et al., 2004).Shvets et al. (2004a) have examined a very short-period of data (March-August, 1997) for two paths (one is the Tsushima-Chofu and another, NWC (Australia)-Chofu) and found that wave-like anomalies in the VLF Omega signal with periods of a few hours (as indicative of the importance of atmospheric gravity wave as suggested by Molchanov et al., 2001) were observed 1-3 days before or on the day of moderately strong earthquakes with magnitudes 5-6.1.Then, Rozhnoi et al. (2004) have extensively studied 2 years data of the subionospheric LF signal along the path Japan (call sign, JJY)-Kamchatka (dis-tance=2300 km), and have found from the statistical study that the LF signal effect is observed only for earthquakes with a magnitude at least greater than 5.5.This paper will be devoted to such a statistical study on the correlation between ionospheric disturbances and seismic activity.A few important distinctions from the previous works by Shvets et al. (2002Shvets et al. ( , 2004a) ) and Rozhnoi et al. (2004) are described.The first point is the use of a much longer period of VLF/LF data (five years long).The second point is that we pay attention to the physical parameters of VLF/LF propagation data in this paper; (1) amplitude (or trend) and (2) dispersion (in amplitude) (or fluctuation).In the previous work by Rozhnoi et al. (2004) they have studied the percentage occurrence of anomalous days, in which an anomalous day is defined as one day during which the difference of amplitude (and/or phase) from the monthly average exceeds one standard deviation (σ ).A superimposed epoch analysis is undertaken in this paper and also the statistical test is performed to estimate the significance of seismo-ionospheric perturbations.
In this paper we pay particular attention to the earthquakes occurring in and around Japan, so that we take a wave path from the Japanese LF transmitter, JJY (40 kHz) (geographic coordinates; 36 • 18 N, 139 • 85 E) and a receiving station of Kochi (33 • 33 N, 133 • 32 E). Figure 1 illustrates the relative location of the LF transmitter, JJY and our receiving station, Kochi, and the distance between the transmitter and receiver is 770 km.

Subionospheric LF data and earthquakes for analysis
The subionospheric LF data for this propagation path is taken over 6 years from June 1999 to June 2005, but we excluded one year of 2004 (January to December, 2004) because of the following reason.As you may know, there was an extremely large earthquake named 2004 Mid Niigata prefecture earthquake which happened on 23 October, with a magnitude=6.8and a depth=10 km (Hayakawa et al., 2006), and the effect of the main shock and also the large aftershocks were so large and so frequent that it could disturb our following statistical result greatly.As a result, we have excluded this year (2004) from our analysis.
We have to define the criterion for choosing the earthquakes.The sensitive area for the wave path, the JJY transmitter to the Kochi receiving station, is defined as follows.First, we adopt the circles with a radius of 200 km just around  the transmitter and receiver, and then the sensitive area is defined by connecting the outer edges of these two circles.This is indicated in Fig. 1, as well.All of the 92 earthquakes with a magnitude (conventional magnitude (M) by Japan Meteorological Agency) greater than 5.0 are plotted in Fig. 1, but the earthquake depth is chosen to be smaller than 100 km (taking into account our previous result that shallow earthquakes can have an effect on the ionosphere by Molchanov and Hayakawa, 1998).We have normally been using the fifth Fresnel zone as the VLF/LF sensitive area (Molchanov and Hayakawa, 1998;Rozhnoi et al., 2004), but we have found that the area just around the transmitter and receiver is also sensitive to VLF perturbation (e.g.Ohta et al., 2000), taking into account the possible size of the seismo-ionospheric perturbation.In this sense the sensitive area we choose here seems to be very reasonable because the width of the sensitive area is very close to the 10th Fresnel zone.
In the following statistical analysis, we undertake the socalled superimposed epoch analysis (Taylor et al., 1994), in order to increase the S/N ratio.Here we define the earthquake magnitude in the following different way.Since we treat the data in the unit of one (a) day (we use UT (rather than LT) to count a day because we stay on the same day even when we pass midnight when we use UT), we first estimate the total energy released from several earthquakes with magnitude M in one day within the sensitive area for the LF wave path, as shown in Fig. 1, by integrating the energy released by a few earthquakes (down to the conventional magnitude M=2.0) and by inverting this into an effective magnitude (Meff) for this particular day.This Meff is much more important than the conventional magnitude for each earthquake, because the LF propagation anomaly on one day is the effect integrated over several earthquakes taking place within the sensitive area on that day. Figure 2

VLF data analysis and superimposed epoch analysis
Diurnal variations of the amplitude and phase of subionospheric VLF/LF signal are known to change significantly from month to month and from day to day.Therefore, following our previous works (Shvets et al., 2002(Shvets et al., , 2004a, b;, b;Rozhnoi et al., 2004;Hayakawa et al., 2006;Horie et al., 2006), we use, for our analysis, a residual signal of amplitude dA as the difference between the observed signal intensity (amplitude) and the average of several days preceding or following the current day: where A(t) is the amplitude for a current day and <A(t)> is the corresponding average for ±15 days (15 days before, 15 days after the earthquake and earthquake day) in this paper.In the paper by Rozhnoi et al. (2004), they have defined an anomalous day when dA exceeds the corresponding standard deviation.Figure 3 shows the occurrence histogram of the recurrent time interval (in units of a day) for a Meff greater than 5.5.The result in Fig. 3 gives us some support to our adoption of ±15 days (total, one month) when estimating dA(t).When we try to perform the superimposed epoch analysis, the event should be well isolated to increase the S/N ratio.In our analysis we have studied the nighttime variation (in the UT range from UT=10 h to 20 ) (or LT 19 h to 05 h).Then, we use two physical parameters: average amplitude (we call it "amplitude") (or trend) and amplitude dispersion (we call it "dispersion) (or fluctuation).We estimate the average amplitude for each day (in terms of UT) by using the Then, we are ready to undertake a superimposed epoch analysis (e.g.Taylor et al., 1994;Oike and Yamada, 1994).For the study of the correlation between ionospheric perturbations in terms of the two parameters (amplitude and dispersion) and seismicity, we choose two characteristics periods; seismically active periods with a Meff greater than 5.5 and greater than 6.0.The number of events with a Meff 5.5 is 19, and that with a Meff 6.0 is 4. The raw (stacked) data of superimposed epoch analysis is illustrated in Fig. 4. Figure 4a refers to the result on amplitude, while Fig. 4b refers to that on dispersion.The upper panels refer to a Meff 6.0, and the lower panels to a Meff 5.5.
Let us look at Figs. 4a for amplitude and 4b for dispersion.The abscissa indicates the lag day; minus means the day before the relevant (earthquake) day, and day 0 indicates the earthquake day.The ordinate indicates the amplitude and dispersion values (in dB) after stacking all the data.A seismically active period is characterized by a Meff either greater than 5.5 or 6.0.Rozhnoi et al. ( 2004) have already studied extensively the effects of solar flares, geomagnetic disturbances, in addition to the seismic effect, and they have shown that the seismic effect is definitely seen for a conventional earthquake magnitude (M) greater than 5.5.Let us look at the result for the Meff greater than 5.5 (as shown in red).The amplitude seems to be depleted before the day zero and also depleted 5 days after the zero day.At any rate, the amplitude seems to be affected before and after the zero day.Also the dispersion is likely to be rather enhanced before the earthquake day.However, it is rather difficult for us to judge whether the variations in amplitude and dispersion have shown so significant changes or not, before going into the statistical test.Then, we go to the results for a very seismically active period with a Meff greater than 6.0, which is given in the top panel in Fig. 4a in blue.It is clearly seen from this figure after stacking over 4 events that the amplitude is showing a significant depletion over a few days (6 to 2 days) before the earthquake as a precursory effect of earthquakes.Further, Fig. 4b also suggests that the dispersion is absolutely extremely enhanced over one week before the earthquake.
Finally, we would like to undertake the statistical test for the results in Fig. 4. When we perform the Fisher's ztransformation to the data in Fig. 4 (amplitude and dispersion, respectively), the z value is known to follow approximately the normal distribution of N (0, 1) with zero average and dispersion of unity (Takeuchi et al., 1982;Bickel and Doksum, 2001).
Figures 5a and b represent the corresponding statistical ztest result for the corresponding Figs.4a and b, respectively.The 2σ (σ : standard deviation over the whole period of five years) line is indicated as the statistical criterion.First of all, we look at Fig. 5a.It is clear that the blue line for the Meff greater than 6.0 exceeds the 2σ line a few days before the earthquake.This suggests that the ionospheric perturba-  tion in terms of amplitude (trend) shows a statistically significant precursory behavior (3 to 5 days before the earthquake).Also, when the Meff becomes a little bit smaller (Meff 5.5), the corresponding figure (in red) is found to be very close to the 2σ line before the earthquake, but is exceeding the 2σ line after the earthquake.This is indicative again of the statistically significant effect of earthquakes on the ionosphere just around the day of high seismicity.
Next we go to Fig. 5b for the dispersion.The enhancement of dispersion (fluctuation) is clearly visible for extremely high seismic activity (Meff 6.0), i.e. the dispersion is found to exceed the 2σ line 2-6 days well before the earthquake day, while, in the case of a smaller Meff (Meff 5.5) the dispersion is seen to be considerably enhanced several to a few days before and just after the earthquake.But the dispersion value does not exceed the 2σ criterion.
Finally, we comment on the corresponding result for M 5.0 (further below Meff=5.5 by 0.5).We have found that the variations in amplitude and dispersion, as in Fig. 5 well inside the ±2σ line for a Meff 5.0, and, together with our previous findings, we say that the seismic effect can only be seen, at least, for a Meff 5.5.

Conclusions
A superimposed epoch analysis has been applied to the subionospheric LF data, in order to find the possible effect of earthquakes on the ionosphere.We define the Meff of the earthquakes for a certain day in such a way that the total energy released by all earthquakes in one day is re-converted to the magnitude formula.The use of this Meff is effective because any earthquakes in the LF sensitive area can affect the ionosphere.Two categories for this Meff are considered; seismically active periods with a Meff 5.5 and a Meff 6.0 We can conclude, from the present statistical analysis, the following main important results: 1.When the Meff is greater than 6.0, the superimposed epoch analysis has yielded that the ionosphere is definitely disturbed in terms of both amplitude and dispersion.The amplitude is depleted about one week to a few days before the earthquake, and also the dispersion is very much enhanced during the same period before the earthquake.
2. By using the z-transformation to the previous amplitude and dispersion data, and by performing the statistical ztest, it is found that such changes in both amplitude and dispersion during the same period before the earthquake exceed the 2σ (σ : standard deviation) criterion, indicating the statistical significance.
3. When the Meff becomes a little bit smaller, but greater than 5.5, the statistical z-test indicates that both amplitude and dispersion exhibit the similar tendency as for the case of a Meff 6.0, but not so significant.
We compare our present statistical result with previous ones (Shvets et al., 2002;Rozhnoi et al., 2004).Shvets et al. (2002) have performed different kinds of correlations, but unfortunately the analysis period was too short, of the order of a few months.Rozhnoi et al. (2004) have used the data of two years, and they have treated an "anomalous" day in the definition that the difference (or residual) (in amplitude or phase) (dA or dP) exceeds the corresponding 1σ .
They have studied the percentage occurrence of such anomalous days for different conventional earthquake magnitudes (Ms).After examing different effects (solar flares, geomagnetic storms, etc.), they have succeeded in detecting the seismic effect in subionospheric VLF/LF propagation only when the earthquake magnitude exceeds 5.5.In our paper, we do not pay attention to the percentage occurrence of anomalous days, as studied by Rozhnoi et al. (2004), but we pay attention to two physical parameters of subionospheric LF prop-agation: (1) amplitude (trend) and ( 2) dispersion (or fluctuation).A superimposed analysis, together with the statistical test, has yielded that both these two physical quantities, amplitude and dispersion, are statistically significantly disturbed for the Meff greater than 6.0, i.e. the amplitude is decreased by ∼3 dB one week to a few days before the earthquake, and at the same time during the same period the dispersion is significantly enhanced, as well.The similar tendency is confirmed for the Meff greater than 5.5 (but not so significant as compared to the case of Meff 6.0, but statistically significant enough).Our result seems to have confirmed and supported our previous result by Rozhnoi et al. (2004) by using the much longer-period data.The present statistical study has given a strong validation of the use of the nighttime fluctuation method to determine seismo-ionospheric perturbations (Hayakawa et al., 2006;Horie et al., 2006) Based on our present statistical result, we are ready to go into the details on the generation mechanism of ionospheric perturbations caused by seismic activity.A few hypotheses on the lithosphere-atmosphere-ionosphere coupling have already been proposed (Hayakawa, 2001(Hayakawa, , 2004;;Molchanov et al., 2001;Hayakawa et al., 2004); (1) chemical channel, (2) acoustic channel and (3) electromagnetic channel (Molchanov et al., 1993).It seems likely that the first two channels are promising.As for the chemical channel, we expect modification of the electric fields and currents in the atmosphere due to the change in atmospheric conductivity (e.g.due to the radon emanation) over an earthquake zone, and corresponding effects in the ionosphere (Grimalsky et al., 2003;Hayakawa, 2004;Pulinets and Boyarchuk, 2004;Sorokin et al., 2005).Then, in the acoustic channel we assume the transfer of the disturbances from a seismic source to the atmosphere and ionosphere by means of acoustic and internal gravity waves (Molchanov and Hayakawa, 1998;Molchanov et al., 2001;Miyaki et al., 2002;Hayakawa et al., 2002;Shvets et al., 2004).Furthermore, it is interesting to investigate the correlation and relationship of this VLF propagation anomaly in the lower ionosphere with the upper ionospheric condition (e.g.Liu et al., 2000).Further studies to elucidate the mechanism are being carried out.

Fig. 1 Fig. 1 .
Fig.1 Fig. 1.Relative location of the LF transmitter, JJY in Fukushima and an observing station, Kochi.The sensitive area for this LF propagation path is also indicated; the circles with radius of 200 km around the transmitter and receiver and by connecting the outer edges of these two circles.Also 92 earthquakes with a conventional magnitude (M) greater than 5.0 are plotted, which took place within the sensitive area.

Fig. 2 .
Fig. 2. Occurrence histogram of effective magnitude (Meff) during the whole five years.The ordinate is the number of days with relevant effective magnitude.

Fig. 3 .
Fig. 3. Occurrence histogram of the recurrence interval (in the unit of a day) for the effective magnitude, Meff 5.5.

Fig. 5 .
Fig. 5. Statistical test result for the amplitude (a) and dispersion (b).The day on the abscissa has the same meaning in Fig. 4. The important 2σ (σ : standard deviation) lines are plotted for the statistical test.