Electromagnetic Precursors to the 2016 Kumamoto Earthquakes

The Kumamoto area of Kyusyu Island was attacked by a series of large earthquakes (EQs) in April, 2016. The first two foreshocks had the magnitudes of 6.5 and 6.4, and about 1 day later there was the main shock on 15 April (UT) with magnitude 7.3. These are fault-type EQs, and so we would expect a variety of electromagnetic precursors to these EQs because we had detected different phenomena for the 1995 Kobe EQ, same fault-type EQ. As for the lithospheric effect, the ULF data at Kanoya observatory (about 150 km from the EQ epicenters) are used, but the simple statistical analysis could not provide us with any clear evidence of ULF radiation from the lithosphere. However, our conventional analyses indicated clear signatures in the atmosphere as ULF/ELF impulsive emissions and also in the ionosphere as observed by means of VLF propagation anomalies and ULF depression. ULF/ELF radiation appeared on 8-11 April (in UT) (maximum on 10 and 11 April (UT)), while ULF depression took place on 8 and 10 April (in UT), so that both atmospheric radiation and ionospheric perturbation took place nearly during the same time period.


Introduction
Even after the disastrous 2011 Tohoku earthquake (EQ) [1], the western part of How to cite this paper: Schekotov, A., Izutsu, J., Asano, T., Potirakis Japan remained seismically rather quiet though some enhancement in seismicity was observed such as an M (magnitude) 6 class EQ in 2013 in the same place as the 1995 Kobe EQ [2].Though the medium-term forecast suggested that the occurrence probability of M7 class EQs for the coming 30 years in the Kumamoto region of Kyushu Island was about 1% or so, there happened a series of serious EQs in the region of Kumamoto.Unlike the subduction-type 2011 Tohoku EQ, the EQs in Kumamoto were of the fault-type, just like the 1995 Kobe EQ.As a general tendency, it seems that the fault-type EQs are rich in electromagnetic precursors such as the case of the 1995 Kobe EQ [3], while electromagnetic phenomena are not as abundant for the sea EQs such as the 2011 March Tohoku EQ [4].
In this paper, we pay attention to the electromagnetic phenomena in possible association with this Kumamoto EQ series.It is needless to say that statistical studies are highly required on the correlation of any seismogenic phenomenon with EQs on the basis of long-term observation (e.g., [5] on the upper ionospheric perturbations, [6] on the VLF lower ionospheric perturbations, [7] on the lithospheric ULF emissions, [8] on the ULF/ELF radiation and ULF depression).

EQs Treated
There happened three successive EQs in the Kumamoto area with an extremely low probability of medium-term forecast: 1) M = 6.

Electromagnetic Phenomena in Possible Association with the Kumamoto EQ Event
There may be a few interesting electromagnetic phenomena in possible association Open Journal of Earthquake Research with the Kumamoto EQ event [23] [24]: 1) ULF (ultra low frequency, f < 3 Hz) magnetic field variations (both lithospheric radiation and ULF depression) observed at Kanoya (KNY) station in Kyushu island, 2) ULF/ELF radiation (in the frequency range of 0.1 -24 Hz) as observed with the use of Chubu University ELF network, and 3) monitoring of ionospheric perturbations with the use of our multi-stationed VLF network all over Japan.JMA has been conducting the regular observations of DC/ULF magnetic field variations at several standard observatories in Japan and KNY is one of such stations, which is closest to the EQ epicenters.The ULF/ELF observation network of Chubu University is composed of three observing stations of Nakatsugawa (NAK), Shinojima (SHI) and Izu (IZU) (see Figure 1) [21].The last, our own VLF/LF (very low frequency/ low frequency) network has been in operation for the last five years, and we currently have a network consisted of eight stations in Japan [22].

