Study on VLF Electric Field Anomalies Caused by Seismic Activity on the Western Coast of the Pacific Rim

Abstract

In order to explore the correlation between earthquakes and ionospheric very low-frequency (VLF) electric field disturbances, this article uses VLF data observed by the China Earthquake Electromagnetic Satellite (CSES) to analyze very low-frequency signals before and after earthquakes from January 2019 to March 2023 in terms of the amplitude and signal-to-noise ratio of electric field power spectrum disturbances. Taking 73 earthquakes with a magnitude of 6.0 or higher occurring in the Circum-Pacific seismic belt as an example, comprehensive research on the VLF electric field disturbance phenomenon caused by strong earthquakes is conducted, considering both the earthquake location and source mechanism. The research results indicate the following: (1) there is a strong correlation between earthquakes with a magnitude of 6.0 or above and abnormal disturbances in the VLF electric field, which often occur within 20 days before the earthquake and within 800 km from the epicenter. (2) From the perspective of earthquake-prone areas, the VLF electric field anomalies observed before earthquakes in the Ryukyu Islands of the Taiwan region exhibit small and concentrated field fluctuations, while the Taiwan Philippines region exhibits larger field fluctuations and more dispersed fluctuations. The discovery of this correlation between seismic ionospheric phenomena and seismic activity provides a new and effective approach to earthquake monitoring, which can be used for earthquake prediction, early warning, and disaster prevention and reduction work.

1. Introduction

China is located at the junction of the Eurasian and Pacific plates, where earthquakes occur frequently. Therefore, studying earthquakes on the western coast of the Pacific Rim is of great significance [1]. Previously, the western Pacific Rim seismic zone was divided into the Ryukyu Islands Taiwan Island seismic zone and the Taiwan Islands Philippines regional seismic zone. Due to the different mechanisms of earthquake occurrences between these two seismic zones, mainly manifested in the different properties of upper and lower stresses, the northeastern Japan and Kuril Kamchatka regions are designated as the Type I double seismic zones, while the northeastern Taiwan and North Island regions of New Zealand are designated as the Type II double seismic zones [2].

Since Gokhberg et al. first reported electromagnetic signal anomalies observed by satellites before earthquakes, research based on ionospheric anomalies observed by satellites began to emerge [3]. Theoretical research has shown that stress changes, microcracks, and microvibrations before earthquakes can cause changes in the ground electromagnetic field. Through induction, the transfer of energy by charged particles in the atmosphere and ionosphere can cause changes and disturbances in the density of particles, such as the spatial ionosphere, known as the coupling effect between lithosphere atmosphere ionosphere coupling (LAIC) [4]. Through coupling effects, the electromagnetic anomalies generated during earthquake preparation can cause ionospheric disturbances and trigger changes in various spatial physical parameters [5,6,7,8,9].

The propagation of very low-frequency and low-frequency electromagnetic signals (VLF/LF) in the lower ionosphere can be used to monitor ionospheric disturbances [10,11]. The advantages of this method lie in its long propagation distance and high sensitivity to the ionosphere, but there are also challenges, such as interference from the Earth’s magnetic field.

Therefore, scientists have deployed multiple ground-based artificial source launch stations worldwide, continuously emitting VLF electromagnetic waves into space at different frequencies. These signals can penetrate the ionosphere and propagate upwards to the top of the ionosphere. At the top of the ionosphere, these signals are reflected back to the ground. Satellites and ground-receiving stations can synchronously receive these signals [12]. The VLF electromagnetic waves emitted by artificial sources have the characteristics of a long transmission distance and low energy loss, which can propagate over long distances in the Earth’s ionosphere waveguide system and have obvious wave-particle interaction effects [13,14,15]. When a satellite is flying over an artificial source launch station, it can receive a stable artificial source signal in a specific frequency range [16]. Electromagnetic wave propagation in the ionosphere has complex properties, which can affect the expression of the electromagnetic response [17]. It is known from the literature [18] that most of the many artificial source-launching stations are located north of the equator, while the North West Cape (NWC) station is in Australia. It is the only artificial source transmitter station south of the equator still in operation and is located on the northwest coast of Australia (geographic coordinates 21.82°S, 114.17°E, transmitter power: 1000 kW, transmitter frequency: 19.8 kHz), and the latitude of its location has also induced a strong effect on the radiation band [19].

The NWC transmitting station provides VLF band radio communication signals to warships and submarines and is by far the largest transmitting station in the southern hemisphere in terms of transmitting power. The NWC transmitting station continuously transmits a fixed-frequency VLF signal into space, which propagates upwards and penetrates the ionosphere [20]. The VLF signal in the ionosphere is recorded when the CSES passes over and near the NWC station, and significant electric field radiation can also be detected in the magnetic conjugation zone [21].

