Statistical analyses on the ionospheric total electron content related to M ≥ 6.0 earthquakes in China during 1998-2015

To verify the pre-earthquake ionospheric anomaly (PEIA), statistical analyses are implemented on the relationship between the total electron content (TEC) of global ionosphere map (GIM) and 62 M ≥ 6.0 earthquakes in China during 1998 2015. A median-based method together with z test is employed to determine the criteria and/or characteristics of TEC anomalies related to earthquakes. It is found that the GIM TEC significantly decreases at 18:00 22:00 UT (universal time, post-midnight to pre-dawn) 4 5 days before 37 6.0 ≤ M < 6.5 earthquakes, at 01:00 04:00 UT (morning) 3 6 days before 18 6.5 ≤ M < 7.0 earthquakes, and 04:00 10:00 UT (pre-noon to afternoon) 3 5 days; but increases 08:00 12:00 UT (late afternoon to early evening), 18 20 days before 7 M ≥ 7.0 earthquakes. The receiver operating characteristic (ROC) curve is used to compare the TEC anomaly-based method with some competitive alternatives for predicting the earthquakes under study. We found, based on possible TEC anomalies, that the observed PEIAs are significantly earthquake-related. Moreover, the results of regression analyses show that the PEIA strength is associated with the magnitude of earthquakes. Article history: Received 29 December 2017 Revised 11 February 2018 Accepted 11 March 2018