ULF Magnetic Field Variations (Lithospheric Radiation and ULF Depression)
There is a ULF observatory at KNY belonging to JMA whose geographic coordinates are 31˚25'27''N, 130˚52'48''E (Figure 1).Three magnetic field components (H (horizontal), D (declination), and Z (vertical)) of ULF magnetic fields have been measured at KNY by a fluxgate magnetometer with a sampling frequency of 1 Hz.
Figure 2 illustrates the temporal evolutions of ULF magnetic field variations to both the geomagnetic activity (K p ) and EQ index with the use of K LS [24], the 2nd and 3rd panels correspond to the horizontal (H) and vertical (Z) magnetic components, 4th to the polarization (i.e., the ratio of vertical to horizontal component) P Z/H [25], while the 5th panel refers to the inverse of horizontal component (Dh) and the bottom panel refers to the relative depression, δDep.The vertical line is the time of the main shock on 15 April (UT).We have to mention that the results on the 2nd to 4th panels refer to the lithospheric radiation, while the latter two refer to the ULF depression as an indicator of the lower ionospheric perturbation.First, we pay attention to the lithospheric radiation (direct radiation from the lithosphere) as seen from the 2nd to 4th panels in Figure 2. The horizontal component (F h ) is found to exhibit an increase around 10 March, about 1 month before the EQ.A comparison between the K p index variation and F h and F z variations, may indicate that there may be some periods of enhanced F h and F z when the K p index showed an increase, however there is a few periods of enhanced F z only.In this sense, it is quite difficult to say that there appeared clear seismogenic ULF radiation.The most significant parameter of polarization (P z/h ) (as proposed by [25]) exhibited a significant peak in the end of February and also we can say that some peaks in P z/h exceeding unity tend to appear mainly before the Open Journal of Earthquake Research EQ (during about half a month before the EQ).But, we cannot come to any definite conclusion on the presence of precursory ULF radiation with this conventional statistical analysis.So that, we have to proceed to the non-conventional analysis such as the application of natural time to these ULF data (e.g., [26] [27] [28]), and those national time results will be published elsewhere [29].
On the other hand, even the conventional statistical analysis on D h and δDep has yielded very clear peaks about one week before the EQ.Especially a peak on 8 April is very notable, with δDep exceeding greatly the value of 10.Nobody can oppose to this presence.The lead time of about one week is found to be very consistent with our previous statistical results [30].

ELF Electromagnetic Radiation
The observatories of ULF/ELF electromagnetic radiation are indicated in Figure 1 as NAK (Nakatsugawa, 35.42˚N, 137.55˚E),SHI (Shinojima, 34.67˚N, 137.01˚E) and IZU (Izu, 34.64˚N, 137.01˚E).Three magnetic field components (H, D and Z components) are observed at each observatory by means of three orthogonal magnetometers in the frequency range of 0.1 -24 Hz (with sampling of 100 Hz) [21].Schekotov et al. (2007)  [31] have suggested in their first paper that seismoatmospheric radiation in the ULF/ELF frequency band seems to provide us with a possibility of predicting an EQ; not only predicting the occurrence time of a forthcoming EQ but also indicating the source location (or its epicentral position).Some later works [24] [32] have confirmed that the direction of radiation source is coincident approximately with the position of the epicenter of a future EQ.
Though the details of signal processing method have been summarized in [30], we repeat only the essential points on how to detect seismo-atmospheric ULF/ELF radiation and then how to determine the direction of arrival of the radiation.

1) Direction finding (arrival direction)
The source direction of seismo-atmospheric ULF/ELF radiation is determined as being perpendicular to the main axis of polarization ellipse.The angle between the main axis of the polarization axis and the d (EW)-component axis is denoted by θ, and its tangent is given by the following equation.For the required accuracy of direction finding, we have chosen the signal Open Journal of Earthquake Research whose amplitude is 5 times the average amplitude in order to have sufficient signal to noise ratio.

