West Pacific Earthquake Forecasting Using NOAA Electron Bursts With Independent L-Shells and Ground-Based Magnetic Correlations

Recent advances in statistical correlations between strong earthquakes and several non-seismic phenomena have opened the possibility of formulating warnings within days or even hours. The retrieved correlations have been discovered for those ionospheric physical observations which lasted a long time and realized using the same instruments, including multi-satellite recordings. One of those regarded the electron burst phenomena detected by NOAA, for which the conditional probability of a seismic event was calculated. Then an earthquake probability greater than its frequency was assigned when a satellite realized such a phenomenological observation. This approach refers to the correlations obtained between high-energy electrons detected using the NOAA POES and strong Indonesian and Philippine earthquakes. It is reformulated here to realize a test of earthquake forecasting. The fundamental step is obtained by using a unique electron L-shell interval of 1.21 ≤ L ≤ 1.31 , which decouples the electron parameters from the earthquake parameters. Then, the optimized correlation was recalculated to be 1.5–3.5 h early, between electron bursts and an increased number of seismic events with M ≥ 6 , therein improving the significance too. Moreover, this methodology is reconnected to the frequency theory, and to Molchan’s error diagram, by the probability gain, where a comparison among the significances of various methods is given. The previously proposed physical link between the crust and the ionosphere through magnetic interaction, presumably operating 4–6 h before strong earthquakes, is examined quantitatively on the basis of recent magnetic pulse measurements. Consequently, the probability gain of earthquake forecasting is hypothetically calculated for both the dependent measurements of electron bursts using NOAA satellites and possible ground-based magnetic pulse detection. This method of combining probability gains for earthquake forecasting is general enough that it can be applied to any pair of observables from space and the ground.


