Electromagnetic attenuation of eight earthquakes registered in Mexico using FFT-based spectrum and t-test statistical analysis for ULF Q-R ratios signals

Abstract

A method to improve the detection of seismo-magnetic signals is presented herein. Eight events registered for periods of 24 hours with seismic activity were analyzed and compared with non-seismic periods of the same duration. The distance between the earthquakes (EQs) and the ultra-low frequency detector is of  ρ = (1.8) 10^0.45M, where M is the magnitude of the EQ reported by the Seismological National Service of Mexico, in a period of three years. An improved fast Fourier transform analysis in the form of the ratio of the vertical magnetic field component to the horizontal one (Q = Bz/Bx) has been developed. There are important differences between the frequencies obtained during the days of seismic activity compared with those with no seismic activity.

1. Introduction

There has been a substantial accumulation of convincing evidence of electromagnetic phenomena associated to earthquakes (EQs) in the last decades (e.g. Hayakawa 1999, 2009; Hayakawa & Molchanov 2002; Molchanov & Hayakawa 2008). We can indicate a few examples of those seismo-electromagnetic effects which might be used for the short-term EQ prediction: (1) ULF (ultra-low frequency) electromagnetic emissions from the lithosphere (e.g. Fraser-Smith et al. 1990; Kopytenko et al. 1993; Hayakawa et al. 1996) and (2) seismo-ionospheric perturbations in the lower ionosphere (D/E regions) (e.g. Hayakawa et al. 2007) and in the upper ionosphere (F layer (or so) (e.g. Pulinets & Boyarnuch 2004; Liu et al. 2006). A recent monograph edited by Hayakawa (2009) includes several review papers on different seismogenic aspects containing recent progress, present knowledge and future work (Biagi 2009; Fraser-Smith 2009; Freund 2009; Molchanov 2009; Parrot 2009; Pulinets 2009).

It is well known that a great quantity of evidence on such seismo-electromagnetic phenomena has already been accumulated, though few papers have been recently published criticizing the presence of seismogenic ULF emissions for earlier events (Campbell 2009; Thomas et al. 2009) (although there is no guarantee that the conclusions of those papers are right). Further evidence in support of the presence (or absence) of seismo-electromagnetic phenomena is still highly required. In this paper, a new possibility is proposed; a study of spectra fluctuation (or modulation) of time series data of seismo-electromagnetic activity, in which attention is given to the modulations of different kinds of physical processes possibly taking place in the lithosphere, atmosphere and ionosphere/magnetosphere. The investigation of these modulations of the seismo-electromagnetic phenomena allows offering further support to the presence of those anomalies and also provides further understanding on their generation mechanisms. Furthermore, the slope of spectra fluctuation provides additional useful information on the self-organized criticality taking place in the lithosphere.

Evidence of electromagnetic emission at several frequencies associated with EQs has been accumulated over the years (e.g. Hayakawa & Molchanov 2002; Molchanov & Hayakawa 2008). The presence of strong interferences hampers the accumulation of the experimental data to test existing ideas on the mechano-electromagnetic conversion in the Earth's crust (Guglielmi et al. 2006). In addition, many authors perform signal processing in a wide range of frequencies for reported geomagnetic perturbations related to EQs. Such anomalies have been reported for several decades (Parrot & Johnston 1989; Johnston 1997; Kushwah et al. 2009). Seismogenic phenomena in the ULF range related to EQs with frequencies ranging from 0.005 to 1 Hz (Kopytenko et al. 1993) have been studied since Fraser-Smith et al. (1990) and a number of papers have been published on this topic since then (Smirnova et al. 2004; Serita et al. 2005; Kushwah et al. 2009). This frequency range is of practical importance for EQ activity, since these waves propagate easily on the Earth's surface where ULF sensors are installed (Ida et al. 2008).

One problem regarding these seismogenic ULF emissions is their detection as their signals are weak. Different methods of analysis have been developed: polarization analysis based on the ratio of the vertical magnetic field component to the horizontal (Hayakawa et al. 1996), fractal analysis (mono- and multi-) (Hayakawa et al. 1999; Gotoh et al. 2004; Kotsarenko et al. 2004, 2005, 2007; Smirnova et al. 2004; Ida et al. 2005; Ida & Hayakawa 2006), principal component analysis (Hattori et al. 2004; Kotsarenko et al. 2005) and singular value decomposition (Hattori et al. 2006) among others.

All the aforementioned methods are applied to improve both the detection of the ULF signals associated with seismogenic phenomena at different frequencies (Hayakawa et al. 2008), and the understanding of electromagnetic phenomena associated with tectonic and volcanic activity. In this paper, a simple analysis involving the ratio between vertical and horizontal components has been used. This analysis is completed with a statistical method applied to the FFT ULF data observed in Mexico during the years 2009, 2010 and 2011. Significant differences were found for the data obtained in periods of 24 hours with seismic activity and random signals of the same duration with no seismic activity for eight EQs recorded in Mexico near the observing station, considering a distance of ρ = (1.8) 10^0.45M (Chavez et al. 2011), where M is the magnitude of the EQ reported by the Seismological National Service of Mexico.

2. ULF geomagnetic data and EQs

ULF geomagnetic data was obtained at the Juriquilla (JU2) observatory, located at the geographic coordinates: –100.45° N, 20.70° E. In this site, three geomagnetic components (Bx; N–S component, By: E–W component and Bz: vertical component) are recorded using fluxgate sensors. The sampling frequency is 1 Hz. The analysis is focused on eight big EQs (with magnitude greater than 6.0) within a distance to the station of ρ as described previously. The events occurred between the years 2009 and 2011. The main characteristics of the EQs are presented in Table 1. The period of time analyzed for each event was of 24 hours. The distance between the EQs and the geomagnetic observatory is consistent with that used in previous studies (Chavez et al. 2010; Chavez et al. 2011). In those studies, the perturbations associated with the magnitude of the EQs could be observed.

Table 1. Earthquakes occurred in Mexico during the year 2008–2010. Magnitude greater than 6 was selected for this analysis (E), and magnitude less than 4.5 was selected for the random analysis (S). Year/month/day/hour/min are: the exact time of the EQ (local time), latitude and longitude: the geographic coordinates of the epicenter, magnitudes and depth: magnitude and depth of the EQ, distance: the distance between the epicenter and Juriquilla station and ρ: the radius of the EQ preparation zone estimated by the proposed equation. The EQ magnitudes are presented in bold (catalogue of National Seismological Service, Mexico).

	Event	Year	Month	Day	Hour	Min	Longitude	Latitude	Magnitude, M	Depth, km	Distance, km	ρ, km	Distance/ρ
Events selected M > 6.0 (E)	1	2008	2	12	5	50	−94.54	16.19	6.6	90	800	1679	0.48
	2	2008	4	14	22	3	−91.04	13.27	6.5	40	1262	1514	0.83
	3	2008	9	23	21	33	−105.16	17.16	6.4	42	630	1365	0.46
	4	2008	10	16	14	41	−92.5	13.87	6.6	23	1130	1679	0.67
	5	2009	8	3	13	0	−112.24	28.48	6.9	10	1473	2292	0.64
	6	2009	9	24	2	16	−107.43	17.72	6.2	21	805	1109	0.73
	7	2010	4	4	17	40	−115.36	32.54	7.2	10	1974	2245	0.87
	8	2010	6	30	2	22	−98.03	16.22	6.0	8	563	684	0.82
	1	2008	1	7	19	50	−99.07	15.73	3.7	13	567	83	6.81
	2	2008	3	3	10	26	−98.77	15.86	3.7	16	562	83	6.75
Days with events of M < 4.5 (S)	3	2008	16	5	11	37	−93.84	14.78	4.2	16	960	139	6.87
	4	2008	30	8	8	39	−98.03	16.28	3.8	27	553	92	5.99
	5	2009	28	7	19	46	−103.4	18.2	3.5	3.5	280	68	4.14
	6	2009	12	10	13	27	−99.95	16.86	3.6	15	424	75	5.65
	7	2010	29	3	16	41	−94.95	15.73	4.1	10	796	126	6.32
	8	2010	30	7	16	19	−95.58	15.71	4.3	39	752	155	4.85

3. FFT ratio Analysis

The polarization method developed by Hayakawa et al. (1996) is based on the measurement of the ratio of spectral power of the vertical magnetic field (Z) to the horizontal magnetic field (Bx and By). This ratio provides information on the variation of the geomagnetic field of ionospheric origin (or solar-terrestrial effect) or seismic origin (Ida et al. 2008). In this paper, the proposed ratios for the analysis are Q = Bz/Bx and R = Bz/By, representing the relation between the vertical and the horizontal components. These ratios are analogue to those previously used in other publications (e.g. Hayakawa et al. 2007).

The time series of the ULF data were processed to obtain the Q and R ratios of the magnetic field. Furthermore, the fast Fourier transform (FFT) was applied to the signals to obtain the power spectrum. Figure 1 shows three ratio signals corresponding to the FFT spectra for each one of the eight seismic events. This treatment was performed in order to analyze the components of the frequencies associated to seismogenic ULF emissions (Hayakawa et al. 2007). It can be observed that the FFT spectrum results in significant differences between the Bz/Bx components of the magnetic field.

Figure 1. FFT-based spectrum of the Bx, By and Bz components from ULF geomagnetic data. The blue lines show non-seismic spectra (XS, YS, ZS) while the red lines are seismic (XE, YE, ZE). FFT-based analysis does not show important differences between seismic and non-seismic data from the Bx, By and Bz components. To view this figure in colour, please see the online version of the journal.

The aforementioned behaviour is common to the eight analyzed seismic events. There are some important attenuation patterns during the days when the EQs occurred in all the analyzed frequency ranges (from 0.0 to 0.5 Hz). Those results can be associated with the distance, location, depth and magnitude of the EQs.

4. Observation and FFT analysis results

The proposed method allows observing some differences between the FFT analyses for the Bx, By and Bz components compared with the Q and R ratio signals. As aforementioned, each event's signal consisted in 24 hour data which was sampled at 1 Hz frequency using a fluxgate magnetometer.

The analysis procedure consisted in performing a 1024-point FFT for each event's signal from directions Bx, By, Bz, Q and R to obtain its spectrum followed by a scaling procedure where the spectrum was converted to dB scale to clarify results. The results for the Bx, By and Bz axes can be observed in Figure 1, where the blue lines correspond to superposed epochs with non-seismic signals while the red lines correspond to superposed epochs associated to seismic signals. As can be observed in Figure 1, there is no significant difference between non-seismic and seismic spectra in the Bx, By and Bz components because seismic and non-seismic spectral amplitude overlaps on most part of the signal. Furthermore, figure 2 shows the same procedure performed for Q and R ratios signals. An interesting attenuation in the seismic signal spectrum for the Q ratio (red line) is observed when compared to the non-seismic spectrum (blue line). On the other hand, the R spectrum shows a similar behaviour to the Bx, By and Bz components with no difference between seismic and non-seismic spectra.

Figure 2. FFT-based spectrum of the Q = Bz/Bx and R = Bz/By ratio signals from ULF geomagnetic data. The blue line indicates epochs with no seismic activity (QS and RS) while the red line represents the seismic activity (QE and RE). FFT-based analysis does not show important differences on R ratio but an interesting attenuation in seismic (red line) data from Q ratio spectra can be observed. To view this figure in colour, please see the online version of the journal.

Accordingly, a post processing analysis based on the minimum square procedure was performed for a certain range of the Bx, By, Bz and R spectra to obtain a quantitative measurement of the peak amplitude of FFT-based analysis. The minimum square range was established between 0.0303 and 0.616 Hz and was determined experimentally. Figures 1 and 2 show the selected frequency range for the minimum square procedure. We know that the geomagnetic field is affected by solar activity as shown by Wen et al. (2012). However, as will be seen in this paper, the FFT analysis results in geomagnetic relations on the Q ratio only, which indicates that the solar activity can be discriminated for this purpose.

5. Square summation-based quantitative analysis

In order to quantify the FFT-based results, a well-known square summation procedure was applied to the resulting FFT spectra by performing a summation of the square of each point in the spectrum in the range b–a as indicated by equation (1)(1) where w is the frequency of the FFT domain, B(w) is the FFT spectrum of the magnetic field B, a is the lower limit, b is the upper limit and Square sum (Sqsum) is the result of the summation procedure. Consequently, a scalar value is obtained from each FFT signal that represents a quantitative measurement of its spectrum amplitude in the aforementioned b–a range. (1)

6. Statistical analysis results

Furthermore, a statistical analysis was performed using t-test and Mann–Whitney algorithms to compare the statistical properties in each axis under calm conditions and seismic activity. It means that the analysis focused on the existence of statistical differences between calm and seismic periods for each axis from the data extracted during FFT analysis and quantified by equation (1)(1) . The statistical analyses groups were organized by choosing an N number of elements per group of eight (N = 8). Therefore, eight signals from non-seismic readings (Bx = XS, By = YS, Bz = ZS components and RS and QS ratios) and eight during seismic conditions (Bx = XE, By = YE, Bz = ZE components, RE and QE ratios) from each axis were used for this analysis as it can be observed in Table 2.

Table 2. Proposed statistical tests.

Test	Statistical tool for comparison
X axis calm (XS) vs. X axis seismic (XE)	t-Test
Y axis calm (YS) vs. Y axis seismic (YE)	t-Test
Z axis calm (ZS) vs. Z axis seismic (ZE)	Mann–Whitney rank sum test
Q ratio calm (QS) vs. Q ratio seismic (QE)	Mann–Whitney rank sum test
R ratio calm (RS) vs. R ratio seismic (RE)	Mann–Whitney rank sum test

Mean values and standard deviations of each group were obtained using Sigma Plot 11 software. Figure 3 shows the results of the aforementioned statistical analysis for each group including the mean value and the +/− standard deviation bar. This plot shows that, when comparing each component or their ratio during seismic and non-seismic conditions, an overlap occurs for axis Bx (XE–XS), By (YE–YS), Bz (ZE–ZS) components Q (QE–QS) and R (RE–RS) ratio. However, there is no overlap between calm and seismic periods for the Q ratio. These results indicate that the Q ratio may be a means to identify seismic events.

Figure 3. Statistical analysis for Bx (XE–XS), By (YE–YS), Bz (ZE–ZS), Q = Bz/Bx (QS–QE) and R = Bz/By (RS–RE). Mean +/− three sigma bars for each statistical analysis group where only QS and QE groups have a good statistical correlation.

In order to identify how different are calm and seismic periods in each axis, statistical tests were performed for the pairs of values shown in Table 2.

Different statistical tools can be used to compare different kinds of experimental data. However, this project requires statistical comparisons between seismic and calm periods for each axis. The t-test has shown to be appropriate to statistically compare two groups (Rice 2006). In that sense, this test was selected as the statistical analysis tool to compare the Bx and By components of the magnetic field. However, the t-test is not recommended for non-parametric or highly scattered data, and therefore, it results inadequate for the Bz component and Q and R ratios as they are of non-parametric nature. Because of this, the Mann–Whitney rank sum test was selected to statistically compare the calm and seismic periods for the Bz component and the Q and R ratios (Fay & Proschan 2010).

On the other hand, the p-value represents a probability of obtaining the same sample results when the null hypothesis is true. Where null hypothesis means that there is no relationship between two phenomena. In order to perform a comparison tests, a constant value for the p-value must be chosen. When the p-value resulting from a comparison test between two groups is less than 0.05, it is considered that there is statistical difference between these groups. However, if the p-value obtained from a comparison test exceeds the 0.05 threshold, it is considered that there is no statistical difference between these groups.

The statistical analysis performed in this study shows p-values during calm and seismic periods exceeding the 0.05 threshold for the Bx (XE–XS), By (YE–YS) and Bz (ZE–ZS) components and the R (RE–RS) ratio (0.934, 0.466, 0.878 and 0.442, respectively) which means that it is impossible to distinguish between seismic and calm signals for these components, as it can be observed in Table 3. However, the Q ratio showed a p = 0.001 which indicates that it can be a promising index to discriminate between seismic QE from non-seismic QS events due to its highly significant statistical difference gathered by the Mann–Whitney rank sum test.

Table 3. Statistical analysis for paired groups Bx (XS–XE), By (YS–YE), Bz (ZS–ZE), Q = Bz/Bx (QS–QE) and R = Bz/By (RS–RE).

Group name	N	Mean	Standard deviation	Comparison test	P value	Significant statistical difference
XS	8	203403.825	35129.963	XS–XE t-test	0.934	No
XE	8	201769.508	42371.540
YS	8	200739.302	12748.461	YS–YE t-test	0.466	No
YE	8	206734.216	18694.726
ZS	8	201770	42372	ZS–ZE Mann–Whitney rank sum test	0.878	No
ZE	8	187260	8171.5
QS	8	173570	48958.0	QS–QE Mann–Whitney rank sum test	0.001	Yes
QE	8	10885	3459.0
RS	8	187310	171610	RS–RE Mann–Whitney rank sum test	0.442	No
RE	8	258590	50834

7. Conclusions

There are seismo-magnetic signals related to EQs disturbing the earth's electromagnetic field. The main problem of the signal detection is the discrimination of the background noise. For that reason, the use of the FFT-based analysis on Q and R ratio signals is proposed. Additionally, an improved statistical ratio analysis method has been applied to the ULF geomagnetic data at Juriquilla station in Mexico. The three geomagnetic field components recorded behave in different ways, for that reason a general spectrum analysis is performed for the Bx, By and Bz components. Also, the spectrum analysis methodology is proposed for the Q and R ratio signals with interesting results. It is well known that the spectra in the Bz (XE–XS), By(YE–YS) and Bz (ZE–ZS) axis do not provide significant information to discriminate between seismic and non-seismic events. The same happens with the R (RE–RS) ratio signal spectrum which does not provide with any interesting information. However, the Q (QE–QS) ratio signal showed an interesting attenuation in its spectrum when seismic events occur. By using the well-known minimum square methodology to quantify the data dispersion in the resulting spectra, a significant difference between the Q ratio for seismic and non-seismic data was confirmed. The minimum square methodology used for the Q ratio spectrum shows to be enough robust to discriminate seismic from non-seismic electromagnetic signals. The p-value obtained from the Mann–Whitney rank sum test for eight different events is far below the threshold that indicates statistical difference between two groups (p = 0.001). This methodology can be useful in future projects to accurately discriminate seismic from calm periods.

Acknowledgements

The authors are grateful to CONACyT for their support and collaborations in this research under the SNI assistance under the project number 54600. The authors would also like to thank PhD Hilda Romero Zepeda for editing the English version and to LAV Rodrigo Beltrán Soria for editing the figures of this paper. The authors thank to the anonymous referees for their observations and recommendations to improve this paper.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

CONACyT: SNI, Project Number 54600 [grant number 209837].