2) Detection of the radiation
In [31], a new following parameter ∆S is proposed in order to detect seismoatmospheric radiation, which is given by, ( ) The numerator refers to the ratio of two horizontal spectral components P hh (NS component of magnetic field) and P dd (EW component).And the denominator is the root mean square (rms) of the deviation of the signal ellipticity.The expression of β is given then by the following equation. ( Here Im means imaginary part.Because Shekotov et al. ( 2007) [31] have compared different wave parameters with respect to their association with EQs, and have found an enhancement in the spectral ratio of P hh /P dd and a reduction in the polarization ellipticity before an EQ.So the parameter introduced by Equation ( 2) is proved to be most sensitive and reproducible to seismic shock.The ellipticity or the ratio of minor axis to major axis is defined by tanβ.The sense of polarization is determined by the sign of β; the polarization is right-hand (RH) when β > 0, and β < 0 means the left-hand (LH) polarization.The most important linear polarization is expressed by β = 0. We have calculated the field component power spectral densities, P hh , P dd and their cross-power spectral densities P hd , P dh were calculated by using Fourier transforms with frequency resolution of about 0.1 Hz.Spectral components in a frequency range from 0.1 to 24 Hz were used, and they were averaged over one-Hz intervals such as 0.1 -1, 1.1 -2, ••••••, 23.1 -24 Hz, so that we have 24 spectral components.
Figure 3 illustrates the result on ∆S at a particular station of SHI during a period of 06 April to 15 April 2016 (EQ day), because we observed a lot of interference at other two stations.In this analysis we have used only the local nighttime data in JST from 0.5 to 5.0 h where we expect the minimum local noise.The bottom panel illustrates the temporal evolution of ∆S, which indicates that there is observed a significant broad peak from 8 to 11 April (in UT) (maximum on 10 and 11 April) (so that 9 to 12 April (in JST)), which is likely to be a precursor to the EQ event.This lead time of about one week is in good agreement with our previous papers [21] [30] [31].
Figure 4 illustrates the temporal evolutions of the azimuthal distribution of ULF/ELF impulsive radiation.Unfortunately, the noise conditions at other two stations of NAK and IZU are not so good enough to perform the direction finding,  so that the results only from SHI are plotted there.The azimuthal directions of ULF/ELF radiation on the active periods of 8 to 11 April are found to be approximately directed to the Kyushu area, being consistent with the speculation that those ULF/ELF emissions are likely to be generated close to the EQ epicenters.

1) VLF monitoring
Based on our network observation of subionospheric VLF/LF signals, Hayakawa and Asano (2016) [33] have published some preliminary analysis and have shown that the VLF propagation anomaly in association with a Japanese transmitter, JJI (frequency = 22.2 kHz) located in Miyazaki prefecture [33] is definitely observed before the EQ.This was based on the VLF data observed at a particular station.
Asano and Hayakawa (2017) [34] have made further detailed analysis on the basis of VLF data observed at all 8 stations in Japan on the assumption that the ionospheric perturbation appears first above the EQ epicenter and it is assumed to develop temporally and spatially (horizontally and vertically), we compute the temporal evolutions of amplitude changes at all stations by changing independently the reflection height of one-hop sky wave and that of the second-hop wave close to the transmitter side of the relevant path, and we compare those theoretical values based on the wave-hop method with the observational results in order to infer extensively the spatio-temporal evolution of the ionospheric perturbation.Details on the results will be published elsewhere, but we indicate the following conclusion that the ionospheric perturbations appears about two weeks before the EQ, then it seems to develop by expanding horizontally and vertically (its reflection height becomes lower), and the most perturbed perturbation happens on 10 -12 April (in UT), followed by a decay.

2) ULF depression
The horizontal component of ULF waves observed at KNY has been investigated to investigate the ULF depression [8] [35].The frequency of 0.01 Hz (10 mHz) and local time data are used, and the results of the ULF depression have already been presented in the bottom panel of Figure 1, but they are plotted again on the second panel of Figure 3 during a shorter period before the main shock (the depression in the Z component is also plotted on the 3rd panel in

Summary
We here summarize the observational results on electromagnetic precursors to the 2016 April Kumamoto EQ event.

1) Lithospheric effect Open Journal of Earthquake Research
As for the direct lithospheric ULF radiation, we cannot come to any definite conclusion on the presence of seismogenic ULF emissions simply on the basis of the statistical analysis.

2) Atmospheric effect
Next, how about the indirect effect of seismogenic atmospheric phenomena?
The characteristic parameter ∆S is clearly increased on April 8 to 11 (in UT), which shows the presence of seismogenic ULF/ELF (f = 0.1 -24 Hz) impulsive emissions.Then, the direction finding result for those ELF emissions indicates that the azimuths are approximately directed to Kyushu Island, but it is not clear that the azimuth is not always coincident with the EQ epicenter.

3) Ionospheric effect
The lower ionosphere is definitely perturbed with the use of our subionospheric VLF/LF network.It is found that the ionosphere is extremely perturbed on 10 -12 April (in UT).In exact consistence with this result, the ULF depression (as a signature of lower ionospheric perturbation) at KNY indicated very clear effects on 8 and 10 April (in UT).