INTRODUCTION
Different phenomena, possibly connected with seismic activity, have been reported in recent years by many authors researching anomalies, both geoscientific (Rikitake, 1976;Rikitake, 1987) and macroscopic (Rikitake, 2003). Their research has reported instrumentally repeated observations occurring with strong earthquakes (EQs), and this has permitted them to identify important statistical behavior (Rikitake, 2003). However, in the absence of a recording network, almost all the results are dependent on individual properties of recording, and this renders it rather difficult to obtain an estimation on reasonable statistics (Molchanov and Hayakawa, 2008). The limited number of observatories on the ground and their punctual observations, even when operative (Console, 2001), reduce the number of considered strong EQs, making it too small to calculate a statistical correlation over several decades. Only when moderate magnitude EQs are considered, a statistical correlation is currently calculated for ULF geomagnetic fluctuations at ground stations (Schekotov et al., 2006;Hattori et al., 2013;Han et al., 2014;Han et al., 2017). In other studies, Pc1 anticipated EQs by 6-7 days (Bortnik et al., 2008), VLF noise by 2 days (Oike and Yamada, 1994), lightning activities 17-19 days before EQs (Liu et al., 2015), and geoelectric fields with lead times from days to weeks . A review of several correlation increases corresponding to 3 days between ELF Q-bursts and the Kanchakta EQs has been reported with the possible associated physical models . A method to predict the time, epicenter, and magnitude of such events has been suggested  based on the works cited above.
Observations made by low-orbit satellites are able to monitor large portions of the ground in a few hours, allowing the monitoring of the area affected by each seismic event (Barnhart et al., 2019), and to consider all, or a large portion, of strong EQs. Satellite detection techniques and communication developments are made using electromagnetic instruments and have therefore been immediately used to monitor electromagnetic fields in the ionosphere. Electromagnetic fields measured in the low Earth orbits were associated with strong EQ occurrence for the first time in the 1980s, with regard to electric and magnetic intensity in the range of 1 Hz-10 kHz, when satellites arrived close to EQ epicentres (Larkina et al., 1983;Parrot and Lefeuvre, 1985;Larkina et al., 1989;Parrot and Mogilevsky, 1989;Mikhaylova et al., 1991;Serebryakova et al., 1992). A space-borne system for short-term EQ warning has been suggested (Pulinets, 1998a;Parrot, 2002;Pulinets, 2006), and as ionospheric perturbations measured using satellites are not only due to EQs and are not found for all EQs, a statistical analysis of a possible influence of the seismic activity on the ionosphere is preferred (Parrot, 2011). Statistical results were obtained in 1993 using the Intercosmos-24 satellite (Molchanov, 1993), where the probability of charged particle burst observations was from 6 to 24 h before the event increased by 50%, and the DE-2 satellite (Henderson et al., 1993), where no significant differences occurred between EQ orbits and control orbits. Moreover, the average wave intensity received on board the AUREOL-3 satellite (Parrot, 1994) increased with seismic activity, resulting in an extension in the latitude direction but not in the longitude, with respect to EQ epicentres. The ISIS 1 and 2 satellites have been used to identify the spectra of electromagnetic radiation of seismic events under control data (Rodger et al., 1996), showing no significant evidence for differences between data sets of the EQs and control orbits. A statistical study of intensity for VLF electromagnetic waves has been realized in the vicinity of EQ epicentres using the micro-satellite DEMETER (Nemec et al., 2008). It has evidenced a significant decrease in the measured wave intensity, 0-4 h before strong EQs. A confirmation of this result, on a longer set of data, was obtained (Nemec et al., 2009), suggesting that a significant decrease is occurring for larger EQs, is stronger for shallower EQs, and does not seem to depend on whether the EQ occurs below an ocean or not.
Anomalies appearing in electron densities of the ionospheric F-region a few days before strong EQs were observed (Pulinets 1998b;Liu et al., 2000;Pulinets and Boyarchuk, 2004;Liu et al., 2006). These anomalies concern the electron densities recorded using local ionosondes, where the critical frequency of the F2peak, foF2, significantly decreased days before several EQs. Moreover, decreasing electron densities, days before strong EQs in Taiwan, had been compared with the total electron content (TEC) calculated using ground-based GPS receivers and satellite transmitters (Liu et al., 2004). Anomalous TEC signals were observed in Southern California, but no statistically significant correlations regarding time and space between these TEC anomalies and the occurrence of seismic events resulted (Thomas et al., 2007). On the other hand, positive results for the correlation of EQ from 2 to 5 days after TEC fluctuations have been obtained (Li and Parrot, 2013). Study results from Taiwan support the result that the equatorial ionization anomaly crest significantly moves equatorward 1-5 days before strong EQs (Liu et al., 2010). A statistical analysis carried out on TEC data from the global ionosphere map evidenced that the largest occurrence rates of anomalies were for those EQs with larger magnitudes and lower depths 1-5 days before the EQs (Liu et al., 2011;Zhu et al., 2014). This was confirmed by Japanese (Kon et al., 2011) and Chinese (Ke et al., 2016) studies. Vertical TEC positive and negative anomalies aligned parallel with the local geomagnetic field were repeatedly observed 20-40 min before three M > 8 Chilean EQs (He and Heki, 2016). Concentrations of electron density and magnetic anomalies for more than two months to some days before the EQ occurrences have been reported worldwide (De Santis et al., 2019). However, a 14-year analysis of data did not reveal any statistically significant changes prior to EQs when considering all of the 1,279 EQs together (Thomas et al., 2017). Using GPS TEC measurements, a statistical analysis and comparison of the temporal and spatial distributions for the pre-EQ ionospheric anomalies before the 1, 339 M ≥ 6.0 EQs, which occurred globally between January 2003 and December 2014, did not provide reliable evidence of pre-EQ changes on the global ionospheric map of TEC data (Zhu et al., 2018).
Ionospheric perturbations with seismic activity have included wave paths of VLF and LF transmitters (Molchanov and Hayakawa, 1998;Biagi et al., 2001). The VLF and LF amplitude observations are connected with the entire path covered by waves, even if they are obtained in a punctual station, thus representing an intermediate type of observation between the purely punctual and the completely diffused realized using satellites. Seismoionospheric effects on long sub-ionospheric paths have been investigated in amplitude variations of signals and have used the VLF terminator time method (Clilverd et al., 1999), indicating that the occurrence rate of successful EQ predictions using it cannot be distinguished with respect to a random one. Statistical results obtained by the superimposed epoch analysis in Japan (Maekawa et al., 2006) yielded that the ionosphere was definitively disturbed in terms of both amplitude and dispersion. For an integrated energy, released within the interested area for the LF wave path, the amplitude is depleted and the dispersion is very much enhanced for about one week to a few days before the EQ. A statistical correlation between EQs and VLF/LF signals over 10 years or so, obtained by means of the Japanese VLF/LF network, revealed perturbations 3-6 days prior to wave paths (Hayakawa et al., 2010).
A correlation analysis between earthquakes and atmospheric temperature variations over several months observed using a portable meteostation obtained a time anticipation of about one day (Molchanov et al., 2003). Consequently, thermal infrared anomalies have also been observed with strong EQs from space. A complete review has reported the main contributions and results achieved over 30 years (Tramutoli et al., 2015). Molchan's error diagram analysis computed for different classes of magnitude and significant sequences of thermal infrared anomalies has suggested a prognostic probability gain when compared to random guess results, both for strong EQs in Greece (Eleftheriou et al., 2016) and in the Sichuan area (Zhang and Meng, 2019).
Sudden variations in high-energy charged particle fluxes near the South Atlantic Anomaly have also been associated with seismic activity (Voronov et al., 1989). In fact, numerous experiments followed the discovery of the Van Allen Belts (Van Allen, 1959) to determine safe conditions for near-Earth space exploration. Further detection of charged particle flux variations associated with strong EQs was obtained using the Intercosmos-24 satellite (Galperin et al., 1992;Boskova et al., 1994), resulting in precipitating particles which escaped the trapped conditions of the geomagnetic field. High-energy precipitating particle fluxes have been statistically analyzed in relation to seismic activity in various near-Earth space experiments such as the MIR orbital station and the METEOR-3, GAMMA, and SAMPEX satellites, which have shown particle bursts 2-5 h before EQs (Aleksandrin et al., 2003). A reanalysis of the more recent and extended SAMPEX database has also shown a 3-4 h correlation with precipitating high-energy electrons anticipating strong EQs (Sgrigna et al., 2005). The NOAA-15 satellite particle database, which has been collecting data since 1998, has been systematically studied (Fidani et al., 2010). Sudden variations in high-energy charged particles have been connected with strong EQs in periods of weak solar activity (Fidani, 2015). This statistical correlation analysis evidenced that exceptional increases in electron fluxes occurring 2-3 h prior to the largest quakes had struck the Indonesian and Philippine investigated areas between 1998 and 2014. The correlated EQs occurred at a depth less than 200 km, independent of sea or land. Precipitating particles have been detected far from EQ geographical positions, so the possible disturbances above the EQ epicenters due to particle drift have been estimated to be in the range of 4-6.5 h before strong seismic events (Fidani, 2018).
Definitely, several studies' results have suggested that ionospheric phenomena appear to statistically precede strong earthquakes by up to a week, and some studies even propose longer anticipation times. However, the demonstration of the physical link between the two phenomena is essential to affirm that one of the two is more than a candidate precursor for the other. Moreover, short-term EQ precursors are thought to precede by 1-2 days to several weeks, and near-seismic precursors are thought to precede by several hours to 1-2 days (Molchanov and Hayakawa, 2008). Therefore, electromagnetic fluctuations detected using satellites both in the ELF and VLF bands, together with particle precipitation, may belong to the class of near-seismic candidate precursors, whereas TEC, together with both VLF and LF path amplitude depletion, and geoelectric ULF fields may belong to the class of short-term candidate precursors. As for ground observations, VLF noise belongs to the class of near-seismic candidate precursors, and Pc1 pulsations and lightning activity may belong to the class of short-term candidate precursors, while ELF Q-bursts, Schumann resonances, and ULF magnetic depressions may belong to both classes. The recent results on EQ observations, carried out using low Earth orbit satellites, can be found in a publication by Ouzounov et al. (2018).
Given the possibility reported above to observe physical phenomena that recurrently, and with statistical significance, anticipate strong EQs, a verification of EQ forecasting has been formulated on the basis of schemes already used, starting from the statistical study of seismicity (Console, 2001). In particular, the previously calculated correlations between NOAA electron precipitations and EQs are examined in the fundamental steps, showing what remains ambiguous for electron identification, in Reviewing NOAA Electrons' Statistical Correlations. A reproduction of the statistical correlation between completely decoupled NOAA electron precipitations and EQs, so as to define the precursor phenomenon unambiguously (Wyss, 1997), is obtained in Unambiguous NOAA Electrons' Statistical Correlation. NOAA correlations are reviewed using classic statistical frequency techniques, and their statistical significance is calculated using recent methods of error diagrams in Forecasting Methodologies and Evaluating Significance, where a complete equivalence among these approaches is demonstrated. The objective is to define a methodology to introduce one or more statistically verified precursors in EQ forecasting using conditional probabilities. Prediction Model is devoted to building a prediction scenario following the work of Console (2001). A discussion of any possible relevance to the improvement of the probability gain derived from dependent precursors used together, observed both on the Earth's surface and from space, is presented in Dependent Observables. Finally, the electron precipitations' possible dependence on magnetic pulses is shown using a physical model in Different Precursors Combined, and a hypothetical experiment demonstrates the probability gain due to these two dependent observations. The conclusion is reported in Conclusion with a sequence of steps that combine interdependent observations for EQ forecasting.

MATERIALS AND METHODS
To systematically test methodologies of precursors, for which a statistical analysis of past cases is feasible, the conditional probability of occurrence (Aki, 1981) is a desirable parameter that should supplement the usual time-location-magnitude parameters of the prediction (Console, 2001). Several physical observations on the Earth's surface and from space have been successfully correlated with EQs. However, a conditional probability has only been obtained from the analysis of the NOAA satellites' precipitating electrons (Fidani, 2018;Fidani, 2019;Fidani, 2020), so it would be the most suitable methodology to be tested.

