Ionospheric Disturbances Possibly Associated with Yangbi Ms6.4 and Maduo Ms7.4 Earthquakes in China from China Seismo Electromagnetic Satellite

: In this paper, the data of electron density, electron temperature and ion compositions detected by the China seismo-electromagnetic satellite (CSES) before Yangbi Ms6.4 and Maduo Ms7.4 earthquakes in China during 21–22 May 2021 are analyzed by using the algorithms of moving average to get the high frequency perturbations, revisited orbit comparison to construct the normal background, and wavelet transform to obtain the spectrum characteristics. Most of the parameters we studied were characterized by synchronous changes and similar frequencies, and even have symmetrical structures in the seismogenic zone, suggesting that they may originate from the same disturbing source. The research results also demonstrated that the electron density and electron temperature forms changed from relatively stable to unstable about 40 days before the earthquakes, and then gradually recovered after the earthquakes. Meanwhile, the disturbances are concentrated near the epicenter. As the earthquakes approach, the number of disturbances gradually becomes more frequent. Based on those observational facts, we suggest that the Yangbi and Maduo earthquakes may have affected the ionosphere through the ionosphere and magnetosphere branch and thermal branch of the lithosphere atmosphere ionosphere coupling model mechanism.


Introduction
In the process of earthquake preparation, energy is dynamically accumulated and transmitted due to crust movement. During an earthquake, the tectonic near the source is destroyed. Such changes before or during earthquakes may affect the physical and chemical processes of the lithosphere, atmosphere and ionosphere, which may be observed by satellites. These variations before the main earthquakes are considered as earthquake precursors [1]. Pulinets and Ouzounov [2] proposed a lithosphere atmosphere ionosphere coupling (LAIC) model based on many observations and experimental results. This model demonstrates that the plasma and magnetic fields in the ionosphere are ultimately affected by thermal, cloud, ionosphere and magnetosphere during earthquake preparation. Among them, the variations in atmosphere electricity induce variations in the ionosphere, including the electron and ion concentration, electron and ion temperature, ion composition and height scale. The research results of many scholars have also confirmed that the electron density (Ne), electron temperature (Te) and ion density (Ni) in the ionosphere will change before the earthquake. These changes may increase or decrease and sometimes may increase and then decrease repeatedly [3][4][5][6][7], and the cumulative number of anomalies increase significantly before the earthquake [5]. In terms of ion compositions, such as H + , He + and O + , different scholars have observed that they will change annually or seasonally before the earthquake by using satellite in situ measurement data and other data [8][9][10]. Akhoondzadeh et al. [1] statistically analyzed the changes of Ne and Ni before strong earthquakes (M > 6.0) by using DEMETER satellite and GPS data, and illustrated that there would be some anomalies in the ionosphere 1-5 days before the earthquakes. By using the SWARM satellites to study 12 strong earthquakes from 2014 to 2016, De Santis et al. [11] illustrated that the magnetic field and Ne observed by the satellites have a statistical relationship with earthquakes.
In order to study the changes in the ionosphere before the Yangbi and Maduo earthquakes, we use the data of CSES to analyze the Ne, Te and ion compositions. We only analyze the nightside data because the nightside is less affected by solar activity and the anomaly is obvious. In addition, according to the 'train radius' formula, ρ = 10 0.43M km , proposed by Dobrovolsky et al. [3], the seismogenic zone of Yangbi earthquake is completely within the seismogenic zone of Maduo earthquake, and the occurrence time of them is only 4 h apart. In such a situation, it is difficult to distinguish whether the possible ionospheric disturbances are from Yangbi or Maduo, or both. Therefore, the two earthquakes are analyzed as the same event in this paper. Nonetheless, we still hope to have new discoveries, because most of the previous studies, except for the statistical studies, are for a specific earthquake, but the two strong earthquakes we studied are close together and the time interval is short. Figure 1 shows the epicenter, seismogenic zone and satellite revisit orbits of the earthquakes. This paper is structured as follows: Section 2 includes the data of CSES and methods implementation. Sections 3 and 4 are results and adequate discussions, respectively. Finally, conclusions are presented in the last section.