The earthquake monitoring satellite Detection of Electromagnetic Emissions Transmitted from Earthquake Regions (DEMETER), launched by the European Space Agency on 1 October 2004, observed many typical ionospheric disturbances during its operation and discovered significant changes in various parameters for studying the ionosphere before earthquakes, such as the 2005 Sumatra 8.6 earthquake, the 2008 Wenchuan 8.0 earthquake, identifying that each parameter of the ionosphere changed significantly before the earthquakes, such as electric field strengths, high-energy particle flux, electron density, etc. [22,23,24]. Moreover, the statistical analysis work of seismic events recorded during the operation period obtained a series of results reflecting the characteristics of seismic ionospheric phenomena. Yan R et al. [25] used the statistics of ion density observed using the DEMETER satellite between 2004 and 2010 to find that the anomalous pre-seismic perturbations mainly appeared in the first 5 days before the earthquake and within 200 km from the epicenter. Zhang et al. [26] collected and analyzed 69 strong earthquakes of magnitude 7 and above that occurred between January 2005 and February 2010; 46% of the strong earthquakes had ULF and ELF electric field disturbances detected prior to the earthquake; 35% of the strong earthquakes had maximum disturbance amplitudes exceeding three times the standard deviation in the 10 days prior to the earthquake. And the earthquakes often occurred in the time interval with the strongest disturbance amplitudes, where the values of disturbance amplitudes dropped back to their normal range after the earthquake.

CSES is the first electromagnetic monitoring test satellite launched by the China Earthquake Administration on 2 February 2018 and is the first star of China’s geophysical field exploration satellite program. CSES is mainly used to obtain dynamic quasi-real-time monitoring data of the global electromagnetic field, ionospheric plasma, energetic particles, and other physical quantities. CSES provides new technical means to explore earthquake precursor information, space environment monitoring and forecasting, and Earth system scientific research and is making technical preparations for the establishment of an operational satellite system for the electromagnetic monitoring of earthquake precursors in the future. The VLF band in CSES is 1.8–20 kHz. CSES carries eight types of payloads, of which the electric field detector (EFD) has four different detection frequency bands; the sampling rate of the VLF band is 50 kHz, and the sampling period is 2.048 s, so there are 2048 sampling points in each working cycle [27,28].

These studies of anomalous changes in electric and magnetic fields before earthquakes cover almost the entire electromagnetic wave band from DC to HF [29,30,31]. This suggests that pre-seismic electric and magnetic field anomalous changes are not only limited to a specific frequency band but may occur in the entire electromagnetic wave band. This provides new ideas and methods for earthquake prediction, which is of great significance for improving the accuracy and reliability of earthquake prediction. This article only studies the VLF electric field. Accompanied by the continuous development of satellite monitoring methods and technologies, more and more researchers have found that low-orbit satellites can monitor VLF/LF signals entering into the ionosphere, which is an important reference value for the study of seismic electromagnetic waves penetrating into the ionosphere, the signal of ionospheric electric field perturbation [32,33], and the signal of ionospheric electric field perturbation. Although the relationship between VLF disturbances and earthquakes has been statistically analyzed in previous studies, these statistics have not been classified according to the different source mechanisms of earthquakes and the region of seismicity.

Domestic and foreign scholars have made some progress in the study of pre-seismic electromagnetic anomalies, providing rich information for people to understand pre-seismic electromagnetic anomalies. CSES differs from DEMETER satellites in terms of the operating altitude and orbital period; therefore, the statistical study of very low-frequency electric fields observed by CSES is innovative. It is absolutely necessary and meaningful to conduct statistical research on the VLF electric field recorded by the CSES, which provides new ideas and methods for earthquake prediction.

There are some significant questions, such as is the anomalous disturbance of the ionosphere caused by earthquakes consistent with the different earthquake source mechanisms? What are the characteristics of changes om various space physics parameters? It is interesting to study and discuss these questions, which are of practical significance for earthquake monitoring and prediction.

2. Data Sources and Selection

Since the launch of the CSES, it has been operating both stably and normally and has now generated a huge amount of space physics data by recording the ionospheric disturbance caused by all the seismic events that have occurred during this period, which provides valuable data information for this research.

2.1. Selection of Study Area and Earthquake Magnitude