Discussion
The conventional statistical analyses used in this paper to explore electromagnetic precursors to the 2016 April Kumamoto EQ event, have indicated clear precursory signatures both in the atmosphere and ionosphere (as ULF/ELF atmospheric radiation and ionospheric perturbation, respectively), but lithospheric ULF radiation is not so evident.The importance of these two phenomena has already been confirmed for the recent 2011 Tohoku EQ, so that the occurrence of these two precursory phenomena among the many, 1) atmospheric ULF/ELF radiation and 2) ionospheric perturbation are likely to be universal for either inland (fault-type) or sea (subduction) EQ, which would be of the first priority for short-term EQ prediction [30].However, we have to mention that the mechanisms for both phenomena are not well understood, and require further investigation [30].Our non-conventional analysis based on critical analysis has yielded the presence of lithospheric ULF radiation for the 2011 Tohoku EQ [26], even though no clear signature in lithospheric ULF radiation was found by the conventional statistical analysis for this EQ.The lithospheric ULF radiation is again not so evident for this Kumamoto EQ event, but our future paper [29] will indicate a possibility of finding out criticality even for this Kumamoto EQ event.
It seems interesting to make a comparison between the present results for the Kumamoto EQ event and those for the 1995 Kobe EQ [3], because both EQs are of the same fault-type (inland EQ) and have the similar magnitude.A lot of electromagnetic precursors have been observed prior to the Kobe EQ.The most definite was the ionospheric perturbation as observed with the use of subionospheric VLF signals [9], which seems to be consistent with the present Kumamoto case.The second important effect for the Kumamoto case is the clear signature of ULF/ELF radiation for the Kumamoto case, but unfortunately there was no ELF measurement at the time of Kobe EQ, but lightning activity (in VLF) seemed to be enhanced for the Kobe EQ.As for lithospheric effect, there was no ULF station close to the EQ epicenter for the Kobe EQ, to be compared with the present Kumamoto case.
Finally it seems to the authors that inland, fault-type EQs are rich in electromagnetic precursors in contrast to sea (subduction) EQs.There are only few publications paying attention to this difference in electromagnetic characteristics between inland and sea EQs.Parrot (2011) [36] is the only one that studied this difference only in the ionospheric parameters.This point requires further extensive studies, because it would be of great importance in studying the whole seismogenic effects and lithosphere-atmosphere-ionosphere coupling mechanism.

Figure 1 .
Figure 1.A map of the magnetic field observatory at Kanoya (KNY) and EQs with M w > 5 (two foreshocks and main shock at Kumamoto as red circles).Also the positions of main fractures and troughs and ULF/ELF stations (NAK, SHI and IZU) are plotted.

Figure 2 .
Figure 2. Evolutions of magnetic field characteristics in the frequency band of 10 -20 mHz during 15 Jan-31 May 2016.The vertical red dashed line shows the time of main EQ.The top panel illustrates the temporal evolutions of geomagnetic (K p ) and seismic (K LS ) activity.Evolutions of horizontal, vertical components and their ratio are shown on 2 -4 panels.Bottom two panels reflect the evolutions of ULF absolute (inverse of F h ) and normalized depression.
Here h A , d A and h ϕ , d ϕ are instantaneous amplitudes and phases of the field component signals.h corresponds to the NS component of magnetic field, while d, the EW component.They are computed from appropriate complex signals which are obtained using from the real signals (U h and U d ) by means of the Hilbert transform.The last ones are extracted from the recorded signals with narrow-band filtration.

Figure 3 .
Figure 3. Temporal evolutions of geomagnetic (K p ) and seismic activity (K LS ) just before the EQ (top panel), 2nd and 3rd panels refer to the ULF depression for H and Z components in the frequency band of 10 -20 mHz at KNY (just for comparison), respectively.The rectangular plot on the bottom panel refers to ∆S at SHI at f = 11 -12 Hz.

Figure 4 .
Figure 4. Small rectangular panels and one large panel on the bottom right illustrate the daily evolutions of main lobes of azimuth distribution of ULF/ELF radiation in the frequency band 11 -12 Hz.Also main EQs and fractures and troughs are shown.Dependence of depth on color and magnitude on dimension of circles are shown to the right of the large panel.

Figure 3
Figure 3 just for comparison).It is clearly recognized from the 2nd plot that the ULF depression (horizontal component) is observed on 8 and 10 April (in UT) without any doubt.