INTRODUCTION
Anomalous changes in the electron density and/or electromagnetic signals in the ionosphere before large earthquakes have been intensively studied (Hayakawa and Fujinawa 1994;Hayakawa 1999;Liu et al. 2000Liu et al. , 2006Liu et al. , 2013aHayakawa and Molchanov 2002;Pulinets and Boyarchuk 2004). For the short-term earthquake prediction or forecast, it needs recognizable and reliable precursors. Liu et al. (2001) pioneer developing the total electron content (TEC) derived from measurements of ground-based GPS receivers and find a significant decrease of the GPS TEC in the afternoon period 1, 3, and 4 days prior to the 21 September 1999 M 7.7 Chi-Chi earthquake. Liu et al. (2004a) further conduct a statistical analysis to examine the GPS TEC and 20 M ≥ 6.0 earthquakes in the Taiwan area from September 1999 to December 2002. Results show that the negative (i.e., de-crease) TEC anomalies appear during the late afternoon and evening period of 18:00 -22:00 LT (local time) 1 -5 days prior to 16 out of the 20 M ≥ 6.0 earthquakes. Moreover, the characteristic of the negative polarity (i.e., the TEC decrease), the appearance local time in the afternoon/evening period, and the 1 -5 lead days of the earthquakes reached by the statistical analyses generally agree with those for the Chi-Chi earthquake. In the statistical results, 16 out of 20 earthquakes being successfully alarmed suggests that the TEC could be useful to detect pre-earthquake ionospheric anomalies (PEIAs) of large earthquakes. In fact, many statistical methods have been applied testing characteristics of PEIAs in the electron density in terms of the ionospheric F2peak plasma frequency foF2 (Liu et al. 2000(Liu et al. , 2004b(Liu et al. , 2006Chen et al. 2004) and of the GPS TEC (Liu et al. 2004a(Liu et al. , b, 2010a to see whether recognizable and reliable precursors exist in Taiwan. All the related statistical studies yield consistent results in the PEIA characteristic that negative E a r l y R e l e a s e Terr. Atmos. Ocean. Sci., Vol. 29, No. 5, 1-14, October 2018 anomalies frequently appear in the afternoon/evening period 1 -5 days before the earthquakes under study. Therefore, the PEIAs of the ionospheric TEC are strongly suggested to be reliable seismo-ionospheric precursors (SIPs) in Taiwan.
Based on Liu et al. (2001Liu et al. ( , 2004a, Liu et al. (2009) first time employ the global ionospheric map (GIM) to find the PEIA of the GPS TEC associated with 35 M ≥ 6.0 earthquakes in China during the 10-year period of 1 May 1998 to 30 April 2008. Anomalies of the earthquakes and the 12 May 2008 M 8.2 (Mw 7.9) Wenchuan earthquake (30. 986°N, 103.364°E) show that the PEIA characteristic in China is the profound decrease of the GIM TEC in the afternoon period 3 -5 days before large earthquakes. Liu et al. (2010b) report that PEIAs of the 26 December 2004 M 9.3 Sumatra-Andaman earthquake meet the characteristic of the statistical results, which the GIM TEC around the epicenter significantly reduces during the afternoon period 1 -5 days before 100 M ≥ 6.0 earthquakes occurring in Indonesia from 1 May 1998 to 31 December 2008. By contrast, Kon et al. (2011) and  statistically examine TEC anomalies associated with M ≥ 6.0 earthquakes in Japan during 1998 -2011, and find the PEIA characteristic to be positive anomalies significantly appearing 1 -5 days before the earthquakes. These statistical results show that the PEIA characteristic (i.e., polarity, appearance local time, duration, lead day, etc.) varies in different locations. Therefore, for any monitoring area, it is essential to find whether the PEIA characteristic reached by statistical analyses is reliable SIPs or not.
Note that  statistically study the GIM TEC associated with 56 M ≥ 6.0 earthquakes in China during 1998 -2012, and find that the PEIA characteristic is the significant decrease of TEC in the afternoon period 2 -9 days before the earthquakes. The characteristic agrees with PEIAs of negative anomalies appearing in the afternoon period 2 -6 days before the Wenchuan earthquake, although Liu et al. (2009 simply count percentages of the earthquakes with positive and negative anomalies, respectively, in their statistical studies. To see if SIPs exist and reliable, Chen et al. (2015) investigate the GIM TEC and 56 M ≥ 6.0 earthquakes in China during 1998 -2013 by means of the z test (Neter et al. 1988) and the ROC (receiver operating characteristic) curve (Swets 1988). In this study, we extend the existing studies, examine 62 M ≥ 6.0 earthquakes reported by CENC (China Earthquake Networks Center, http://www.csndmc. ac.cn/newweb/catalog _ direct _ link.htm) and the GIM TEC at a fixed location (32.5°N, 95°E, the location center of those earthquakes) retreated from CODE (Center for Orbit Determination in Europe, ftp://ftp.unibe.ch/aiub/CODE/) in China during 1998 -2015 (Figs. 1, 2, and Table 1). The statistical significance of the PEIAs associated with the earthquakes is investigated by the z test, and ROC curve. To avoid possible confounded effects, we subdivide the earthquakes into three non-overlapping groups, 6.0 ≤ M < 6.5, 6.5 ≤ M < 7.0, and M ≥ 7.0, and find the associated characteristic of the observed    E a r l y R e l e a s e PEIA in each group. Meanwhile, we randomize the observed anomalous days to verify the significance of the PEIAs. Finally, a logistic regression model (Hosmer and Lemeshow 1989) is applied to find the relationship between earthquake parameters and the associated PEIA strength.

Detection of PEIAs
There are 62 M ≥ 6.0 earthquakes occurring in China during the 18-year period of 1998 -2015 ( Fig. 1 and Table 1). To detect PEIAs that are possibly associated with the earthquakes, we extract and examine the TEC at a fixed location (32.5°N, 95°E), where is the center of these earthquakes, from the CODE GIM with 2.5°-latitude × 5°-longitude spatial resolution. The GIM TEC, routinely published with 2-hr time resolution, is smoothed by interpolating with a cubic spline function (De Boor 1978) to obtain a value every 15 min. Therefore, it has 96 TEC observation points daily. To identify anomalous signals at each time point, we calculate the median (second quartile; 50%), lower (first quartile; 25%) and upper (third quartile; 75%) quartiles of TEC, denoted by M O , LQ and UQ, respectively, based on previous 15-day TECs. Then the lower bound and upper bound are given by where k is a threshold constant. Based on Liu et al. (2009), to have a stringent criterion, we set k = 1.5, the lower bound, O , and hence the probability of a new GPS TEC in the interval (LB, UB) is approximately 65% (Kotz et al. 1983). If an observed TEC value falls outside the associated LB or UB, a negative or positive anomaly is declared. To minimize and/or remove seasonal and solar activity effects as well as to find the relative difference of TEC to its reference value, we compute the standardized TEC Figure 2 depicts the TEC, δTEC, and the earthquakes for the entire 18-year period. It can be seen that the TEC experiences prominent solar activity (cycle) and seasonal effects, while those effects on δTEC are drastically reduced, especially the solar cycle one. Figure  When the TEC locates out of the lower or upper bound, which means the negative or positive anomaly being detected, the anomaly strength is marked, and δTEC is denoted for anomalous decrease or increase accordingly. It can be seen that negative anomalies frequently appear during 00:00 -09:40 UT [06:20 -16:00 LT, local time is universal time (UT) plus 06 hrs and 20 mins at 95°E] 0 -6 days before the Wenchuan earthquake. Note that positive anomalies prominently appear during 08:00 -11:40 UT (14:20 -18:00 LT) on 9 May, day 3 before the earthquake. For a detailed study, Figs. 4a and b illustrate diurnal variations of δTEC with and without denoting its negative and positive anomaly (k = 1.5) 30 days before the Wenchuan earthquake, respectively. It can be seen that δTEC frequently decreases anomalously during 00:00 -09:40 UT, from the early morning to afternoon 06:20 -16:00 LT, 1 -6 days before the earthquake (Fig. 4b).

z Test
To see whether the PEIA exists or not, the GIM TEC and the 62 M ≥ 6.0 earthquakes occurring in Southwest China during the 18-year period are examined (Fig. 1). Figure 5 illustrates the median of δTEC 30 days before three group of 37 6.0 ≤ M < 6.5, 18 6.5 ≤ M < 7.0, and 7 M ≥ 7.0 earthquakes. A z test is then applied to find if δTEC is a statistically significant increase or decrease. Let r be the observed proportion of earthquake-related anomalies and π 0 the background proportion of anomalies in the entire 18year of 6488 (28 March 1998 to 31 December 2015) days. The z value is then given by where n = 62 is the number of earthquakes. If z > 1.96, we claim, at significant level 0.05, that > 0 r r . Note that the z test is conducted for negative and positive anomalies separately.
Since the common characteristic of the three group earthquakes is that negative anomalies appear 3 -6 days before the earthquakes, if defining the negative or positive anomaly day in which more than one third of negative or positive anomalies appearing within the detected anomalous period, it will allow us studying the overall negative or positive anomaly days appearing before and after the earthquakes. Figure 6 illustrates the percentage of the negative and positive anomaly days, as well as their difference 30 days before and after the earthquakes. The percentage is the number of the negative (or positive) anomaly days dividing by the total earthquake number of each group. It can be seen that the negative anomalies prominently appear day 4 -5 before the 6.0 ≤ M < 6.5 earthquakes, day 3 -6 before the 6.5 ≤ M < 7.0 earthquakes, and day 3 -5 before the M ≥ 7.0 earthquakes, which generally agree with significantly negative anomalies in the three zones of the three groups shown in Fig. 5. A detailed study of Figs. 5 and 6 gives 59.5% (22 out of 37) of 6.0 ≤ M < 6.5, 72.2% (13 out of 18) of 6.5 ≤ M < 7.0, and 85.7% (6 out of 7) of M ≥ 7.0 earthquakes being preceded by the PEIAs of negative TEC anomalies (also Table 1) and shows that a larger earthquake has a better chance to experience the PEIA.
The negative anomalies appear during 00:00 -09:40 UT (06:20 -16:00 LT) within 1 -6 days before the Wenchuan earthquake (Fig. 4b), which yields the best agreement in the anomaly appearance time of Zone C (i.e., M ≥ 7.0 earthquakes) among the three negative zones (Fig. 5c). This agreement suggests that an observed PEIA could shed some light on both the magnitude and the lead day of a possible forthcoming earthquake. On the other hand, for the positive TEC anomalies and the threshold k = 1.5, 100% (7 out of 7) of M ≥ 7.0 earthquakes are preceded by the positive PEIA (Table 1, Fig. 5c). Thus, the M ≥ 7.0 earthquakes most likely experience both the negative and positive PEIAs.

ROC Curve
To further verify if the PEIA is a candidate for a reliable precursor of SIP, we treat the PEIAs as alarms for earthquakes in Zones A, B, C, and D (Fig. 5), and construct the ROC curve to evaluate the reliability of the earthquake alarming. Taking Zone B as an example, based on its PEIA characteristic for k = 1.5, when negative anomalies appear E a r l y R e l e a s e E a r l y R e l e a s e more than one third of the period of 01:00 -04:00 UT, we can issue an alarm for an earthquake with the magnitude of 6.5 ≤ M < 7.0 occurring in following 3 -6 days. Note that k = 1.5 is simply an imperial threshold in the previous studies (Liu et al. 2008(Liu et al. , 2009(Liu et al. , 2010a(Liu et al. , b, 2011(Liu et al. , 2013a. In fact, for small k values, many alarming days are issued, and hence more false alarms are obtained. However, for large k values, limited alarming days are issued and both the false alarm and successful rates are drastically reduced. Hence, to test the preference of the PEIAs, we consider constructing the ROC curve for various k values. For each k value, we examine four different conditions, an alarm day being followed by earthquakes or no earthquake, and a non-alarm day being followed by earthquakes or no earthquake within a certain lead day period (Table 2). Let TP(k) and FN(k) stand for numbers of earthquake days with and without being led by alarm days, respectively. Moreover, denote FP(k) and TN(k) to be numbers of non-earthquake days with and without being alarmed, respectively. Then, we have a 2 × 2 contingency table (Table 2) and the true positive rate TPR(k) and false positive rate FPR(k) as given by where TPR(k) is the probability that an earthquake is successfully alarmed, and FPR(k) is the probability to make a false alarm. Hence, the ROC curve with FPR(k) as the xaxis and TPR(k) as the y-axis can be constructed. The area under the ROC curve (AUC) is further used for assessing the effectiveness of the PEIAs (Bradley 1997). Note that, when TPR(k) = FPR(k) for all k, an equal chance to alarm earthquake day and non-earthquake day, we have AUC = 0.5. Therefore, a reliable precursor should have AUC > 0.5. Note that the value of k varies from 0 to 10 by increasing 0.1, and therefore there are 101 k values under investigation. Figure 7 shows that the ROC curves (red curves) of the observed PEIAs, a TPR is generally greater than the associated FPR, and the AUCs are greater than 0.5 in the three groups of earthquakes.
To find the significance of the observed ROC curve, we perform a simulation-based statistical test for the null hypothesis H 0 (SIPs are randomly observed) versus the alternative hypothesis H A (SIPs are practically observed). We investigate the efficiency of any possible sequence of PEIAs that raising an earthquake alarm for those earthquakes under study via the criteria determined by z test. Again, taking Zone B as an example, we find that with k = 2.1, there are 1014 PEIA days, more than one third of negative in the period of 01:00 -04:00 UT being detected, which results in an appearance distribution of 1013 (= 1014 -1) inter-PEIA times during 1998 -2015. Here, for each k, we obtain associated distribution of PEIA appearances for earthquake alarming. We then simulate 1000 sequences of alarming days under the distribution for each of 101 k values and construct the E a r l y R e l e a s e related ROC curves. Figure 7 reveals 1000 simulated ROC curves (gray curves) with the associated TPR(k) * , FPR(k) * , and AUC * . The simulation results allow us further compute the p value, which is the proportion of the simulated AUC * larger than the observed AUC. Therefore, small p values (< 0.05) lead to the rejection of the null hypothesis H 0 . Table 3 shows that the AUC of the three earthquake groups are all greater than 0.5 and the resulted p values are all zero. Therefore, taking the AUC as an overall performance, the statistical analysis of ROC curve confirms that the PEIA is a reliable precursor. Meanwhile, Table 3 depicts that the ob-served AUC is proportional to the earthquake magnitude, which indicates that a larger earthquake has a better chance to be led by significant SIPs.
Note that the blue curve in Fig. 7 is the simulated 95% upper confidence bound for the ROC curve and most part of the observed ROC curve with k either not too small or not too large is above the blue curve. Therefore, it is essential to find an optimal k value for identifying a useful PEIA as a reliable earthquake precursor of SIPs. In fact, the observed AUC curve above the blue line locates the suitable range of k value and therefore, k = 1.1 -1.6 for 6.0 ≤ M < 6.5

Earthquake
No earthquake

Alarm TP(k) FP(k)
No alarm FN(k) TN(k)  Fig. 7. ROC curves for alarming earthquakes based on precursory information from the four zones indicated in Fig. 5. The red, gray, and blue curves denote the ROC curves of the observations, 1000 simulations, and the 95% line of simulations respectively. The red arrow denotes the best point yielding the maximum R score (= TPR -FPR), which is called the Youden index (Youden 1950).
Please describe (a) to (d)! E a r l y R e l e a s e earthquakes (Fig. 6a), k = 1.1 -2.2 for 6.0 ≤ M < 7.0 earthquakes (Fig. 7b), and k = 1 -2.2 for M ≥ 7.0 earthquakes (Fig. 7c). However, to search for an optimal k value, we further compute the Youden index (Youden 1950), which is the maximum of the k-dependent R scores (Shi et al. 2001;Chen et al. 2004) with R as the difference between TPR(k) and FPR(k), or TPR(k)-FPR(k). For three negative polarity groups of earthquakes, the Youden index of 0.2139, 0.3412, and 0.3725 give the optimal value of k = 1.5, k = 2.1, and k = 1.4 respectively. We find that the optimal k = 1.5 of the 6.0 ≤ M < 6.5 earthquakes is coincident with the empirical value of k = 1.5 (Liu et al. 2009) with TPR(1.5) = 0.5946 and FPR(1.5) = 0.3807. Once again, taking the 6.5 ≤ M < 7.0 earthquakes as an example, the optimal k = 2.1 produces TPR(2.1) = 0.7222 and FPR(2.1) = 0.3810, while k = 1.5 yields TPR(1.5) = 0.8333 and FPR(1.5) = 0.5221. Obviously, both TPR and FPR, as a penalty, are increasing with k, but the R score may decrease. For M ≥ 7.0 earthquakes, the optimal k = 1.4 gives TPR(1.4) = 0.8751 and FPR(1.4) = 0.4100. It is worthy to mention that although the three groups of earthquakes have different optimal k values, the FPRs are somewhat similar but the TPRs are proportional to the earthquake magnitude. Figure 5c reveals that positive anomalies significantly appear in the afternoon period 14:20 -18:20 LT (08:00 -12:00 UT) 18 -20 days before the 7 M ≥ 7.0 earthquakes. Figure 7d shows that the optimal k = 2.0 produces TPR(2.0) = 1 and FPR(2.0) = 0.4155, and the AUC = 0.8026. This confirms that both the negative and positive PEIAs appear before the 7 M ≥ 7.0 earthquakes. Table 1 and Fig. 6 depict that the greater earthquakes have a higher chance to be preceded by PEIAs, while Fig. 5 reveals that a greater earthquake tend to be led by a larger PEIA of δTEC. Therefore, we further examine the relationship between occurrences of the PEIAs and parameters of the related earthquakes. The logistic regression model (Hosmer and Lemeshow 1989) is used to investigate how the occurrence of the PEIA is related to the associated magnitude or distance between the epicenter of each earthquake and the TEC monitoring location. On the basis of the 62 earthquake days, 42 with and 20 without PEIA, the fitted logistic regression models for the logarithm of the odds of PEIA occurrence against the corresponding magnitude and distance, respectively, are obtained as

Logistic Regression
and where P is the probability that an earthquake experiences the PEIA, which is the number of earthquakes with the PEIA divided by the total earthquake number, and M is the magnitude of corresponding earthquake and D is the associated distance in 100 km. Results show the larger the earthquake, the better chance for the earthquake to be recognized by the PEIA (Fig. 8a). Moreover, the odds are greater than 2.0 only for earthquakes with magnitude at least 6.4. Here, 73.3%, 22 out of 30, M ≥ 6.4 earthquakes experience the PEIA. However, the slope in Eq. (7) is not significantly different from zero (also see Fig. 8b). This result seems to indicate that the more nearby earthquakes may not have a better chance to be recognized. For example, we observe from Fig. 1 that there are some more nearby small earthquakes without PEIA but some large earthquakes far away with PEIAs. In other words, the magnitude and distance of earthquakes under study are confounded. Nevertheless, these findings demonstrate that the PEIA is related to the energy with a rather complex geometry released from the associated earthquake (Båth 1966;Lay and Wallace 1995). Figures 5 and 8a show that a larger earthquake has a better chance to be preceded by PEIAs (also see Table 3), and it would be interesting to find whether the PEIA strength, which includes both the duration and the δTEC value, is proportional to the magnitude of possible forthcoming earthquakes. Because only the M ≥ 7.0 earthquakes can issue the positive anomalies, we simply deal with the negative anomalies. For each alarmed earthquake, we define the PEIA strength, S, as the sum of the absolute value of δTEC over the PEIA time period and days. Thus, the strength S E a r l y R e l e a s e simultaneously takes δTEC and the persistency of the PEIA into account. We fit the 42 alarmed earthquakes' magnitude M versus the associated strength S with the linear regression model (Devore 1988), which can be written as (the dashed red line in Fig. 9), .
. M S 0 0420 6 2489 = + However, there are some earthquakes with the same magnitude but different strengths and therefore, a median strength S K is calculated. In total, 13 medians are obtained. The linear regression model can be expressed as (the solid black line in Fig. 9), Hereafter, at (32.5°N, 95°E), when δTEC negatively polarity anomalies meet the PEIA characteristic, based on Eq. (9), the strength of PEIA could shed some lights on the magnitude for possible forthcoming earthquakes in China. Liu et al. (2009) observed GIM TEC above the forthcoming epicenter anomalously decreases in the afternoon period of day 4 -6 and in the late evening period of day 3 before the earthquake, but enhances in the afternoon of day 3 before the 12 May 2008 Mw 7.9 Wenchuan earthquake (30.986°N, 103.364°E). Their result generally agrees with the negative anomalies in δTEC frequently appear over the monitoring point (32.5°N, 95°E) during the morning-afternoon period of 06:20 -16:00 LT 0 -6 days before the Wenchuan earthquake ( Figs. 2 and 3). This agreement and the distance between the epicenter and the monitoring point indicate that the SIP area of the GIM TEC of the Wenchuan earthquake is rather huge. Liu et al. (2009) also examined the GIM TEC over the epicenter associated with 35 M ≥ 6.0 earthquakes in China during the 10-year period of 1 May 1998 to 30 April 2008. They find that the GIM TEC above the epicenter often profoundly decreases in the afternoon period on day 3 -5 before 17 M ≥ 6.3 earthquakes, which generally agree with that the GIM TEC frequently decreases in the morning period day 3 -6 before the 6.5 ≤ M < 7.0 earthquakes, and in the pre-noon to afternoon period day 3 -5 before the M ≥ 7.0.earthquakes shown in Figs. 4 and 5. The slight discrepancies in the above two studies suggest that the appearance time and the lead day are a function of the earthquake magnitude. Figure 8a shows that the PEIAs associated M ≥ 6.4 earthquakes are statistical significance, while Liu et al. (2009) reports that the GIM TEC over the epicenter significantly decreases before M ≥ 6.3 earthquakes. Liu et al. (2009) reporting PEIAs associate with M ≥ 6.3 earthquakes is because they examine the GIM TEC over each epicenter and do not take the distance into account, which certainly is not realistic. Therefore, PEIAs associated M ≥ 6.4 earthquakes are statistical significance in China. Chen et al. (2015) conduct rigorous statistical analyses of using z test and ROC curve on the GIM TEC over the fixed location (35°N, 90°E) associated with the 56 M ≥ 6.0, 27 M ≥ 6.5, and 6 M ≥ 7.0 earthquakes. They find that the PEIAs of the three groups appearing in the morning period 06:00 -10:00 LT (00:00 -04:00 UT) 1 -6 days before the earthquakes are significant, which are reliable SIPs. By contrast, we find Zone A (18:00 -22:00 UT 4 -5 days before E a r l y R e l e a s e 37 6.0 ≤ M < 6.5 earthquakes), Zone B (01:00 -04:00 UT 3 -6 days before 18 6.5 ≤ M < 7.0 earthquakes), and Zone C (04:00 -10:00 UT 3 -5 days before M ≥ 7.0 earthquakes) being significant. The main discrepancy between the two studies is that PEIAs reported by Chen et al. (2015) appear at, their Zone B, 00:00 -04:00 UT 1 -6 days before the 6 M ≥ 7.0 earthquakes, while this study finds PEIAs occurring at, Zone C, 04:00 -10:00 UT 3 -5 days before the 7 M ≥ 7.0 earthquakes. Note that to conduct tests statistical significance for various PEIA zones, Chen et al. (2015) simply based on the PEIA characteristic of the 27 M ≥ 6.5 earthquakes set the same single zone for their three earthquake groups. Their Fig. 3c shows for z test that the PEIAs appearing in afternoon period of 12:00 -16:00 LT (06:00 -10:00 UT) 3 -5 days before the 6 M ≥ 7.0 earthquakes are significant, which in fact well agrees with time period of 04:00 -10:00 UT 3 -5 days before the 7 M ≥ 7.0 earthquakes reached by this study. Nevertheless, owing to Chen et al. (2015) overlapping the classifications of M ≥ 6.0, M ≥ 6.5, and M ≥ 7.0 earthquakes, their PEIA characteristic is contaminated by the zone with the greater earthquakes, and in turn is difficult to find the relationship between the PEIA strength and the earthquake magnitude. Since different magnitudes have different PEIA characteristics, our study shows the appearance times and lead days could provide additional information for possible forthcoming earthquakes. In addition, Chen et al. (2015) randomly generate earthquake days with the observed ones during the study period, and then calculate the ROC curve and the corresponding AUC for each simulation. Based on 500 simulations, the p value of the significance AUC was obtained. Their results show that almost all the AUCs are greater than 0.5, and the p values are less than 0.05 when negative PEIAs appear in the morning 06:00 -10:00 LT 1 -6 day before their 3-range earthquakes, which somewhat agrees with our Zone A (00:20 -04:20 LT, 4 -5 days before the 37 6.0 ≤ M < 6.5 earthquakes). However, they reported that the p values for negative PEIAs appearing in the afternoon 13:00 -17:00 LT, 1 -6 days before their 3-range earthquakes are greater than 0.05 which are insignificant, and all their p values lie between 0.004 and 0.568. In contrast, our p values all equal to zero for the 3 negative (also 1 positive) polarity zone tests. The discrepancy results from the earthquake groups being overlapped, the local time and lead day being not properly assigned, and the simulations being conducted by randomizing earthquake days in Chen et al. (2015). Nevertheless, after carefully grouping the earthquakes, specifying the appearance time and lead day based the characteristic by z test, and carrying out correct simulations, the ROC results show that the SIPs of the ionospheric TEC are reliable and useful. Liu et al. (2006) on the basis of 93 earthquake days in Taiwan during 1994 -1999, fit the logistic regression models (Hosmer and Lemeshow 1989), and report that the logodds is a linear function of the earthquake magnitude with a positive slope of 1.5 but inversely proportional to the square of the distance. Equations (6) and (7) show that the logodds is a linear function of the earthquake magnitude with a positive slope of 0.67, and however that of the distance is not significantly different from zero (Fig. 8b). No clear relationship can be found in the distance might result from E a r l y R e l e a s e that we set the monitoring location at the center of the study earthquakes (32.5°N, 95°E), which enhances the chance of detecting PEIAs but provides no the distance information because the monitoring center is nearby the epicenter (Liu et al. 2013). Nevertheless, it is important to choose the monitoring location to enhance the PEIA detectability, where usually is the center of previous earthquakes. All positive slopes of the two linear regression models of the strength and the magnitude of the 42 alarmed earthquakes show that the larger PEIA strength follows by the greater earthquake. Therefore, the strength of PEIA could be useful to estimate the magnitude of possible forthcoming earthquakes.

DISCUSSION
Meanwhile, 7 and 6 out of the 7 M ≥ 7.0 earthquakes simultaneously are leaded by the positive and negative PE-IAs, respectively. In addition to the characteristic reached by the statistical analyses, the dual PEIA features might be useful in turn to suggest possible forthcoming M ≥ 7.0 earthquakes in China.

CONCLUSION
In conclusion, z test shows that to identify PEIAs at a certain area, the characteristic of the polarity, local time, duration, lead day etc. has to be statistically constructed first. Since different magnitude earthquakes inhabit different characteristics, in turn the observed polarity, local time, duration, lead day, etc., could leak some information of the magnitude of forthcoming earthquakes. The ROC curve confirms that the SIPs of the GIM TEC in China are statistically significant and reliable. The logistic regression proves that the larger earthquake has a better chance to be preceded by PEIAs, while the regression shows that the PEIA strength is proportional to the associated upcoming earthquake magnitude, which implies that larger earthquakes have the greater preparation energy, and in turn release a stronger PEIA strength. It also means larger amplitude of PEIA related to larger earthquake is more reliable. The SIPs in China are negative PEIAs appearing post-midnight to pre-dawn (00:20 -04:20 LT) 4 -5 days before 6.0 ≤ M < 6.5 earthquakes. The SIPs in China are negative PEIAs appearing post-midnight to pre-dawn (00:20 -04:20 LT) 4 -5 days before 6.0 ≤ M < 6.5 earthquakes), morning (07:20 -11:20 LT) 3 -6 days before 6.5 ≤ M < 7.0 earthquakes, and prenoon to afternoon (10:20 -16:20 LT) 3 -5 days before M ≥ 7.0 earthquakes. Positive PEIAs also significantly appear in the afternoon (14:20 -18:20 LT) 18 -20 days (Zone D) before M ≥ 7.0 earthquakes.