The West Coast of the Pacific Rim seismic zone described earlier has different mechanisms for generating earthquakes and is an ideal region to study whether the ionospheric disturbances caused by different earthquake source mechanisms are consistent. Earthquakes have a very rapid impact on the ionosphere, but not all earthquakes can produce obvious ionospheric anomalies. The larger the magnitude, the easier it is to observe the ionospheric anomalies before the earthquake. For this reason, this paper adopts the data of earthquakes of a magnitude 6.0 and above on the Richter scale within the study area from January 2019 to March 2023, published on the website of the United States Geological Survey, USGS (https://earthquake.usgs.gov/earthquakes/search/, accessed on 1 May 2023), and discards earthquakes with an epicenter depth greater than 70 km because the electromagnetic radiation from rock ruptures is difficult to radiate from a deeper subsurface [34,35]. Meanwhile, considering that ionospheric disturbances in high latitudes are more complicated [36], only earthquakes that occurred in middle and low latitudes were selected, i.e., earthquakes with latitudes outside the 50° north and south latitudes were discarded. In order to improve the authenticity of the statistical analysis of earthquake cases, this paper only studies the main shock of each earthquake, i.e., only the earthquake with the largest magnitude within the spatial range of latitude and longitude ± 2° and time ± 15 days were counted. After this restriction, the final seismic examples included 73 seismic events, the geographical distribution of which is shown in Figure 1.

The time range was chosen in accordance with the period from January 2019 to March 2023 after the normal operation of the CSES to the present time. According to the earthquake data released by the USGS, there were 73 cases of magnitude 6 or above earthquakes during this period. For each earthquake case, we used data from 65 to 30 days before the earthquake to construct the background field and analyze the signal-to-noise ratio 30 days before and 5 days after the earthquake.

2.2. Spatio-Temporal Selection of Data

The extent of ionospheric anomaly disturbances caused by seismic events depends on the magnitude of the earthquake and the time of occurrence. The greater the seismic magnitude, the greater the extent of the impact on ionospheric anomalies. Therefore, when determining the spatial range, it is necessary to divide different ranges according to the earthquake’s magnitude. In order to ensure the comparability of these statistical results, the range of anomalous disturbances of seismic events of the same magnitude needs to be unified when determining the spatial extent. The size of the preparatory seismic zone can be calculated according to the preparatory seismic Formula (1):

$$R={10}^{0.43M}$$

(1)

where $R$ represents the radius of the preparatory zone in km, and $M$ represents the magnitude [37]. For an earthquake of magnitude 6 on the Richter scale, the diameter of the preparatory zone is about 380 km; for an earthquake of magnitude 7, the diameter reaches 1023 km. Based on the experience of the previous researchers [34], the study area was chosen as a detection range with a radius of 1000 km centered on the epicenter of the earthquake.

In the time range selection, we refer to the scope of previous research time selection studied more than 30 days before the earthquake to 5 days after the earthquake [22,23].

Therefore, the range from 30 days before the earthquake to 5 days after the earthquake was chosen as the research period.

In order to exclude the disturbance of the ionosphere caused by solar activity and other factors, the following disposition was made: if the geomagnetic activity index is Kp > 3 on the day when the disturbance occurs, all the orbital data of that day are deleted so as to avoid anomalies due to the excessive geomagnetic activity interfering with the statistical results. Because of the weak geomagnetic activity at night side, nighttime data of the electric field recorded by the CSES were selected for this study.

3. Extraction Anomaly Methods for VLF Data

3.1. Signal-to-Noise Ratio Variation Characterization Method

The signal-to-noise ratio (SNR) of a signal reflects the relationship between the strength of the propagated signal, and a larger SNR indicates a stronger useful signal and better signal quality. According to the calculation methods proposed in the relevant literature [23,38], the formula for calculating the SNR of the VLF artificial source transmitting station in the preparatory seismic range is shown in Formula (2):

$$SNR=\frac{2A\left(f_0\right)}{A\left(f_{+}\right)+A\left(f_{-}\right)}$$

(2)

where $A\left(f_0\right)$ is the power spectral density value corresponding to the fixed transmitter frequency of the artificial source transmitter station, $A\left(f_{+}\right)$, $A\left(f_{-}\right)$ are the power spectral density values corresponding to the upper and lower limit frequencies of the transmitter station, and the bandwidth of the transmitter station’s frequency is determined by its transmitter power, as shown in Figure 2.

In this research, we analyzed the frequency bands of NWC stations. The analysis steps are as follows:

• The E-ab component observation data of the three components of the electric field are obtained within a range of 10° around the epicenter from 30 days before the earthquake to 5 days after the earthquake.
• For each track data obtained in the first step (current track data C0), the revisited tracks within the previous 35 days are obtained, which means that C0 data for each track correspond to a set of seven revisited tracks.
• The VLF artificial station used is the NWC station located in Australia. Formula (2) is used to calculate the SNR of the VLF frequency band of the space electric field recorded by the satellite when passing through the earthquake preparation area.
• At the same time, based on the revisited period of CSES, the SNR changes are obtained for a total of 7 periods of 5 days, including 30 days before and 5 days after the earthquake.

3.2. Methods for Analyzing the Ionospheric Perturbation Amplitude

According to the relevant literature, the perturbation amplitude of the electric field power spectrum relative to the background can be obtained by subtracting the mean value of the background field power spectrum from the mean value of the electric field power spectrum and then dividing it by the standard variance matrix of the background field power spectrum [39]. The formula for calculating the perturbation amplitude of the ionosphere is shown in Formula (3):

$$\theta =\frac{\alpha -\beta }{\gamma }$$

(3)

where $\beta$ is the mean matrix of the background field constructed using data from 65 to 30 days before the earthquake, $\gamma$ is the standard variance matrix of the background field, matrix $\alpha$ calculates the disturbance amplitude value for each earthquake preparation period (five days), and $\theta$ represents the perturbation amplitude matrix of the seismic spatial electric field power spectrum with respect to the background field, each value of which is the ratio of the number of the corresponding position of the matrix. We used each value in the $\theta$ matrix as the ionospheric disturbance amplitude value at the corresponding location.

The steps of this statistical study are similar to the methodological steps described above, differing only in the methodological formulas applied and no difference in the data used, so the presentation of the methodological steps is not repeated.

4. Statistical Results and Analysis

The process of generating earthquakes involves a variety of factors, such as the rock structure, geophysics, and geochemistry of the Earth’s interior, and, therefore, the abnormalities that result from them are inevitably affected by a variety of factors that manifest themselves in a complex and varied manner.

In the following, the relationship between different regional types of seismic events and VLF electric field anomalies are statistically analyzed from the perspectives of spatial distribution and the spatial-temporal range of anomalies disturbance. We expected to obtain a common feature of pre-seismic anomalous disturbances from the results of the statistical analysis of multiple seismic events.

4.1. Characterization of the Spatial SNR Triggered by Earthquakes of Different Regional Types

Previous studies have shown that during the seismic gestation period, the medium and magnetic field has an effect on the propagation of electromagnetic signals during the upward propagation of the seismic source, and different geological formations and seismic source mechanisms between the two regions give rise to different ionospheric perturbation phenomena. In order to compare the statistical results with those of the different regions in the later section, in this paper, we first conduct the statistical analysis on the pre-seismic SNR of all the selected seismic cases, and the results of the statistical analysis are shown in Figure 3.

The discontinuity between 360 degrees and 0 degrees here is because we did not perform cyclic interpolation on the matrix when processing the data. Azimuth is the angle between the north direction and the direction to a point of interest, measured in a clockwise direction from the north. It is used to describe the direction of an object and is typically expressed in degrees or kilometers.

Figure 3 gives the spatial distribution of pre-seismic SNR relative to the background in the total earthquake case preparation area. It can be seen that the pre-seismic SNR intensity is mostly concentrated in the north within 200 km from the epicenter, widely distributed within 400 km from the epicenter, and rarely distributed within 600 km–1000 km from the epicenter. The maximum SNR value reaches 3.2σ, and the minimum value is −1.92σ. The variation in SNR within the entire preparation seismogenic area is basically within the range of 600 km from the epicenter, and the numerical range is mostly between 2.5 × 10⁻³ and 3.3 × 10⁻³. The strongest variation is located between 315 degrees and 0 degrees of the epicenter azimuth angle.

4.1.1. Spatial Comparison Analysis of the SNR in Preparatory Seismic Regions

In order to analyze the similarities and differences in the spatial distribution of pre-seismic SNR within the preparatory seismic region of the two types of regions, this study collected data on the variation in the SNR and analyzed them using the superposed epoch analysis and the results of the statistical are shown in Figure 4 and Figure 5.

In the Ryukyu Islands—Taiwan Island region the pre-seismic SNR in the preparatory space is strong on the right side within 200 km of the epicenter, on the left side at 500 km from the epicenter, and the SNR is lower than that of the Taiwan Island—Philippine Islands region. Whereas the pre-seismic SNR in the preparatory space of the Taiwan Island—Philippine Islands region is centered in the area within 400 km, the SNR is present at the epicenter in the area of 600–800 km but is lower. In the Ryukyu Islands—Taiwan Island region, the maximum SNR value reaches 3.4σ, and the minimum value reaches −1.42σ. In the Taiwan Island-Philippine Islands region, the maximum SNR value reaches 2.75σ, and the minimum value reaches −2.0σ. The degree of variation in the SNR ratio before earthquakes and the differences in regions could be related to the different geological structures in these two regions.

4.1.2. Spatial and Temporal Comparison Analysis of the SNR Ratio in Preparatory Seismic Area

The VLF signals generated by earthquakes interfere with the VLF artificial source signals in the process of upward propagation, thus causing the phenomenon of SNR ratio reduction. Different source mechanisms may produce different SNR ratio anomalous perturbation phenomena, so it is necessary to count regional earthquakes with different source mechanisms separately. The seismic examples used in this experiment include 38 earthquakes in the Ryukyu Islands—Taiwan Island region and 35 earthquakes in the Taiwan Island—Philippine Islands region. The results of the superposed epoch analysis of the statistics of the seismic cases of the two regional earthquakes are shown in Figure 6. The color bar here shows the level of SNR. The cycle here refers to the revisit cycle of CSES, and each revisit cycle is 5 days. The research time frame is 35 days, so there are 7 periods in total with one picture for each period, so 7 pictures are shown here.

From Figure 6, it can be seen that there is a significant change in the SNR in the fifth period before the earthquake. The change in the signal-to-noise ratio in the fourth period before the earthquake is very small, ranging from 0.4 × 10⁻³ to 0.8 × 10⁻³. The change in the SNR in the third period before the earthquake is generated near the epicenter, with a maximum value of 1.4 × 10⁻³, concentrated at a distance of 200 km from the epicenter. The spatial variation range of SNR in the second period before the earthquake is expanded compared to the previous period and concentrated at a distance of 600 km from the epicenter. The change in the SNR in the first period before the earthquake still exists, but it is relatively small, with an average value of around 0.8 × 10⁻³. During the earthquake, a significant change in the SNR is observed at the epicenter, with a maximum value of 1.5 × 10⁻³ and a diffuse distribution in space. The change in the SNR in the first period occurred after the earthquake slowed down.

The results of the superposed epoch analysis of earthquake case statistics for the Ryukyu Islands—Taiwan Island region are shown in Figure 7.

As shown in Figure 7, the changes in the SNR in the fifth period before the earthquake exist and are concentrated within a range of 200 km and 500 km from the epicenter, with a maximum value of 1.1 × 10⁻³. The changes in the SNR in the fourth period before the earthquake decreased to 500 km from the epicenter compared to the previous period. The changes in SNR near the epicenter in the third period before the earthquake ranged from 0.6 × 10⁻³ to 0.9 × 10⁻³. And the spatial range was wide, up to 600 km away from the epicenter. The variation in the SNR in the second cycle before the earthquake ranged from 0.2 × 10⁻³ to 0.4 × 10⁻³ within the range of 400 to 200 km from the epicenter. The maximum change in the SNR in the first cycle before the earthquake reached 1.1 × 10⁻³. The change in the SNR during the earthquake period was not significant and mostly concentrated within a range of 200 km from the epicenter. The SNR in the first period after the earthquake still varied, which we believe is related to aftershocks or other additional geological activities after the earthquake.

The results of the superposed epoch analysis of the earthquake case statistics for the Taiwan Island—Philippine Islands region are shown in Figure 8.

It can be seen from Figure 8 that the SNR of the earthquake preparation area in the Philippine region of Taiwan is greater than that of the Ryukyu Taiwan region. The SNR variation in the fifth cycle before the earthquakes is concentrated within a range of 500 km, with a maximum value of 1.8 × 10⁻³ and an average of 0.7 × 10⁻³. The spatial distribution characteristics of the changes in the SNR in the fourth period before the earthquake are similar to the previous period, with a maximum value of 1.5 × 10⁻³ and an average of 0.64 × 10⁻³. The variation in the SNR in the third period before the earthquake is concentrated in the northeast within 200 km of the epicenter, but its maximum value is smaller than that in the fifth period before the earthquake with an average of 0.7 × 10⁻³. The spatial distribution of SNR changes in the second period before the earthquake is diffuse but still concentrated within a range of 600 km from the epicenter. The maximum SNR in the first period before the earthquake reached 1.8 × 10⁻³, with an average of 0.7 × 10⁻³, and was spatially distributed within 400 km of the epicenter. This variation in the periodic SNR during earthquakes is distributed within the azimuth range of 225 degrees to 45 degrees. The first cycle changed after the earthquake weakened. Through the comparison between Figure 7 and Figure 8, it was found that the spatial distribution characteristics and temporal trends of the SNR changes were different, which could be related to different geological structures in the two study areas.

4.2. Analysis of the Amplitude of Ionospheric Disturbances Triggered by Earthquakes of Different Regional Types

Previous studies have shown that electromagnetic signals propagate upwards from the source of an earthquake during its gestation period and at the time of its occurrence. In the propagation process, electromagnetic signals are interfered with by the medium and magnetic field. The dielectric constant and magnetic permeability of the medium affect the propagation speed and direction of electromagnetic waves, while the magnetic field affects the propagation direction of electromagnetic waves. Therefore, the medium and magnetic field affect the propagation distance, propagation time, and propagation direction of electromagnetic signals during the upward propagation of the seismic source. Moreover, different geological structures and earthquake source mechanisms between the two regions cause different ionospheric disturbance phenomena. Therefore, this study divides the seismic examples based on different source mechanisms and superimposes the results of the superposed epoch analysis, as shown in Figure 9.

As shown in Figure 9, the spatial evolution characteristics of the ionospheric disturbance amplitude before the total earthquake compared to the background. From Figure 9, it can be seen that the ionospheric disturbances in the pre-seismic area are the strongest near the epicenter, concentrated in the northeast region within 200 km from the epicenter, with a maximum value of 13 × 10⁻⁴. These disturbances are distributed in the pre-seismic area in the form of ripples and diffusion. Ionospheric disturbances exist within the range of 800 to 1000 km from the epicenter, and the disturbance values are in the range of 7 × 10⁻⁴ to 8 × 10⁻⁴. At 500 km from the epicenter, changes in ionospheric disturbances can also be observed, with a maximum value of 10 × 10⁻⁴. And this indicates that θ, the maximum value, reaches 4.0σ.

4.2.1. Spatial Comparison Analysis of the Amplitude of Ionospheric Disturbances in the Preparatory Seismic Region

In order to analyze similarities and differences in the spatial distribution of pre-seismic ionospheric disturbances within the preparatory seismic region of these two regions, this article collects ionospheric change data and analyzes them using the superposed epoch analysis method, as shown in Figure 10 and Figure 11, which show the spatial distribution of pre-seismic ionospheric disturbances in the two regions, respectively.

From Figure 10, it can be seen that the ionospheric disturbances in the earthquake preparation period region are the strongest at the epicenter. The ionospheric disturbances in the Ryukyu Islands—Taiwan Island region are within 200 km of the epicenter, with azimuth angles ranging from 0 to 180 degrees. There is a circle of disturbances around the epicenter and 500 km away from the epicenter, with disturbance amplitude values ranging from 12 × 10⁻⁵ to 13 × 10⁻⁵. There is a circle of weaker disturbances at 800 km from the epicenter, with disturbance amplitude values ranging from 9 × 10⁻⁵ to 11 × 10⁻⁵. The rest of the region is relatively calm, with disturbance levels ranging from 10 to 5. And it is shown that the maximum value of θ reaches 3.3σ. The average amplitude of ionospheric disturbances is 2.8 × 10⁻⁵, and the maximum value is 15 × 10⁻⁵.

As seen in Figure 11, the regional ionospheric disturbances in the Taiwan Island—Philippine Islands region are stronger in the northeast direction within 200 km of the epicenter, and the disturbances at 400 km from the epicenter are mostly located in the northwest direction of the epicenter. And the ionospheric background disturbances are more frequent, with the disturbance magnitude in the order of 10⁻³. It shows that the maximum θ value reaches 4.48σ. According to the comparison between Figure 10 and Figure 11, we found that the spatial distribution of ionospheric disturbances is different between the two, and the disturbance values are also significantly different, which could be related to their different geological structures.

4.2.2. Temporal and Spatial Comparative Analysis of the Amplitude of Pre-Seismic Ionospheric Disturbances in the Preparatory Seismic Region

In order to explore the similarities and differences of the ionospheric perturbations in the preparatory seismic zone of earthquakes in these two regions with respect to time, this article collects the ionospheric perturbation data of the preparatory seismic zone in different study cycles and restricts the color maps of the two regions to be in the same range. The following comparative temporal and spatial analyses of the amplitude of pre-seismic ionospheric disturbances in the preparatory seismic region are shown in Figure 12, Figure 13 and Figure 14.

As shown in Figure 12, the strongest value of the ionospheric disturbance amplitude in the fifth period before the earthquake within a range of 200 km from the epicenter is within the azimuth range of 0 to 90 degrees, and the maximum value is 8 × 10⁻⁴. The average amplitude of the disturbance in the fourth period before the earthquake increases within 500 km of the epicenter. The strongest value of the ionospheric disturbance amplitude in the third period before the earthquake, within the range of 200 to 400 km from the epicenter, is located in the azimuth angle range of 270 to 360 degrees with a maximum value of 8 × 10⁻⁴. The disturbance amplitude in the second period before the earthquake fills the space within 400 km of the epicenter, and the strongest value is distributed 400 km away from the epicenter. A very high amplitude of disturbance in the first period before the earthquake can be observed within 1000 km of the epicenter, ranging from 6 × 10⁻⁴ to 8 × 10⁻⁴. The amplitude of periodic disturbances during earthquakes ranges from 5 × 10⁻⁴ to 6 × 10⁻⁴ at 800 km from the epicenter, and a disturbance amplitude of 8 × 10⁻⁴ can be observed at 400 km from the epicenter. The disturbance amplitude of the first period after the earthquake ranges from 6 × 10⁻⁴ to 8 × 10⁻⁴ within the range of 800 to 1000 km from the epicenter, but the disturbance amplitude is almost zero within the range of 600 to 800 km from the epicenter.

As shown in Figure 13, the amplitude of ionospheric disturbances in the fifth period before the earthquake is concentrated within 600 km of the epicenter, with a spatial mean of 9.1 × 10⁻⁴. These disturbances are mostly distributed at 200 and 400 km of the epicenter and exhibit a banded distribution feature. The ionospheric disturbance in the fourth period before the earthquake is distributed within 1000 km of the epicenter, with the maximum value within 200 km of the epicenter and an azimuth angle of 0 to 90 degrees. The ionospheric disturbance in the third period before the earthquake is distributed within a range of 400 km from the epicenter and presents a three-segment distribution. The first segment is located at an azimuth angle of 0 to 90 degrees, the second segment is located at an azimuth angle of 130 to 230 degrees, and the third segment is located at 270 to 360 degrees. The spatial distribution of the ionospheric disturbance amplitude in the second period before the earthquake is similar to that in the previous period, with a maximum disturbance amplitude of 3 × 10⁻⁴. The spatial distribution of the ionospheric disturbance in the first period before the earthquake is similar to that in the previous period, with the maximum amplitude of the disturbance located within 200 km of the epicenter. The amplitude of periodic ionospheric disturbances during earthquakes ranges from 100 to 200 azimuth angles at 900 km from the epicenter, with a disturbance value of 2.3 × 10⁻⁴. The average amplitude of ionospheric disturbance in the first period after the earthquake is 2.8 × 10⁻⁴ within a range of 1000 km from the epicenter, indicating a small overall disturbance amplitude.

From Figure 14, it can be seen that the ionospheric disturbance in the fifth period before the earthquake has a disturbance amplitude of 7 × 10⁻⁴ at a distance of 500 km from the epicenter, and there are disturbances ranging from 3 to 6 within a range of 600 to 800 km from the epicenter. The ionospheric disturbance in the fourth period before the earthquake is within a range of 5 × 10⁻⁴ to 8 × 10⁻⁴ within 400 km of the epicenter, with an average of 3.3 × 10⁻⁴. The amplitude of ionospheric disturbances in the third period before the earthquake is concentrated within a range of 400 km from the epicenter, with values ranging from 3 × 10⁻⁴ to 4 × 10⁻⁴ within a range of 600 to 1000 km from the epicenter. The amplitude of ionospheric disturbances in the second period before the earthquake reached a maximum of eight within a range of 400 km from the epicenter. The amplitude of the first cycle disturbance before the earthquake is within 600 to 1000 km of the epicenter, with values ranging from 6 × 10⁻⁴ to 8 × 10⁻⁴ and widely distributed in space. The periodic ionospheric disturbances during earthquakes are mostly concentrated within a range of 400 km from the epicenter, with a maximum value of 8 × 10⁻⁴. The disturbance amplitude values at 800 km from the epicenter range from 5 × 10⁻⁴ to 6 × 10⁻⁴. The amplitude of ionospheric disturbances in the first period after the earthquake reaches a range of 6.5 × 10⁻⁴ to 7.5 × 10⁻⁴ at a distance of 800 km from the epicenter. From the comparison between Figure 13 and Figure 14, it can be seen that there are significant differences in the spatial distribution and temporal distribution of the disturbance amplitude between these two regions, which could be related to different geological structures in the two regions.

5. Discussion

The statistical results of this paper show that the time of occurrence of anomalies of earthquakes with a magnitude of 6.0 or above in the study area is within the range of about 30 days before the earthquake to 5 days after the earthquake, and these anomalies occur within 800 km of the epicenter, which is basically in line with the results obtained from the previous studies [23].

In the studied earthquake cases, the seismic zone division identifier given by other scholars is divided into two regions for study [1], and the results show that the anomalies in the Ryukyu Islands—Taiwan Island region have a smaller magnitude and range, mostly occurring within 400 km, while the anomalies in the Taiwan Island—Philippine Islands region have a larger magnitude and range, mostly occurring within 800 km.

According to Formula (1), it is reasonable that the anomalies generated in the process of seismicity are concentrated in the range of 800 km, with the most obvious ones in the range of 400 km. These two regions are located in the same longitude, with a large difference in latitude, and the electric fields in the breeding zone exhibited have their own characteristics. From the perspective of SNR analysis, the left side of these two regions is stronger than the right side, and it is interesting to note that in the 30 days before the earthquake, as the time of the earthquake approaches, the distribution of the SNR shows a ripple-like signal distribution centered on the epicenter and spreading in all directions.

This may be an expressive feature of the ionosphere’s absorption of electromagnetic signals triggered by earthquakes, and it is conjectured that, based on the dynamic change in this ripple phenomenon over time, we can then calculate at which stage the current earthquake gestation is in. Analyzing from the perspective of ionospheric disturbance, the responses to earthquakes in these two regions are very different, with a slight disturbance in the Ryukyu Islands—Taiwan Islands region. By contrast, the Taiwan Island—Philippine Islands region showed a significant disturbance amplitude and time-dependent ripples within 30 days before the earthquake. Anomalous disturbances were also found after the earthquake, which could be due to the tectonic plates still colliding with each other after the main shock. Such disturbances may be caused by tectonic plates still colliding with each other after the main earthquake or by secondary hazards such as tsunamis and volcanic eruptions caused by the earthquake.

Specifically, after the mainshock occurs, the tectonic plates are still colliding with each other to generate a large amount of energy, which is transmitted to the ionosphere in the form of seismic waves and electromagnetic waves, thereby causing ionospheric disturbances. And secondary disasters such as tsunamis and volcanic eruptions caused after the earthquake also generate a large amount of energy, which affects the ionosphere, thus causing ionospheric disturbances. According to data in this study, 70% of the earthquake cases occur in the sea. Earthquakes occurring in the sea area are more likely to cause secondary disasters such as tsunamis and volcanic eruptions, so these secondary disasters also cause ionospheric disturbances.

From the perspective of the earthquake source mechanism, the Ryukyu—Taiwan region belongs to the upper-stress zone, which is a Type I double seismic zone dominated by compressive stress. Abnormal seismic disturbances often occur within a range of about 400 km near the epicenter, and this amplitude is relatively small. By contrast, the Taiwan Philippines archipelago region belongs to the Type II dual seismic zone, with the upper stress mainly being tensile stress; thus, abnormal disturbances often occur around 800 km from the epicenter and have a large amplitude.

In summary, the causes of post-seismic ionospheric disturbances may include the collision between tectonic plates after the main earthquake and secondary disasters such as tsunamis and volcanic eruptions caused by the earthquake [40,41]. Earthquakes are formed by plates colliding and crushing each other, and different seismogenic mechanisms may cause different ionospheric perturbation phenomena. And displacements and faults caused by earthquakes may also cause ionospheric perturbation in the process of post-earthquake crustal recovery. Therefore, how to screen out these ionospheric disturbances in relation to earthquakes needs further research.

6. Conclusions

In this paper, the space electric field VLF data observed by CSES for a total of 51 months from January 2019 to March 2023 were used for analysis, and from 73 earthquake cases of magnitude 6.0 and above, data were processed using the SNR ratio method with the ionospheric perturbation method, and then statistically analyzed via superposed epoch analysis. This study was carried out in terms of both spatial and temporal scopes and spatial distribution characteristics and statistically analyzed according to the location of the earthquakes, where the following two conclusions were obtained:

(1)

In terms of the number of anomalies that appeared, the electric field anomalies in the Ryukyu Islands—Taiwan Island seismic region are small and concentrated, and the electric field anomalies in the Taiwan Island—Philippine Islands seismic region are large and scattered;

(2)

In terms of the timing of the anomalies, the anomalies in the Ryukyu—Taiwan seismic region occur in the fourth pre-seismic cycle and the seismic time cycle, while the anomalies in the Taiwan—Philippine Islands seismic region occur in the fourth pre-seismic cycle to the first pre-seismic cycle.