Reviewing NOAA Electrons' Statistical Correlations
NOAA polar satellites use particle detectors which monitor fluxes of protons and electrons in polar orbits at altitudes between 807 and 854 km (Davis, 2007). The particle detectors (Space Environment Monitor SEM-2) consist of the total energy detector (TED) and the medium energy proton and electron detector (MEPED). The MEPED is composed of eight solid-state detectors measuring proton and electron fluxes from 30 keV to 200 MeV (Evans and Greer, 2004) which include the radiation belt populations, energetic solar particle events, and the lowenergy portion of the galactic cosmic ray population. Data can be downloaded at the link http://www.ngdc.noaa.gov/stp/satellite/ poes/dataaccess.html. As all of the sets of orbital parameters are provided every 8 s, this value was chosen as the basic time step for our study (Fidani and Battiston, 2008). Consequently, all other variables were defined with respect to the 8-s step. Thus, 8-s averages of the counting rates (CRs), latitude, longitude, MEPED, and omnidirectional data were calculated. Unreliable CRs with negative values were labeled and excluded from the analysis.
To systematically test the methodology proposed for NOAA data, a quantitative and rigorous definition of the concerned precursor (Console, 2001) was established. The daily averages of particle CRs exiting the entrapment in the geomagnetic field (precipitations) were calculated, and then the condition for which a CR fluctuation was not likely, due to possible statistical fluctuations, was set. This calculation was formulated with a probability larger than 99% (Fidani et al., 2010). The sudden increase in particle flux that satisfies this condition was named particle burst (Sgrigna et al., 2005), and for electrons, the same condition was named electron burst (EB). According to a previous work (Aleksandrin et al., 2003), the daily averages of CRs were calculated in the invariant coordinate space. Together with the L-shell and the pitch angle, it was necessary to take into account the CR amplitudes and their variations versus geomagnetic coordinates, since the spatial gradient of particle fluxes near the South Atlantic Anomaly (SAA) is too large (Fidani and Battiston, 2008). CR distributions inside invariant areas are compatible with a Poisson distribution. Being so, an amplitude threshold was introduced for the CRs to define the conditions for which a CR is a non-Poissonian fluctuation with 99% probability. Furthermore, NOAA satellites measure variations of ionospheric parameters not only due to EQs; indeed, they are principally due to solar activity (Sgrigna et al., 2005;Parrot, 2011). To reduce the effects of solar activity, both low values in Dst variations (http:// wdc.kugi.kyoto-u.ac.jp/dst_final/index.html) and geomagnetic Ap indexes (https://www.ngdc.noaa.gov/geomag/data.shtml) were chosen to exclude CR data corresponding to the Sun's influence (Fidani, 2015).
To conclude the quantification of the precursor (Console, 2001), the L-shell invariant parameter was considered to define the magnetic line where a physical interaction, whatever it is, can connect the seismic and ionospheric activities. Following the works by Aleksandrin et al. (2003) and Sgrigna et al. (2005), particle bursts were considered only when their L-shell values referred to M ≥ 6 EQ epicenter projections that occurred on the magnetic lines. This is equivalent to imagining that the physical interaction can occur in the same region near the vertical. The correlation was calculated by filling a histogram with the time differences Δt between the EQs and EBs, which was indicated by {EQ;EB} (EQ × EB) following the work of Fidani (2015). This approach was performed by considering only EB on magnetic lines identified as projecting EQ coordinates at different altitudes with respect to EQ epicenters, from −600 km up to 3,200 km in increments of 100 km. A correlation peak at Δt 2−3 h started to be significant only for 30-100-keV EBs when considering magnetic line altitudes above 1,400 km and was maximized for 2,200 km (Fidani, 2015) (see Figure 1). Correlations were maximized by using EQs with magnitudes M ≥ 6 downloaded at the link https://earthquake.usgs.gov/earthquakes/search/, located in both the Indonesian and the Philippine regions, having 90°-150°longitudes, with few events in South America. The EBs were detected high off the shore of the United States and the west coasts of South America, at longitudes between 200°and 280°. These different electron positions were associated in a causal way, due to the fact that electrons drift eastward and the EQ positions were located west of the EB detection positions (Fidani, 2018;Fidani, 2020). Being so, if the disturbances which caused electron precipitations from inner radiation belts occurred above the EQ epicenters in the ionosphere, they most likely anticipated the EQ times by 4-6.5 h.
However, to consider EB L-shells around the L-shells corresponding to the EQ epicenter projected at several altitudes constitutes an ambiguity in defining the phenomena that preceded EQs. In fact, the L-shell parameter to choose EB involves both EB and EQ events and each EQ event with a multiplicity of altitude projections. So, the request for the L-shell parameter, which is only related to EBs, is lost and EB cannot be chosen unambiguously. Furthermore, if the condition of L-shell similarity between considered EB detection and EQ projections is not satisfied, the correlation cannot be found. In conclusion, the steps used to calculate any statistical correlation between EBs and EQs are ambiguous, so they cannot be used to define a phenomenon that anticipates strong EQs. In particular, the step involving the L-shell of EBs needs to be modified for the EQ forecasting.

Unambiguous NOAA Electrons' Statistical Correlation
The NOAA-15 MEPED telescope used to monitor the electron flux coming from the zenith in three energy bands in the range of 30 keV-2.5 MeV (Evans and Greer, 2004) has been used. The energy detected for the electrons is a cumulative sum over three thresholds: E 1 30 keV, E 2 100 keV, and E 3 300 keV. Since different energies determine different behaviors in particle dynamics, new energy channels were derived from the difference of the energy thresholds to obtain electrons detected in the intervals 30-100 keV, 100-300 keV, and 300-2.5 MeV. CRs were then corrected for proton contamination (Rodger et al., 2010) from the lower energy range, based on both observations (Asikainen and Mursula, 2008) and simulations (Yando et al., 2011), and using software downloaded from the Virtual Radiation Belt Observatory (http://virbo.org/POES#Processing). Furthermore, the escaping conditions from trapped electrons were determined, thus selecting particles perturbed from the inner Van Allen Belts. These precipitating electrons were identified by calculating their minimum mirror point altitudes, h min , through the UNILIB libraries (Krunglanski, 2003). In fact, if h min < 100 km along the drift period, the electrons having energies between 30 keV and 3 MeV are ensured to be absorbed in the residual atmosphere. This occurs at the SAA longitudes due to the geomagnetic field asymmetry. Then, electrons drifting eastward and escaping the trapped conditions can be found by enforcing the condition h min < 100 km on detected CRs. Such electrons cross the NOAA altitudes and are thus able to be detected, up to the 80°longitude from the westward edge of the SAA.
The dynamics of electrons were described using adiabatic invariants such as the geomagnetic field at mirror points B m B/cos 2 α, where B is the geomagnetic field, the pitch angle α is the difference between the electron velocity and geomagnetic field directions and the L-shell parameter. CRs of precipitating electrons were thus represented in a 4-dimensional matrix (t; L; α; B) including time. The introduction of B was useful for describing the strong spatial variability of the CRs when the satellite entered the SAA (Asikainen and Mursula, 2008). B covered the range of 16-47 µT which was divided into nine nonidentical intervals, shorter where the CR was higher, to better describe the CR spatial variations, and larger where the CR was less frequent, to have a greater number of samples (Fidani and Battiston, 2008). The considered intervals in B were as follows: 16.0-17.5, 17.5-19.0, 19.0-20.5, 20.5-22.0, 22.0-25.0, 25.0-28.0, 28.0-32.5, 32.5-37.0, and 37.0-47.0 µT. Having the SEM-2 detectors with a finite aperture of 30°, the α interval was chosen to be of 15°, dividing the complete excursion into 12 equal intervals. The B-field and the L-shell were re-evaluated on the NOAA-15 orbit using the International Geomagnetic Reference Field (IGRF-12) model (http://www.ngdc.noaa.gov/IAGA/vmod/ igrf.html) (Finlay et al., 2015). To concentrate the analysis in the inner Van Allen Belts, L was restricted at the interval of 1.0-2.2 with L-shell steps of 0.1 defining 12 equal intervals. To realize a measure of precipitating electrons which are disturbed from an action coming from the Earth's surface, external Van Allen Belts are excluded by limiting L under 2.2, and the SAA is also excluded by a minimum value in B 20.5 µT.
CRs were summed on 8-s time intervals and associated to each adiabatic interval in B, α, and L. CR histograms were created, and their distributions at all the adiabatic regions resulted as Poissonian. Therefore, to obtain less than 1% probability that a CR fluctuation was of a statistical origin, the condition Poisson (CR) < 0.01 had to be satisfied for the average value corresponding to the same adiabatic coordinates of that CR. Thus, such a CR was considered to be a significant fluctuation with a probability greater than 99%, corresponding to the same adiabatic intervals. The small geometric acceptance ( ∼ 0.1 cm 2 sr) of NOAA detectors required a long time and large adiabatic intervals to obtain sufficient statistics for daily averages. However, in order to obtain a more accurate reading of the particle dynamics, small cell dimensions of adiabatic invariants should be preferred. Being so, an interpolated average value for each adiabatic interval was used to map L and B continuously from the centers of their intervals. In this way, cell sizes were not reduced, and daily averages were accepted only for those cells having at least 20 satellite passages a day. Furthermore, as CRs strongly increase near the SAA, a cubic nonlinear algorithm was used to better interpolate the averages. Starting from averages and variances, it was possible to verify if the NOAA MEPED detected any significant CR fluctuations along the entire satellite orbit.
Electron loss is primarily induced by solar activity; thus, time intervals, when the solar activity influences electron motion inside the internal Van Allen Belts, are excluded from the analysis. The exclusions are obtained by excluding time intervals when the Ap index overcame a threshold which is variable with seasons and years due to the solar cycle. Sun's rotation. Moreover, being that the electron flux was related to substorm activity , CRs were not considered for the analysis when the Dst index was lower than −27 nT satisfying these conditions. The sudden increase of electron CRs were considered EBs influenced by the Earth' surface, and more EBs detected along the same semiorbit were considered as one EB.
The correlation between EBs and EQs was calculated after defining L-shells for an EQ (L EQ ) by projecting the EQ coordinates to different altitudes and then requiring the condition L EQ − L EB ≤ 0.1. This was discussed in Reviewing NOAA Electrons' Statistical Correlations as a problem for discriminating EBs unambiguously, and another criterion is needed. Moreover, the L-shell condition is equivalent to setting a position above the future EQ epicenter where the detected EB passed. So, it might be the position where electrons escaped the trapped conditions following some remote interaction with the EQ preparation zone in the ground (Fidani, 2020). Given that each L-shell is associated with a well-defined altitude at each geographical point and a physical link is reasonably able to reach a certain maximal altitude, it is plausible that for each L-shell, there corresponds a more-or-less defined interval of geographical coordinates. A plot of EQ latitudes with respect to the L EB indicates this correspondence in Figure 2 (left). The plot shows a quadratic dependence of L EB on the EQ latitude with a minimum around 10°. This depends on the shape of the internal Van Allen Belts above the EQ epicenters (see Figure 2; right) that crosses the altitude around 2,000 km with an increasing L-shell, as the latitude moves away from 10°. This asymmetry around the equator is produced by the inclination of the geomagnetic field with respect to the rotational axes. Being so, it is enough to select only EBs with L EB in a well-defined interval to guarantee that they correlate with strong EQs in Indonesia and the Philippines.
To verify the validity of the new EB condition, a correlation was recalculated between EQs and EBs, which now includes only EB parameters, with the L EB being in a restricted interval. After a complete study to maximize the correlation with respect to many EQ and EB parameters, the validity of the new condition was FIGURE 2 | L EB dependence by the EQ latitudes which are considered for the correlation calculus when the L EQ at some altitude projection of the EQ epicenter is near L EB , on the left. The parable gives a quadratic dependence. The cause of this dependence is shown on the right where the L EB covers the expected interval, around 2,000 km of altitude above the EQ latitudes. The satellite altitude at the corresponding longitudes is well under the Van Allen Belts; the satellite altitude will approach the L EB interval at least 60°further east.

Fidani
NOAA's West Pacific Earthquake Forecasting confirmed by choosing EBs with the following: 1.21 ≤ L EB ≤ 1.31, pitch angles 56°≤ α ≤ 74°and 108°≤ α ≤ 126°, and positions −35°to 15°in latitudes and 230°−280°in longitudes. For EQs, the depth must be less than 200 km, the latitude in the −6°to 26°interval, and the longitude in the 90°−170°interval. The correlation was defined by filling the histogram {EQ;EB} (EQ × EB) with the differences Δt T EQ − T EB between the EQ time T EQ and the EB time T EB , only for those EBs with L in the interval 1.21-1.31. The correlation was also optimized according to both the time binning and the time shift. The optimization corresponded to the Δt 1.5−3.5 h interval. Here, the time difference interval, used as the binning, was suggested as the time necessary for the EB to cover the EQ longitude interval of 80°for a 60-keV electron drifting eastward at L 1.26, which was found to be about 2 h. After it, the number of correlation events is increased to 44, thus improving the correlation significance. The updated correlation is shown in Figure 3 (left), and the geographical distribution of correlated EQs is shown on the right. L EB and the EB geographic position were found to be the critical parameters to reveal EB true alarms, and the Δt 1.5−3.5 h interval is used from here on out. It should be noted that the peaks around 48 h and around 0 h are consistent with the results of the study by Anagnostopoulos et al. (2012), even if with a low significance.

Forecasting Methodologies
Following a work by Fidani (2018), the conditional probability of a strong EQ, given the EB measurement, can be calculated using the relation between the covariance and cross correlation (Billingsley, 1995). Moreover, this discussion is valid for binary events using any physical observation other than the EB, which is correlated with the EQ. If applied to the EQ and EB unitary events, the binary correlation is as follows: where the cov(EQ, EB) [P(EQ ∩ EB) − P(EQ)P(EB)] can be explicated throughout the histogram of the EQ to EB coincidences {EQ;EB} (EQ × EB) in the following population formula: where N EQ and N EB are the number of EQs and EBs that participated in the correlation, respectively, while N h is the number of total hours divided by two considered for the correlation. Being so, the probabilities of single events are P(EQ) N EQ /N h and P(EB) N EB /N h . The binary correlation histogram is then calculated as follows: The conditional probability P(EQ|EB) P(EQ ∩ EB)/P(EB) can be rewritten as follows: which means that if a correlation exists between EQs and EBs which is max Δt [ {EQ;EB} (EQ × EB)] and the time difference Δt between EQ and EB is chosen to be that of correlations, the probability of an EQ with M ≥ 6 is increased by a term proportional to the correlation. An equivalent approach to test the results obtained using NOAA particle data refers to the work by Console (2001). Here, a simple definition of an EQ forecasting hypothesis was suggested FIGURE 3 | Pearson's cross correlation between EBs and EQs recalculated using the new condition on L EB for 16.5 years of data is shown on the left; the 1 σ, 2 σ, and 3 σ thresholds are indicated by yellow, orange, and red dotted lines, respectively. The Indonesian and the Philippine strong EQs producing the 1.5-3.5-h correlation peak are shown on the right.

Fidani
NOAA's West Pacific Earthquake Forecasting with a particular sub-volume of the total time-space volume, called the alarm volume, within which the probability of occurrence of strong EQs is higher than the average. Following the work by Console, the prediction related to the detection of a precursor consists in the occurrence of an EQ event of minimal magnitude in the alarm volume. In this framework, the performance of a specific method is carried out through the statistical parameters that can be evaluated in this example, such as the success rate N S /N A , the false alarm rate (N A − N S )/N A , the alarm rate N S /N E , the failure rate (N E − N S )/N E , and the probability gain (Console, 2001), as shown below: where N S is the success number, N A is the alarm number, N E is the EQ number, V A is the alarm volume, and V T is the total volume. It should be noted that this description is completely equivalent to the previous being where A is the Indonesian and Philippine areas and V A 2A is for the alarm duration of 2 h. V A based on NOAA particle data is constant for all the alarms. Thus, the success rate is exactly P(EQ|EB), the false alarm rate is 1 − P(EQ|EB), the alarm rate is max Δt [ {EQ;EB} (EQ × EB)]/N EQ , the failure rate is 1 − max Δt [ {EQ;EB} (EQ × EB)]/N EQ , and the probability gain is as follows: where corr Δt (EQ, EB) is the particular correlation and P Δt (EQ|EB) is the particular conditional probability, both corresponding to a Δt of 1.5-3.5 h.

Evaluating Significance
A criterion for considering one model more valid than another can be made through the log-likelihood of observing that particular realization of the EQ process: under the hypothesis defining the probabilities of occurrence in P sub-volumes p i , and c i being the digital occurrence of at least one event in the subvolume, the following is the case (Console, 2001): The geographical regions of Indonesia and the Philippines are considered, where strong (M ≥ 6) EQs occurred over 16.5 years from July 1998, which were correlated with the NOAA EBs. The space-time alarm sub-volumes are in this case disjointed and separated in time only, each completely covering both areas for 2h time intervals from 1.5 to 3.5 h after the EB observations. If so, the complete volume covers N h hours, of which N EQ are those where an EQ occurred with c i 1, N EB are those where an alarm occurred with p i P(EQ|EB), N h − N EB are those where no alarm occurred with p i P(EQ), and {EQ;EB} (EQ × EB) are those where an EQ followed an alarm and where p i P(EQ|EB). The log-likelihood histogram is as follows: and it is possible to compare this forecasting hypothesis with a simpler model, called the null hypothesis, that is, the Poisson hypothesis. The success rate of this model is constantly P(EQ) and the probability gain is always G 1.0, so that the loglikelihood is calculated as follows: The log-likelihood ratio (Martin, 1971) can now be calculated to test the near-Earth space influence hypothesis against the null hypothesis for Δt of the maximum correlation, evaluating the following: Moreover, the statistical link between EQs and EBs was tested for its significance, starting from their correlation distribution (Fidani, 2015). In the work by Fidani et al. (2010), it was reported that the EQ-to-EB correlation histogram (Eq. 3), obtained collecting {EQ;EB} (EQ × EB) for many time differences, can satisfy a Poissonian process when only main shock EQs and semi-orbit EBs are used. Then, indicating the average correlation histogram with Ave, the standard deviation histogram is σ Ave √ . Being so, the number of standard deviations N σ relative to Δt, which is shown in Figure 4 of the work by Fidani (2015) for a Δt lasting 2-3 h with respect to altitude projections, can be evaluated by calculating Ave Δt and σ Δt for the same Δt as follows: The significance in terms of N σ is shown in Figure 4 for the entire interval of altitude projections corresponding to Figure 1. The significance α corr in terms of N σ can be obtained using tables of Poisson probabilities. The significance of the new correlation was also evaluated using N σ . The maximum obtained N σ 5.4 corresponded to a probability < 0.1%, not being a statistical fluctuation. A summary of the correlation calculated using Eq. 3, time interval Δt and probability gains calculated using Eq. 7, and number of events and N σ calculated using Eq. 12, corresponding to the altitude projections and to the new model based on L EB only, are reported in columns 2, 3, 4, 5, and 6 of Table 1, respectively. The values of significance calculated using tables starting from Eq. 12 are also reported in column 7 of Table 1.
Precursory information can be evaluated using Molchan's error diagram (Molchan, 2003). The quantities needed to characterize the predictive properties of a strategy in an interval (0, T) are the relative number of failures to predict for an EQ magnitude greater than M, as follows: and the relative alert time is as follows: where t a is the alert time, N(t) is the number of seismic events in the interval (0, t), dN(t) N(t + t a ) − N(t) {0 or 1} is the number of events in the alert time interval, and D e is a decision alert which can be {0 or 1} in the interval (t, t + t a ). In the NOAA electron statistical results, t a 2 h, N(T) N EQ , the relative number of failures is the failure rate in the study of Console (2001), being ] EB 1 − max Δt [ {EQ;EB} (EQ × EB)]/N EQ , and the relative alert time is τ EB P(EB). The statistical significance α Molch of a given point (] EB , τ EB ) on the Molchan's error diagram can be tested using the random probability of the Poisson model, which is the diagonal of the Molchan's diagram, and is given by the binomial distribution (Kossobokov, 2006) as follows: In light of this, the probability of obtaining N EQ (1 − ] EB ) or more hits by chance, as there have been N EQ observed target EQs, is described by the following: which produces the confidence bounds and where the index m max Δt [ {EQ;EB} (EQ × EB)]. G on the Molchan's error diagram is the slope of the line connecting (0, 1) to (τ EB , ] EB ) (Zechar and Jordan, 2008), and it is simply calculated as G (1 − ] EB )/τ EB (Molchan, 1991), which is identical to Eq. 7.

RESULTS
The EQ prediction model analyzing the EBs detected using NOAA satellites can now be represented following the more general work by Console (2001). Thus, the recent forecasting results obtained in the work by Fidani (2020) using annual averages must be redefined, in order to fit this more general representation. Subsequently, the final assessment of the hypothesis validity should be carried out via a test on a new and independent set of observations (Console, 2001).

Prediction Model
The scenario representing the model of EQ prediction needs to define volumes where EQs occur, where the target volume V T is 2-d space + 1-d time-space. In this volume, the points of EQ occurrence can be identified, together with alarm volumes V A , as success (S) and failure of prediction (F) events that are EQs occurring inside or outside V A , respectively. In this case, a precursor volume V P containing the alarm events must be defined, which is generally different from V T ; V P is the volume of the area where EB detection using NOAA satellites occurs, multiplied by the time of EB observations. An EB detection in V P is an alarm event which defines V A . With regard to the correlation mentioned above, for the Indonesian and Philippine latitudes and longitudes, V T is obtained by multiplying this area by the time spanned by the EQ observations. In this scenario, the occurrence of an EQ event is considered only with M above a magnitude threshold M 0 , where M 0 6. Unlike the models that consider EQs as precursors themselves, in this model, the EB precursor events are detected at different latitudes and longitudes, with respect to those of EQs. Correlations between EBs and EQs occurred for EB detection in the area to the west of the South Atlantic Anomaly. Thus, V T concerns the longitude interval of 90°-170°and the latitude interval of −6°-26°multiplied by the time interval of the analysis, whereas V P concerns the area of 230°-280°in longitude and the area of −35°-15°in latitude multiplied by the time interval of the analysis. V A is generated by an observation of one EB in V P . It has the same area as V T multiplied by a duration T 2 h, which results in the width of the correlation peak, occurring over the next 1.5-3.5 h, which is found to be the time position of the correlation peak. In this case, the alarm volume V A is constant for all alarms. A success is added if an EQ with M ≥ 6 occurs in the V A . A failure is added if an EQ occurs out of the V A , which means not included in the time intervals, and any alarm EB not followed by EQ is classified as a false alarm. The described model can be represented by the three-dimensional space of geographical coordinates and time reported in Figure 5. Performance of the NOAA EB detection is conditioned by a V P which is not continuous. In fact, as reported above, only a low number of days are magnetically calm enough to be used for the analysis. Moreover, the NOAA satellites go into the detection area west of the South Atlantic Anomaly intermittently, thus further reducing the total time of observations. FIGURE 5 | Volume representation where the forecasting model can be tested is delimited by the product among the geographical coordinates of EQs and EBs and the time of observations. V T , V A , and V P are all discontinuous volumes in time as the EB analysis is suitable for EQ forecasting when the solar activity is very low. Moreover, the alarm duration is 2 h. V T and V A cover the entire West Pacific area, and V P covers the different geographical areas on the western North American and South American coasts. The causal link between EQ disturbance and EB measurement events is represented by green dashed arrows. The possible precursor volume due to a hypothetical physical action on the ionosphere above the epicenters of the earthquakes is represented in red. This is the cause of a noncontinuous V T , where V A appears to fill the same geographical area as V T for a time interval of 2 h. A V A within V T is generated 1.5 h after an EB is observed in V P . This is different from the model based exclusively on seismic activity, where the causal link between V A and the precursor is near the vertical, given the seismic properties to cluster. The causal link between V A and the EB observation event is represented near the horizontal. This causal link is the eastward electron drift, according to the electron energy, which is represented in Figure 5 by green thick arrows moving according to the longitude. The vertical distance between the starting point of arrows above the EQ epicenters and the base of V A is the time anticipation of a possible physical interaction relating the EQ preparation volume to the ionosphere. In the analysis of 16.5 years of NOAA data, V T concerned only the Indonesian and the Philippine areas multiplied by N h 6, 953 h. It should be noted that this value and the following are different from those reported in past publications (Fidani, 2020), as the past reports were rough estimates. In this volume, the following occurred: a total number N EQ N E 600 of EQs with M ≥ 6, a total number N EB N A 1, 892 of alarms, and a total number N S 44 of success. Being so, the success rate of this model is N S /N A 0.023, the false alarm rate is 1 − (N S /N A ) 0.977, the alarm rate is N S /N E 0.073, the failure rate is 1 − (N S /N E ) 0.927, and the G 3.0. Here, the target volume is subdivided into nonoverlapping sub-volumes with time intervals of a day that fill V T completely. For each day sub-volume, the probability of occurrence of at least one target event is estimated to be equal to P(EQ) with no EB observed and P(EQ|EB) with one EB. Analyzing the data, days with more than one burst can be found with a frequency of about 20%. These bursts can be far away in time when the time difference is more than 10,000 s (∼2.8 h) or neighbors when the time difference is 5,000-7,000 s (∼1.4-2.0 h); in the latter case, they belong to successive orbits according to NOAA POES orbit parameters. It should be remembered that all EBs in a semiorbit were considered as only one EB to be counted for the correlation calculus (Fidani, 2015). When two bursts are far away, the time alarm of the first ends before the beginning of the second detection, so two disjoint V A with the same increasing of the conditional probabilities occur. When detected bursts belong to successive orbits, the V A time interval of 2 h is greater than the range of 5,000-7,000 s, and a partial overlapping between two consecutive alarms will also occur. The few cases of overlapping for NOAA alarms are not considered here and will be presented in a future publication.

Dependent Observables
Bayes' theorem allows us to compute the probability that a hypothesis is true, provided that one knows the probable truth of all supporting arguments. It reverses the conditional probabilities and defines the probability of the hypothesis given the evidence. It shows that there is a significant probability gain in using precursors for prediction, even if a phenomenological occurrence is not the proof of a precursor. However, starting from the correlation results (Fidani, 2015), Bayes' theorem, as shown below, P(EB|EQ) P(EB) P(EQ|EB) P(EQ) , can be employed using statistical bases (Eq. 5). Being so, the alarm rate can be rewritten as P(EB|EQ), the failure rate as 1 − P(EB|EQ), and the probability gain as follows: Bayes theorem allows us to compute the probability of an EB given the measurement of an EQ. If calculated, it appears surprisingly high, equal to 0.8 for the correlation. However, this result must not be misinterpreted. In fact, the probability gain remains the same, suggesting that the high probability of detecting an EB when an EQ is observed is due to its greater frequency of EB occurrence. Bayes' theorem tests hypotheses and can be updated on the basis of new information. If used under the condition of many independent precursors EB, EC, ED . . ., the EQ conditional probability in a small time interval of a given area after the simultaneous detection of one EB, one EC, one ED ... can be approximated by (Aki, 1981) the following: It is the product among the unconditional probability and the probability gains of each precursor; those are the ratios between the conditional probability of each precursor and unconditional probability. However, the condition of independent precursors is difficult to prove, and from the studies reported in the Introduction, a set of physical links for only a part of them is suspected (see, for example, the study by Pulinets et al. (2015)). Therefore, it seems that dependent candidate precursors represent the most common occurrence. Therein, the conditional probability of an EQ with a magnitude greater than M 0 is not increased by a product of a further probability gain of another detected precursor, if this has a certain degree of dependence on the first to be detected. Thus, the conditional probability cannot be approximated by Eq. 19; it must be recalculated. Starting with only two dependent precursors EB and EC that generate alarms in the same V T , the conditional probability on EQ, given the observations of both precursors EB and EC, can be expressed using the relations (Eq. 7) and (Eq. 1) as follows: where the covariance can be explicated throughout the histogram of EQ to EB∩EC coincidences {EQ;EB∩EC} [EQ × (EB∧EC)], and by considering the total number N EB∩EC of correlated precursor events. Finally, P(EB|EC) can be calculated through the corr(EB, EC) of the relation (Eq. 5), as with all the other conditional probabilities. However, a more interesting question might arise upon using two observation networks whose observables are dependent: what is the overall probability gain upon using observations without differentiating them? In this case, the warning corresponds to a detection from the set of Frontiers in Earth Science | www.frontiersin.org August 2021 | Volume 9 | Article 673105 11

Fidani
NOAA's West Pacific Earthquake Forecasting dependent observables EB∪EC, which means that the probability, as shown below, is conditioned by the observation of an EB or an EC or of both an EB and an EC, indifferently. Then, using simple algebra with P(EQ∩(EB∪EC)) P((EQ∩EB))∪(EQ∩EC)), we obtain the following: where P(EQ ∩ EB) P(EQ|EB)P(EB), P(EQ ∩ EC) P(EQ|EC)P(EC), P(EB ∩ EC) P(EB|EC)P(EC), P(EQ ∩ EB ∩ EC) P(EQ|EB ∩ EC)P(EB ∩ EC). Considering the relation (Eq. 7), the probability gain due to the observation of an event in EB∪EC is as follows: where G EB and G EC are the probability gains of the single precursor EB and the single precursor EC, respectively. G EB∩EC is the probability gain (Eq. 20) due to the observation of both EB and EC events that are correlated between them. Alternatively, using the relations (Eq. 7) and (Eq. 1), we obtain the following: where the covariance can be explicated throughout the histogram of EQ to EB∪EC coincidences {EQ;EB∪EC} [EQ × (EB∨EC)], and by considering the total number N EB∪EC N EB + N EC − N EB∩EC .

Different Precursors Combined
It needs to be highlighted that a statistical correlation between the two time series does not generally mean that they are physically related (Aldrich, 1995). A causal link has been hypothesized between EB measurements and EQ occurrence (Fidani, 2018;2020), which supports the validity of the hypothesis. Although the Δt T EQ − T EB time difference of the correlation is in agreement with the physical migration of electrons eastward, this migration has not yet been observed for EBs correlated with EQs. Furthermore, a physical link between the EQ preparation zone and the ionosphere above the future epicenter, separated by about 2,000 km, has not been demonstrated. Finally, the existence of some physical phenomena occurring at the future EQ epicenter, which is enough to influence the ionosphere, remains only a hypothesis until all of these passages have been fully demonstrated. To discover EBs with the correct times at different longitudes, which would satisfy their physical migration in the ionosphere, more satellites are needed. This verification is currently possible for EBs, as different NOAA satellites fly together, even though cases having suitable satellite positions must be found. Till date, this has not been calculated, as no correlations have been found with EBs selected using other satellite databases during the same periods.
Regarding a physical link, able to cover more than 2,000 km between the Earth's crust and the ionosphere throughout the atmosphere, magnetic pulses have been hypothesized to influence the high-energy charged particles' motion by pitch angle diffusion (Galper et al., 1995). Even if other processes have been proposed more recently, such as the injection of radioactive substances and charged aerosols into the atmosphere, leading to a change in the vertical electric current in the atmosphere and to a modification of the electrical field in the ionosphere (Sorokin et al., 2001), from a solidstate physics perspective (Freund, 2011), and radioactive ionization to model the lithosphere-atmosphere-ionosphere-magnetosphere coupling by the latent heat flux (Pulinets et al., 2015), by AGW (Yang et al., 2019;Piersanti et al., 2020). Recent measurements of magnetic pulses on the occasions of strong EQs (Bleier et al., 2009;Orsini, 2011;Nenovski, 2015) have suggested an efficient process which may influence the ionosphere . In this latter work, the magnetic data analysis at the Chieti Station of the Central Italy Electromagnetic Network, performed using two independent sample systems of the same signal, showed that a large number of pulses were recorded in the ELF band below 10 Hz with amplitudes mostly in the range of 2.5-80 nT. Specifically, the model proposed for analyzing magnetic pulses consisted of diffused underground electrical currents throughout a conductive strip between the Apennines and the Adriatic Sea. The current required to induce detected pulses is greater than i 40 kA for the strongest pulses. Unipolar magnetic pulses can be, for example, generated deep in the rock column by peroxy defects when rocks are subjected to increasing deviatoric stresses (Freund et al., 2021). The proposed strip of diffuse current constitutes the electromagnetic source of ULF waves, which is able to produce low intensities of magnetic inductions on the Earth's surface, even if it is measurable both near and far from the EQ epicenter. The strip of diffuse current is able to produce significant fields, also far from the ground by integration. To demonstrate this and calculate magnetic induction in the ionosphere, two simple models can be considered. An infinitely long strip of about 6-km-thick and 150-km-large conductive soil has been considered to calculate the magnetic induction very close to the larger strip surface. In a real case, a finite-length strip should be utilized, which gives a correct result even for a coil magnetometer very close to the strip. Referring to Figure 6, on the left side, where the lines of induction are generated by a not-finite strip and wire lengths are compared, it is possible to see that the B intensities near the strip are lower than those near the wire. This is due to the current density differences flowing near the observational point. Relations (Eq. 2) and (C4) of the study by Fidani et al. (2020) were used for the wire and the strip, respectively. Moving away from the current density, magnetic inductions generated by the wire and the strip become equal, when the distance overcomes the strip width. However, currents which are not of finite length deviate too much from the real case at large distances. So, a finite length l of the wire is considered and the magnetic induction is simply calculated as follows: where d is the distance. Figure 6, on the right side, depicts a comparison of B intensities of the two models with the distance. The B intensity results are limited to 0.2 nT at distances of 1,000 km with a total current of 40 kA. This is in agreement with the intensities obtained from the theoretical calculations, which have shown that only a magnetic type source with frequency < 10-20 Hz can be effective for the penetration of fields into the upper ionosphere and magnetosphere (Molchanov et al., 1995). Moreover, magnetic disturbances were observed above moderate EQ epicenters for the frequency band of 0.1-10 Hz (Strakhov et al., 1994). From the upper ionosphere, these waves travel as Alfven waves further along the geomagnetic field line and reach the inner boundary of the inner Van Allen Belt. Alfven waves are thought to resonantly interact with trapped charged particles; this process is most intensive in the equatorial part of the magnetosphere, for L-shell values equal to or less than 2 (Molchanov and Mazhaeva, 1993). In fact, the measured pulse frequency under 10 Hz is around the bouncing resonance of electrons with energies between 30 and 100 keV (Walt, 1994); it is exactly as measured on board NOAA satellites. This means that energy can be efficiently transferred from magnetic pulses producing electron pitch angle disturbances (Galper et al., 1995). Furthermore, the B intensity due to an Earth surface strip of current is a very low value, compared to an example of geomagnetic activity. However, the frequency range of geomagnetic activity is about 0.0001-0.01 Hz (Francia and Villante, 1997), so that the B rate is out of resonance. Finally, currents produced by lightning are generally lower in intensity, around i 10 kA, with frequency emissions in the upper part of the ELF band and the VLF band still out of bouncing resonance (Inan et al., 2010). Indeed, magnetic signals have been correlated with strong EQs in different regions of the world such as Japan (Ohta et al., 2013), Kamchactka , and California (Kappler et al., 2019). The Japanese and Kamchactka studies obtained correlations with time differences of 2-5 days before seismic events, and probability gains of about 1.6 (Han et al., 2014) to 2.6  were reported. Magnetic pulses, from here on identified with ECs, have been hypothesized to induce EBs, and therefore, they may be considered dependent events for hypothesis. The relation (Eq. 23) can be used to calculate the probability gain of EQ probability due to possible observations of both EBs and ECs. This possibility can be useful as the observation of ECs on the Earth's surface can occur when a NOAA satellite is not in a suitable position to reveal the possibly induced EBs. Moreover, magnetic detectors cannot be installed at some positions offshore or are not able to detect ECs where there are EQs. In other words, EBs and ECs can compensate for the gaps on each other in a forecasting experiment, where V T is very large, as for the West Pacific. The magnetic pulse precursor can be represented by the red volume of Figure 5. Being so, the mutual interrelation of the sets of EBs, ECs, and EQs is a completely real case, which can be studied for a total probability gain with the expression (Eq. 23). To show the advantage in using two dependent observables, a forecasting experiment can be imagined where both EB and EC data are available. Here, one can imagine the presence of a magnetometer network distributed on several islands of the West Pacific. Although not still existent, this network is currently achievable. We suppose that it exists and is able to obtain a G EC 1.6 for EQs within a certain distance around the stations, with a time advance of ECs of 4-6 h with respect to EQs. P(EC) is necessary for Eq. 23 and is supposed to be 0.05, which means N EC 324 magnetic pulses, or sets of magnetic pulses, considered as magnetic alarms on the same time interval as for EBs. For what concerns the NOAA observable, a G EB 3 was found from the correlation analysis. However, even if ECs from days with the Ap index above the previously defined threshold were excluded as for EBs, the daily NOAA observation time was always half a day (Fidani, 2020) due to the satellite's orbit crossing the EB detection region. Being so, the probability gain should be evaluated on a double time interval and being able to detect the same number of EBs. The probability gain is FIGURE 6 | Comparison of magnetic induction lines produced by the rectangular and circular sections in the center on the left; continuous lines are generated by the rectangular strip (s), whereas the dotted lines are generated by the wire (w). B dependence, from the distance at half of the strip width of the currents, is shown on the right for both conductors that are infinitely long or of length l using the relation (Eq. 25) and the relation (Eq. 2) of the study by Fidani et al. (2020 . It would also be necessary in this case to consider a range of possibilities from 0, no correlation between EQ and EB∩EC, to 10 common events. The probability gain improvements that are obtainable thanks to a network of EC measurements added to the NOAA satellite EB detection are reported in the contour plot of Figure 7. The maximum value of G EB∪EC 2.05 is obtained for 314 ECs which are correlated with EBs, but none of these EBs are correlated with EQs, shown by point (314, 0) of Figure 7, which means that 10 ECs are correlated with EQs and another part of EBs not correlated with ECs is correlated with EQs; thus, the total number of different observations correlated with EQs is increased by 10. G EB∪EC is a minimum when the correlation events between EBs and ECs are max Δtp [ {EB;EC} (EB × EC)] 0, shown by the point (0, 0) of Figure 7. That is, when there are no ECs correlated with EQs nor with EBs, and G EB∪EC is slightly less than 1.6 because it is a weighted average between G EB and G EC . When the EB-to-EC correlation increases to 10 events, a maximum of 10 correlations with EQs can happen, as shown by point (10, 10) of Figure 7, where G EB∪EC is still lower due to redundancy. G EB∪EC does not exist to the left of the white dotted line and on the right of the black dotted line. When the EB-to-EC correlation is at the maximum and the correlation of such ECs with EQs reaches 10 events, as shown by point (324, 10) of Figure 7, G EB∪EC returns as near 1.5 due to the low number of correlated EQ events.

CONCLUSION
A complete EQ forecasting methodology has been considered in this study based on NOAA satellite high-energy electron detection. It utilizes recently discovered correlations existing between EBs selected from the NOAA database and strong EQs collected at the USGS (Fidani, 2015), which can be classified as an electromagnetic phenomenon and a nearseismic precursor. This correlation, concerning EBs anticipating main shocks, has been obtained also thanks to the geomagnetic disturbance database of the Ap and Dst indexes, as well as the IGRF model. The methodology has been represented by the volumes of target, alarm, and precursor, for testing EQ forecast hypotheses (Console, 2001).
To systematically test this methodology, a quantitative and rigorous definition of the anomaly is given according to a statistical criterion with respect to the Poisson distribution of electron CRs. Here, electrons must be escaping the trapping conditions, that is, precipitating, probably due to a disturbance. Finally, novel in relation to previous publications (Fidani, 2015;Fidani, 2018;Fidani, 2020), the parameter L-shell of the electrons has been disentangled from apparent L-shells associated to EQs. In fact, the condition of the difference L EB − L EQ ≤ 0.1 in previous publications was substituted with the 1.21 ≤ L EB ≤ 1.31, which was deduced from the correspondence between electron L-shells at around 2,000 km and EQ latitudes. This disentangling between anomalies and EQs is essential to provide a precise definition of the observed phenomenon (Wyss, 1997) in order to carry out a forecast.
The statistical correlation between EBs and EQs has been extended to a time difference of 1.5-3.5 h, thanks to the hypothesized causal connection of drifting electrons. Maximizing the significance of this correlation indicated that EQs still belong to the Indonesian and the Philippine areas, collecting more seismic events from the West Pacific. These EQs occurred as in the previous analysis, mainly in the sea, with a depth up to 200 km to correlate with EBs. Regarding EBs, the pitch FIGURE 7 | Probability gain (Eq. 23) due to the contribution of EB and EC detection, where in this case ECs are ULF magnetic pulses. This gain is calculated with respect to the maximum event correlation between EBs and ECs on the abscissa, which comes from P(EB|EC), and to the maximum event correlation between EQ and EB∩EC on the ordinate, which is derived from P(EQ|EB∩EC). Δtp T EB − T EC is thought to be the drifting time of electrons from the ionosphere above the epicenter to the cross with the NOAA satellite. Therefore, G EB∪EC is described by the contours if Δt of EB-to-EC correlation is shifted by Δtp, that is, the dependence between EBs and ECs. The white dotted line on the left and the black one on the right represent G EB∪EC limits of validity. angle intervals were restricted, even if the number of EBs increased. To demonstrate that the correlation calculus is completely equivalent to the frequency approach of Console (2001), the probability gain was recalculated in terms of conditional probability and correlation for digital events.
After having stressed that a statistical correlation between two time series does not generally mean that they are physically related (Aldrich, 1995), the hypothesis concerning a physical link between EBs and EQs was studied. A recent observation of magnetic pulses recorded before strong EQs in Central Italy provided recorded magnetic intensities with a diffuse current model . This model can push back to the hypocenter region the causal connection of physical events, even if the deduced magnetic inductions in the ionosphere must be demonstrated to be able to modify the electron pitch angles. Being magnetic pulses measurable on the Earth's surface, they might be precursors, and indeed, statistical correlation of magnetic pulses was found (Han et al., 2014;Hayakawa et al., 2019). However, given their possible causal link with EBs, magnetic pulses and EBs could not be independent precursors. Thus, starting from Bayes' theorem, a more general relationship of the probability gain due to the combination of two precursors is expressed in terms of single probability gains of each precursor and the correlation between the precursors. An example of improvement in the probability gain due to a couple of digital dependent precursors is tentatively calculated for the first time. A dependence between the precursors is introduced in the probability gain (Eq. 23) by a time shift which correlates the precursors. The best probability gain is obtained for the maximum correlation between precursors not correlated with EQs.
Finally, this methodology is general enough that it could be adapted to the combinations of observations from both the Earth's surface and space, such as electromagnetic, seismic, or other physical observables. To do this, a series of steps must be performed: 1) collect data from the same instrument(s) with the same environmental conditions for many years, 2) search for anomalies of a physical observable with a statistical rigor, following a physical idea of possible equilibrium disturbances, 3) calculate a statistical correlation between EQs and anomalies by selecting physical parameters disentangled from EQ parameters, following a physical idea of possible interaction, 4) calculate the correlation significance, or the likelihood, or the Molchan's error diagram and optimize it with respect to the physical parameters, 5) use the more relevant parameters to determine the correlation significance, or the likelihood, or the Molchan's error diagram for a physical model refinement, 6) demonstrate the Δt and Δtp agreement of the two observables using a unified physical model, and 7) calculate the probability gain to one or more precursors and to their combinations, and verify the results in a target volume of future times or different databases. If step 6 is not obtained, the probability gain (Eq. 23) can be maximized with respect to Δt and Δtp to suggest a probable unified physical model. Moreover, an experiment for the EQ forecasting test in Indonesia and the Philippines is currently feasible using the NOAA-15 satellite, given the presence of the United States West Coast antennas (Fidani, 2020). This could be concluded over a few years with a reasonable response, due to the high frequency of strong seismicity in Indonesia and the Philippines.

AUTHOR CONTRIBUTIONS
CF: scientific analysis and manuscript writing.