Space Environment and CSES Data
In this paper, to discriminate the seismo-ionospheric perturbations from geomagnetic disturbances, the geomagnetic indices Ap (http://www.sepc.ac.cn, accessed on 29 December 2021) and Dst (https://www.nssdc.ac.cn, accessed on 29 December 2021) were checked, and we will only analysis the data in quiet geomagnetic conditions (Ap < 20; Dst > −30 nt). In addition, the F107 (http://www.sepc.ac.cn, accessed on 29 December 2021) is below 100 and does not change much. The Ap, Dst and F107 indexes from five weeks before the earthquake to one week after the earthquake, that is, from 15 April to 31 May, are shown in Figure 2. CSES-1 is Chinese first satellite that is dedicated to monitoring ionospheric disturbance caused by earthquakes. It transits in a solar synchronous orbit with an altitude of 507 km and a descending node at local time 14:00. In order to more effectively realize the coverage monitoring of Chinese domestic seismic belts, the satellite orbits has been specially optimized for the Chinese seismic belts, so that the satellite orbits can cover the North-South seismic belts at the zenith angle to the greatest extent. The payload of the satellite consists of eight kinds of scientific detection instruments, including High-Precision Magnetometer (HPM), Electric Field Detector (EFD), Search Coil Magnetometer (SCM), GNSS Radio Occultation (GRO), High-Energy Particle package (HEP), Langmuir Probe (LAP), Plasma Analyzer Package (PAP) and Tri-Band Beacon (TBB). At the same time, it is also equipped with a High-Energy Particle Detector (HEPD) developed by Italy [14,15].
According to the LAIC model, the plasma and magnetic field in the ionosphere are obviously affected by earthquakes. However, the HPM of CSES-1 was seriously interfered with by TBB due to design reasons. Therefore, in this paper, we obtained the level 02 data of the LAP and PAP data of CSES from 1 January to 15 June 2021. The LAP mainly detects electron density in the range of 5 ×10 2~1 ×10 7 cm − 3 with an error of less than 2‰ and electron temperature in the range of 500~10,000 K with an error of less than 1‰. The PAP has three ion component data: H + , He + and O + ; the detection range of ion density and temperature is consistent with the LAP [16]. More information and data extraction of CSES can be accessed through http://www.leos.ac.cn, accessed on 29 December 2021.

Extract the Disturbance
As the disturbing signals are hidden in the original data, we used the moving average algorithm to extract these disturbing signals. This algorithm adopts the average value in the window as the center point's value, and dynamically adjusts the window size at the extreme point; the size of the window is 50. The disturbing signal can be obtained by subtracting the original data from the smooth data. Then, in order to distinguish the abnormal data, we calculate the root mean square deviation (RMSD) of the window and compare it with the RMSD of all data. Data exceeding the limit will be treated as abnormal data. The RMSD is calculated as follows: where x i is the ith data, x is the average value of all data and n is the number of data, while to calculate the RMSD of all data, n is the number of all data and to calculate the RMSD of windows, n is 50. As the plasma changes drastically in high latitude regions, we do not consider regions with latitudes greater than 50 • . In order to make an in-depth comparison with past data and reflect the trend of Ne and other parameters, according to the characteristics of the revisit orbits specially designed for CSES, the revisit orbits in a period of time before the earthquakes are taken as the reference background field.
As the satellite have some discrete points that exceed the normal data range by one to two orders of magnitude during the measurement process, these discrete points may be measurement errors caused by various reasons and cannot reflect the actual ionospheric conditions. Therefore, to avoid the influence of these discrete points on the results, according to their characteristics, the data points more than six times the average value of the whole track data (from 50 • S to 50 • N) are regarded as discrete points and eliminated. Then, these data are averaged according to the revisit orbits, and the averaged value is used as the background field. Please note that according to the definition in Section 2.1, the data with poor space environment will not be used to avoid affecting the results [17].

Wavelet Analysis
Wavelet transform is a signal processing method developed based on the Fourier transform, and there are already many mature algorithms. It can realize arbitrary "amplification" or "reduction" of partial signal through translation and expansion parameters, so as to obtain the corresponding time-frequency resolution of different frequency components of a partial signal [18]. The disturbance caused by earthquakes contains the time-varying unsteady signal in which we are interested. In order to study these signals, we introduce the wavelet transform. Due to the particularity of satellite measurement, its measurement time corresponds to latitude. Therefore, we convert time into latitude to reflect the spatial information of the disturbing signal.

Trend over Date
In order to reflect the long-term trend, we plot the Ne ( Figure 3) and Te ( Figure 4) data from January 1 to June 15 in chronological order. In each graph, the first subgraph is the original data, the second is the moving average data, the third is the extracted disturbance, and the red line in the graphs is the latitude of Yangbi and Maduo earthquakes. All the data are displayed into two parts: the whole orbits (from 50 • S to 50 • N, Figures 3a and 4a) and the orbits in the seismogenic zone (Figures 3b and 4b). . Seismogenic zone orbits; the seismogenic zone is 21 • N to 49 • N, and we actually use 20 • N to 50 • N. In each figure, the first subfigure represents the original data, the second represents the smooth data, and the third represents the abnormal data. Red lines indicate the latitude of the earthquakes. We take the date of the first earthquake (Yangbi earthquake, 21 May) as the date of these two earthquakes (in red). . Seismogenic zone orbits; the seismogenic zone is 21 • N to 49 • N, and we actually use 20 • N to 50 • N. In each figure, the first subfigure represents the original data, the second represents the smooth data, and the third represents the abnormal data. Red lines indicate the latitude of the earthquakes. We take the date of the first earthquake (Yangbi earthquake, 21 May) as the date of these two earthquakes (in red).
Whether Ne or Te, it can be found that from January to March, the original curves and moving average curves are roughly similar in form. However, since about 40 days before the earthquakes, that is, since 10 April, the form has changed. Ne even began to interweave about 20 days before the earthquakes. This situation changed obviously after the earthquakes; they began to become similar again. In addition, these violent changes are mainly in the seismogenic zone or the magnetic conjugation area of the two earthquakes.
It should be noted that the signal near 45 • N in the above results is caused by CSES when switching the observation mode, and it is not processed separately in this paper [19].

Trend with Space
We further studied the relationship between disturbance and space. The extracted disturbances are plotted on the same graph according to longitude. The result of Te is shown in Figure 5. Obviously, these disturbances are concentrated in the orbits passing through the seismogenic zone and the surrounding orbits, indicating that these disturbances may be correlated with the two earthquakes in space. On May 17 and 19, Te disturbances all appeared on the equatorial side of the seismogenic zone, but there was no obvious disturbance in the seismogenic zone. According to previous studies, the ionospheric anomaly does not completely coincide with the epicenter in space, which may be due to the additional electric field generated in the seismogenic zone. Under the action of this additional electric field, the plasma anomaly deviates from the epicenter [20,21].
After removing the background field of Ne, the extracted disturbance is as shown in Figure 6. Considering that the ionosphere has seasonal changes, we use the data from January to March 2021 as the reference background field. It can be found from the figure that Ne has a relatively obvious increase in the epicenter area. In addition, we treat the orbits that the RMSD of moving windows exceeds the RMSD of all data by more than three times as anomalies; the accumulation of anomalies show in Figure 7. Obviously, these disturbances gradually become more frequent as the earthquake approaches, and their intensity also increases. After the earthquakes, it gradually becomes calm again.

Ion Compositions
According to the research results of Parrot, before a strong earthquake, the ion compositions in the ionosphere will also change [22]. Therefore, in order to corroborate the anomalies of Ne and Te observed by CSES before these earthquakes, we analyzed the ion compositions of these orbits, as shown in Figure 8.
Cyan, green, blue and purple represent the ion density of H + , He + , O + and Ne, respectively, and the two red lines represent the latitude position of the epicenter. Although the data above are not from the same day, it can also be found that in the seismogenic zone, the density of these ion compositions and Ne all began to change drastically, especially between the two earthquakes, in which the magnitude of the change was more drastic. Their change characteristics are consistent, which indicates they may be affected by the same source of disturbance.
Lots of ionosphere structures, such as equatorial plasma bubbles (EPBs), can frequently happen at low and middle latitudes. They can also manifest the obvious ionospheric disturbances similar to those caused by earthquakes. On the basis of Stolle et al. and Luhr et al., Li Bo studied the characteristics of EPBs and provided a method to distinguish them [23][24][25]. When an EPB occurs, its Ne drops sharply and the magnetic field increases in the background. As shown in Figure 9, where the Ne drops sharply at 2 • N and 20-30 • N, the magnetic field has no changes; therefore, this does not fit the characteristics of EPBs.
It should be noted that the HPM of CSES-1 was seriously interfered with while the TBB turns on and off, so the magnetic field near 6-10 • N have been set empty.  The medium scale traveling ionospheric disturbances (MSTIDs) can also generate similar structures in the ionosphere. Saito et al. study the total electron content (TEC) by GPS; the results demonstrate that GPS-TEC can clearly reveal the evolution and propagation process of MSTIDs, and MSTIDs propagate in the southwest direction in the northern hemisphere [26]. Our team performed the same work by GPS-TEC, but did not find similar structures; thus, the anomalies we studied are not related to MSTID. It should be noted that this conclusion was drawn by our colleagues in the study of these two earthquakes. However, these results are needed for their other pending papers. Considering academic ethics, we only use their conclusions with their authorization, but do not show their pictures in this article.

Wavelet Analysis
The above results are the manifestations of earthquake precursor anomalies in time and space. In order to further analyze the characteristics of disturbances from the frequency and amplitude, we use the bump wavelet packets to perform wavelet transformation on the data of several orbits with disturbances, as shown in Figure 10.
It can be observed that there is a significant increase at 0.05 Hz and below, which may be caused by the long-period trend of the data. The frequency of ion compositions and electron density is significantly enhanced around 0.1 Hz, 0.25 Hz and 0.35 Hz, especially the performance of ion compositions around 0.35 Hz, which has a symmetry characteristic in latitude. This symmetrical structure appears in different parameters and on different dates, but it is not similar to the phenomenon of magnetic conjugation. It appears to have intensified around 30 • N, the midpoint of the two earthquakes. Right now, we have no way to determining whether it was the result of a combination of two earthquakes, or just one of them.
These results with similar characteristics on different dates indicate that these disturbances are likely to be from the same source. In order to determine this source, we carefully checked the geomagnetic environment, solar activity, meteorological activity, and even solar and lunar eclipses on these dates, but nothing was found. Therefore, this source may be related to the two earthquakes.

Discussions
In this paper, the integrated analysis is performed in the form of Ne, Te and ion compositions variations in date, space and frequency for Yangbi and Maduo earthquakes. The seismogenic zone of Yangbi earthquake is completely within Maduo's, and the occurrence time of them is only 4 h apart. In such a situation, it is difficult to distinguish whether the possible ionospheric disturbances are from Yangbi or Maduo, or both. Therefore, the two earthquakes are analyzed as the same event. The data in poor space environment is discarded when studying them.
We applied the moving average algorithm to extract disturbance signals from the Ne and Te data, and then arranged these results to reflect the changes of these disturbances over the date. The results of Ne and Te indicate that the disturbance began to appear about 40 days before the earthquakes, and the original form of Ne also changed very drastically about 20 days before the earthquakes, but Te was not so obvious. After the earthquakes, their original form quickly calmed down again.
In addition, our anomaly accumulation curve also shows that the extracted anomalies increased suddenly and rapidly before the earthquakes, and the increase rate of anomalies slowed down again after the earthquakes. These results indicate the fact that these disturbances are temporally related to these two earthquakes. Before the earthquakes, some changes caused these disturbances, and after the earthquakes, the stresses were released, and the environment that produced these changes in the lithosphere was destroyed and these disturbances could no longer be generated.
Further research found that, in terms of spatial distribution, these disturbances are concentrated near the seismogenic zone, and as there are no other possible sources, such as solar activity, solar and lunar eclipses, on the date of disturbance, the EPBs and MSTIDs are also excluded. Therefore, the sources of these disturbances may be related to earthquakes in space. What is curious are the ion compositions and Ne. They always have very similar changing trends. We have used wavelet transform to study the time-frequency characteristics of these ion compositions and Ne. Obviously, they all have strong signals around 30 • N. These signals have similar frequency characteristics, and they all have similar symmetrical structures around 0.35 Hz. These show that their source should be the same, and the location of this source is around 30 • N.
According to the lithosphere atmosphere ionosphere coupling model, a series of complex changes will occur before the earthquake. These changes, on the one hand, caused the atmospheric electric field growth, and finally changed the electromagnetic field in the ionosphere and trapped particle precipitation (ionosphere and magnetosphere branch). On the other hand, it leads to an increase in atmospheric temperature, and finally produces the outgoing long-wave radiation anomalies (thermal branch) [2]. Oyama et al. observed the Te disturbance a few days before the main shock, around the epicenters, through several large earthquakes [27]. Especially for the Maduo earthquake, Qi et al. have already observed the abnormal microwave brightness temperature before it [28]. Our Te results also reflect this disturbance. Combining our results, we suggest that the disturbances observed by CSES are likely to be a manifestation of the above process. However, we have not conducted ground observations, so we cannot verify whether the source of the disturbances are related to radon, methane, etc.
Our investigation based on this multi-parametric approach re-established some facts of the LAIC mechanism. First of all, the electron density, electron temperature, and ion compositions in the ionosphere are easily affected by earthquakes, which indicates the arrival of earthquakes. They demonstrate similar changes in many aspects-for instance, the electron density and temperature began to appear abnormal about 40 days before the earthquakes, and the disturbance increased about 20 days before the earthquakes. In space, they all appear around the epicenter. Further research on disturbances found that these disturbances have similar frequency characteristics. Combining the LAIC model and the others' research results of the Maduo earthquake, we suggest that these earthquakes may affect the ionosphere through two branches (ionosphere and magnetosphere and thermal branch), and be observed.
The deficiency is that although we treat the two earthquakes as one, we still expect new findings. However, we do not have any new findings, or the disturbances in the ionosphere caused by the two earthquakes may be different, but we cannot distinguish them.

Conclusions
This manuscript deals with pre-seismic anomalies during the 2021 Yangbi and Maduo earthquakes, as observed from China seismo-electromagnetic satellite. We applied the moving average algorithm, removing background field and conducted a wavelet analysis on the data of electron density, electron temperature and ion compositions for our study. We discarded the data on the poor spatial environment date to avoid the influence of other factors on the results as much as possible. The following are the main outcomes from this study: The anomalies of electron density and temperature began to appear around 40 days and were enhanced around 20 days before the earthquakes. This suggests that some matter and energy have already begun to affect the ionosphere before the earthquakes.
The disturbances of LAP and PAP sensors are mainly distributed around the seismogenic zone, which indicates that there have been some disturbance sources at the epicenter location before the earthquakes.
Wavelet analysis was performed on the extracted disturbances, and the results demonstrated that they were relatively similar in frequency, especially the ion compositions, which seemed to have a symmetrical structure near 0.35 Hz. These results indicate that the source of the disturbances may be the same. However, we have no way of determining whether it was the result of a combination of two earthquakes, or just one of them.
Our anomalous cumulative curve and other researchers' results about the earthquakes corroborated the disturbances of the ionosphere during the seismogenic stage of the earthquakes.
Finally, based on our results, we suggest that the impact of the earthquakes on the ionosphere may be realized through two major branches viz. (a) ionosphere and magnetosphere branch and (b) thermal branch of the LAIC mechanism.
Author Contributions: Conceptualization, X.Z.; Data curation, X.Z.; investigation, X.Z. and X.D.; methodology, X.D. and X.Z.; supervision, X.Z. All authors have read and agreed to the published version of the manuscript.