© S. Aminalragia-Giamini et al., Published by EDP Sciences 2020

1 Introduction

Solar Energetic Particle (SEP) events and their proton component (Solar Proton Events, SPEs) are some of the most hazardous phenomena of space weather, they contain populations of ionizing radiation and as such they can present dangers for the increasingly complex electronics on board spacecrafts such as those of commercial satellites and sensitive instruments of scientific missions. This is even more important for missions outside the Earth’s protective magnetosphere for satellites in the Sun–Earth Lagrangian Points as well as missions to Mars since SPEs can reach into the interplanetary space and still retain much of their hazardous potential (Jiggens et al., 2019). Moreover, SPEs can cause harm to humans either through extreme exposure, e.g., while in extravehicular activities, or through accumulated exposure over the course of a mission; the latter being increasingly important for manned missions to the Moon and Mars as those envisioned by Space Agencies and private organizations. At the same time SPEs offer a fascinating insight into the complex heliospheric processes (Anastasiadis et al., 2019) and while mechanisms such as the generation, acceleration (Afanasiev et al., 2018) and propagation (Wijsen et al., 2019) of solar energetic protons are not fully understood their ongoing study can uncover characteristics that could help in the in-depth analysis and potential physical modelling of such mechanisms (e.g., Werner et al., 2019). Due to their hazardous consequences the space weather scientific community has engaged over the years in forecasting SPEs in various ways, usually entailing at least partially a statistical or probabilistic approach (Camporeale et al., 2019). State of the art forecasting schemes include the University of Malaga Solar Energetic Particle (UMASEP) system (Núñez, 2011, 2015), the FORecasting Solar Particle Events and Flares (FORSPEF) system (Papaioannou et al., 2015, 2018), the Warning System for Aviation Exposure to SEP (WASAVIES) (Kataoka et al., 2014) and WASAVIES – Earth Orbit (WASAVIES-EO) (Sato et al., 2019) systems, the Relativistic Electron Alert System for Exploration (RELeASE) (Posner, 2007; Malandraki & Crosby, 2018), the Empirical model for Solar Proton Events Real Time Alert (ESPERTA) (Laurenza et al., 2009; Alberti et al., 2017) and the Proton Prediction System (PPS) (Smart & Shea, 1989, 1992; Kahler et al., 2007). UMASEP predicts both when the integral proton flux intensity (E ≥ 10 MeV) will reach values of 10 p cm⁻²sr⁻¹s⁻¹ as well as the flux intensity-time profile in the first few hours of an SPE occurrence. FORSPEF provides forecasts for solar flares and Coronal Mass Ejections, as well as the occurrence of SPEs. WASAVIES-EO outputs nowcasts and forecasts for SEP fluxes in the magnetosphere and ionosphere during Ground-Level Enhancements (GLE). RELeASE uses measurements of relativistic electrons to provide short-term forecasting on the occurrence and intensity of SPEs. The ESPERTA model uses as input parameters the flare location, Soft X-Ray (SXR) fluence and ∼1 MHz radio fluence to provide a warning following the SXR peak for ≥M2 solar flares. PPS provides predictions of solar proton intensities at 1 AU following solar flares based on average observed SEP profiles, peak intensities, and event durations using as input parameters the solar flare X-ray peak or fluence, the time of onset and maximum, as well as the solar flare location. Finally, the SOLar Particle ENvironment Code (SOLPENCO; Aran et al., 2006, 2008) was the first tool developed to provide a prediction of the shape of intensity-time profiles of SEP events associated with interplanetary shocks. SOLPENCO and its more recent version, SOLPENCO2, are derived from physics-based models (e.g., Pomoell et al., 2015) and the latter was used in the interplanetary statistical SEP model of SEPEM (Crosby et al., 2015) to provide predictions away from 1 AU. Overall, forecasting models and systems typically provide a probability of Event occurrence, along with an accompanying warning, and an estimation of important SPE characteristics, usually the Peak flux over a range of proton energies.

Recently Kahler & Ling (2018) studied the relationships of SPE peak flux and fluence values using the data detailed in Papaioannou et al. (2016). They were the first to show that linear fits in the log–log space (henceforth called log-linear fits) can be made from SPE Peak flux values to Fluence values at separate energy channels with varying degrees of scattering. In this work, we have expanded upon these findings and show that robust relationships exist between the whole differential Peak and Fluence spectra of each individual SPE across a wide range of proton energies, 5–200 MeV. Furthermore, it is shown that when fitted with appropriate analytical functions the Peak and Fluence spectra of each SPE show closely similar spectral indexes as well as strongly correlated intensity factors. Additionally, the offsets from the analytical fitting exhibit similarities and can be translated from the Peak spectrum to the Fluence spectrum. To our knowledge, these are relationships not shown or quantified before and offer further insights into two very important characteristics of SPEs which can be valuable for the analysis, categorization and modelling of such Events. Furthermore, using these relationships we have devised a simple yet robust methodology for the direct and accurate prediction of the total Fluence spectrum from a known or forecasted Peak spectrum of an SPE. Typically in systems such as those discussed above, Peak flux spectra are forecasted and predicted ahead of the actual SPE onset or before a detected SPE reaches its Peak fluxes. With this methodology a forecasted Peak spectrum can be used to directly derive a self-consistent prediction of a Fluence spectrum or conversely if both Peak and Fluence spectra are forecasted our findings provide a way to check their internal cross-consistency. Due to its simplicity, this method can function as an easy to implement standalone module in existing or future forecasting systems. Results are shown for Fluence predictions and they are validated using two different Event lists, the SEPEM Reference Event list and the Event list defined within the context of the FORSPEF system in order to verify that list-specific attributes do not influence the results. Finally, we investigate and demonstrate the applicability of this approach even in nowcasting applications as a way to obtain a good prediction on the Fluence spectrum at the moment an SPE reaches its maximum proton flux intensity.

2 Methodology

2.1 SEPEM Reference Event List and data

In this work, we present results using the SEPEM Reference Event List (REL) http://sepem.eu/help/event_ref.html which contains 266 SPEs from 1974 to 2015. In SEPEM REL a proton enhancement is defined as an Event if the proton flux intensity of the second channel (E = [7.23–10.46 MeV]) rises above 0.01 p cm⁻²sr⁻¹s⁻¹, remains for at least 24 h above this threshold, and has a peak value of at least 0.5 p cm⁻²sr⁻¹s⁻¹. The Event ends when the flux intensity drops below the 0.01 threshold, and if two qualifying enhancements occur within 24 h they are compounded into one Event, more details can be found in http://sepem.eu/help/event_lists.html. Using the SEPEM REL the corresponding events were extracted from the SEPEM Reference Dataset v2.0 (RDSv2) http://sepem.eu/help/SEPEM_RDS_v2-00.zip. SEPEM RDSv2 contains flux intensity time-series based on NOAA GOES measurements. They are continuous, cleaned, and cross-calibrated with IMP-8 (Sandberg et al., 2014) differential solar proton flux intensity time-series at 5 min intervals with no missing values or time gaps. The fluxes are rebinned at 10 energy bins exponentially spaced within 5–200 MeV and we use the proton time-series as provided at the 10 energy bins and their respective effective energies as shown in Table 1. We note that integral products of RDSv2 have also been validated using STEREO data (Rodriguez et al., 2017).

2.2 Log-linear relationships of Event Peak Flux and Fluence spectra

For each Event in the SEPEM REL we have derived the Peak flux spectrum as the maximum values of the differential proton flux intensity time-series at each energy bin and the Fluence spectrum as the integrals of the differential proton flux intensity time-series at each energy bin over the duration of the Event. We note that the Peak flux times differ for each energy bin. The resulting 266 pairs of Peak and Fluence spectra show linear relationships in the logarithmic space. Figure 1a and b shows two examples, one with a very strong log–linear relationship and another with larger residuals. Figure 1c presents all the Peak-Fluence pairs and as it can be seen the same trend appears between Peaks and Fluences for all Events. In order to quantify these observations linear fits y = kx + m were performed fitting the natural logarithm of the Fluence spectrum to the natural logarithm of the Peak spectrum for each Event, using the least-squares method. It is found that the slopes k of these log-linear fits are clustered close together having a very small range around a mean value of 0.892 with a standard deviation of 0.104 and few outliers. The main differences exist in the constants m, which are dependent on the actual values of the Fluences and Peaks. Figure 2 shows the values of the slopes k versus the constants m derived from the log-linear fits for all Events. The values are calculated at 95% confidence intervals with resulting ±0.04 and ±0.167 mean margins of error for k and m, respectively. We note the k and m parameters show a weak anti-correlation of −0.281. Below we discuss how this relationship can be quantified and used for the prediction of the Fluence spectrum when the Peak spectrum is available, being measured either in near-real-time or forecasted.

Fig. 1 (a) Fluence-Peak pair showing strong log-linearity, (b) Fluence-Peak pair with larger residuals and (c) all Fluence-Peak pairs from the SEPEM list.

Fig. 2 Slopes k versus constants m of the y = kx + m log-linear fits from all pairs of Peak and Fluence spectra.

2.3 Spectral parameters of differential Fluence and Peak spectra

The relationship between all pairs of Peak and Fluence spectra as well the narrow range of slopes shown above implies the existence of a global relationship of spectral indexes. We therefore perform analytical fits of the Fluence and Peak spectra separately in all 266 SPEs across the 10 differential energy bins. Exponential cut-off power-law functions are used, shown below in equation (1). We note the fits are performed in the log-log space:

(1)here S is the differential spectrum (Fluence or Peak flux), a is the (logarithmic) intensity of the spectrum, b is the power-law spectral index and c the exponential spectral index. An example of analytical fits of the Peak and Fluence spectra of an Event can be seen in Figure 3. While power-law, which is typically used for SPE spectra fitting and is associated with shock-acceleration (Ellison & Ramaty, 1985), can provide adequate fits we have used the exponential cut-off power-law function because it is more versatile. Due to the added free parameter this analytical function can capture well spectra with power-law behavior as well as those that diverge from it showing exponential drop-off, especially in higher energies. The use of exponential cut-off power-law analytical function for SPE spectra as introduced by Ellison & Ramaty (1985) has been used with good fitting results demonstrated in Grimani et al. (2013) and Adriani et al. (2011) found SPE flux spectra which were exponential or quasi-exponential.

Fig. 3 (a) Peak spectrum derived from SEPEM RDS for an SPE Event and the analytical fit made with exponential cut-off power-law. (b) Same for respective Fluence spectrum. The fitted spectral parameters can be seen in the bottom left corners.

The fittings of all the spectra result in 266 triplets of [ a , b , c ] parameters for Peaks and the same for Fluences. The values are calculated at 95% confidence intervals with resulting ±1.19, ±0.47 and ±0.01 mean margins of error in [ a , b , c ] respectively for Peak spectra and ±1.04, ±0.41 and ±0.009 for Fluence spectra. Figure 4 shows cross-plots between the respective spectral parameters of Fluences and Peaks where it is evident there are strong one-to-one relations (see Table 2).

Fig. 4 Cross-plots of the spectral parameters of Fluences and Peaks showing the relationships between them. (a) Spectrum intensity a , (b) power-law index b and (c) exponential index c .

Furthermore, the spectral parameters, a , b , c , show interdependencies among them for both Peak and Fluence spectra separately. The strongest correlations appear for the a – b and b – c pairs. Figure 5 shows plots for these two pairs of parameters from Peak spectra with anti-correlations of −0.75 and −0.81, respectively.

Fig. 5 Cross-plots of Peak spectral parameters for all 266 Events. (a) Spectrum intensity a P with power-law index b P and (b) power-law index b P with exponential index c P .

Since the parameters exhibit one-to-one (Peak to Fluence) relationships as well as internal interdependencies it stands to reason that the relations of all Peak parameters to all Fluence parameters can have statistical significance as well. Indeed, as shown in Table 2 the correlations are statistically non-negligible. As expected, the highest correlations appear for the counterpart parameters in the diagonal of the table (as demonstrated in Fig 4); however, non-negligible values appear for all pairs except for [ c _P a _F ].

2.4 Estimation methodology of Fluence spectral parameters

2.4.1 Parameters of analytical Fluence spectrum

Based on the findings and analysis discussed above we employ the [ a _P , b _P , c _P ] parameters of the Peak spectrum as estimators for the three [ a _F , b _F , c _F ] parameters of the Fluence spectrum. This is done by deriving three fitting curves, one for each Fluence parameter. These three curves serve as the estimating functions; they are quadratic hypersurfaces in the 4-dimensional space defined by the three estimators [ a _P , b _P , c _P ] and one predicted Fluence parameter. The quadratic form of the hypersurfaces is shown in equation (2):

(2)where m _j j = 1:6 are constant terms, d is the intercept and y is 1 of the 3 Fluence parameters [ a _F , b _F , c _F ], where each one is separately resolved. We use the three fitted curves to estimate (reconstruct) the Fluence parameters [ a _F , b _F , c _F ] for each Peak flux spectrum in the REL. Figure 6 shows these estimations where the real fluence spectral values (fitted from the data) are plotted against the values recovered from the curves. The estimated (reconstructed) parameters agree well with the values from the data falling along the y = x equality lines with small scattering and few outliers.

Fig. 6 Cross-plots of Fluence spectral parameters derived from data and estimated from the estimator hypersurfaces, y = x equality lines are in blue. (a) Fluence spectrum intensity a F , (b) Fluence power-law index b F and (c) Fluence exponential index c F .

Using these hypersurfaces the Fluence intensity and spectral indexes can be directly estimated and thus a full spectrum can be calculated assuming an exponential cut-off power-law (cf. Eq. (1)). The estimated Fluence spectrum can be calculated for an arbitrary binning within the 5–200 MeV range used in the SEPEM RDS. Additionally, using the analytical function it can be also estimated for energies outside this range. We note that for extrapolations beyond the data energy range, special care needs to be taken when the spectral index c is estimated to be positive. Such positive values are almost certainly due to the effect of background in the flux intensity time-series which skew the spectra at higher proton energies and locally change the slope. However, background removal is not necessarily a straight-forward task and we have opted to use the time-series as they are provided. Therefore, we note that since at higher proton energies the exponential term tends to overwhelm the power-law behavior, positive exponential terms can lead to unphysical estimations outside the fitted 5–200 MeV range. Despite this, in many Fluence spectra positive c values do provide a much better fit to the data and therefore can offer good estimations inside the used energy range.

2.4.2 Residuals from analytical function fitting

Following the Fluence spectral parameter estimation, an additional step can be implemented that is based on the Peak values but limited to the 10 energy bins. As is the case for every type of analytical fitting, the use of the exponential cut-off power-law function to fit the Peak and Fluence spectra is bound to result in residuals from the actual data values. These residuals are just positive or negative logarithmic offsets from the strict analytical form. The residuals are calculated according to equation (3):

(3)where S is the data spectrum (Peak or Fluence), AS is the analytical spectrum calculated from the fitted spectral parameters, E is the energy and i is the index of the energy bin. We show that the residuals from the Peak spectrum can also be used to estimate the respective expected residuals for the estimated Fluence spectrum. The estimation of these residuals offers added precision and can adjust the estimated Fluence spectrum to match even more closely the data. The residuals from the Peak fits (ΔPeaks) show strong correlations or anticorrelations with the residuals from the Fluence fits (ΔFluences) across all energies. This is not very surprising given the shown similarities in the spectral indexes between Peaks and Fluences, as one would expect similar deviations from the analytical form between spectra with closely similar profiles. We note that this of course means that the residuals across energies within Peaks and Fluences separately also show correlations among them. Figure 7a and b shows cross-plots of the ΔPeaks and respective ΔFluences from all Events at 6 MeV and 115 MeV. Both show strong correlations above 0.8. Figure 7c presents a mapping of the correlation coefficients between ΔPeaks and ΔFluences at all energies. In total 100 values were derived by calculating the correlation coefficients between all 10 Peak residuals and all 10 Fluence residuals. The map is characterized mainly by high correlations (light copper color) transitioning into areas of high anticorrelations (black and dark brown colors), in total 65 of the 100 pairs show correlations with absolute values above 0.5.

Fig. 7 (a) and (b) cross-plots of fitting residuals of all Peak and Fluence spectra at 6 MeV and 115 MeV respectively, and (c) mapping of the correlation coefficients of all Peak residuals to all Fluence residuals at all energies.

The estimation of the ΔFluences is based solely on the fitting residuals of the Peak spectrum and the approach employed here is similar as before for the spectral parameters. Specifically, the ΔPeaks are used to calculate 10 quadratic hypersurfaces, one for the estimation of each of the ΔFluences at the 10 energy bins. Only the ΔPeaks at energies with absolute correlation values above 0.5 are used for the estimation of the ΔFluence at each energy bin. The form of the hypersurfaces is shown in equation (4):

(4)where y is the estimated Fluence residual at a specific energy bin, x _i are the Peak residuals, N is the number of energy bins where the Peak residuals have absolute correlations above 0.5 with the Fluence residuals at that specific energy bin, m and k are constants and d the intercept. The hypersurface for each energy is resolved separately. Figure 8 shows cross-plots of the data ΔFluences and the estimated ΔFluences at the first, third, sixth and ninth energy bins at 6 MeV, 12.6 MeV, 38 MeV and 115 MeV respectively. Figures for all energies can be seen in the Supplementary materials.

Fig. 8 Cross-plots of Fluence residuals derived from data and estimated from the hypersurfaces using the Peak residuals, y = x equality lines are in blue. (a) 6 MeV, (b) 12.6 MeV, (c) 38 MeV and (d) 115 MeV.

As mentioned above, the estimated residuals adjust the initial Fluence spectrum derived from the estimated spectral parameters [ a _F , b _F , c _F ]. This is done by using equation (5):

(5)where eFS is the estimated final fluence spectrum, eAS is the estimated analytical fluence spectrum, eRes is the estimated fluence residual, E is the energy and i is the energy bin index. Figure 9 shows such an example comparing the Fluence spectrum derived from the data, the Fluence spectrum calculated from the estimated [ a _F , b _F , c _F ] spectral parameters as well as the latter adjusted with the estimated residuals. It is seen how the addition of the residuals can improve the final profile of the spectrum estimation so that it shows better agreement with the actual values from the data.

These results illustrate that the Fluence spectral parameters as well as the fitting residuals can be well reproduced using as estimators the respective Peak spectral parameters and fitting residuals. However, these results do not guarantee that apart from estimation, this approach can be used for actual predictions, i.e., to provide reliable results from the parameters of Peak spectra that have not been included in the calculations of the hypersurfaces if they are used as predictive functions. Below, we validate our approach and demonstrate that it is indeed applicable for the actual prediction of Fluence spectra. Finally, we also note that very similar relationships where found for integral Peak fluxes and Fluences meaning that the methodology described here can also be applied to integral proton flux data.

3 Results

3.1 Fluence prediction and validation

The 266 Events in the SEPEM REL with their respective Peaks and Fluences are used to validate our approach. An iterative “leave-one-out” validation process is employed in which the Peak and Fluence parameters of the Event whose Fluence is to be predicted are not taken into account in the calculation of the predictive hypersurfaces. The process is detailed below and depicted in Figure 10.

1. The Peak and Fluence spectra of the 266 Events are fitted with the analytical functions to derive the [a b c ] triplets along with the respective fitting residuals.

2. All the variables of the Event whose Fluence is considered unknown are set apart and the variables from the rest 265 Events are used to calculate the hypersurfaces. Three hypersurfaces are calculated for the prediction of each of the Fluence spectral parameters [ a _F b _F c _F ] from the Peak spectral parameters, and 10 for the prediction of each of the Fluence residuals from the Peak residuals.

3. After the hypersurface calculation, the Peak parameters and Peak residuals of the Event that was set apart (with the unknown Fluence) are input to the hypersurfaces and a Fluence spectrum is predicted.

4. The process repeats iteratively for all 266 Events and in each iteration a different Event is set apart and its Fluence is considered unknown.

Fig. 9 Comparison of Fluence spectra, derived from data in black, calculated from the estimated spectral parameters in blue (analytical spectrum) and in red is the final estimation which is the analytical spectrum adjusted with the estimated residuals.

This process goes through the whole list and allows predictions to be made for each and every one Event. Since in every iteration the parameters from the Event under consideration do not participate in the hypersurface calculations the final Fluence outputs are true predictions, essentially as if the Event was the next one to occur. The Peak values are given only as input to the pre-derived hypersurfaces (they do not participate in their calculation) to predict the Fluence [ a _F b _F c _F ] parameters and residuals. Thus, the Peak spectrum is treated as if it was known or accurately predicted by some other means. Finally, the predicted Fluence spectra are compared to the actual measured Fluences. Figure 11 shows the mean and geometric mean Fluence spectrum from all the Events in the SEPEM REL along with the predictions in this work where the overall good agreement can be seen. The insets show the percentage error for the mean spectrum which is below 25% while the geometric mean spectra are virtually indistinguishable with the error being below 0.2%.

Fig. 10 The process used for the “leave-one-out” validation of the Fluence predictions.

In more detailed comparisons, Figure 12 shows cross-plots of measured and predicted Fluences at the first, fifth and ninth energy bins at 6 MeV, 26.3 MeV and 115 MeV, respectively. Additionally, the bar plots show the distribution of multiplicative factors by which the predictions differ from the measurements, being either higher or lower. Figures for all energies can be seen in the Supplementary materials.

Fig. 11 Comparison of (a) mean spectra and (b) geometric mean spectra of Fluences from all Events in the SEPEM REL and all predicted Fluences in this work. Insets show percentage error at all energy bins.

The cross-plots show the predicted Fluences agree well with the data mostly clustering along the y = x equality lines in blue and remaining almost always within a factor of 4 (red lines above and below). The bar plots show the distribution of the fraction of Fluences over predictions where 70% or more of the predictions are within a factor of 2 of the data at all energies and more than 95% are within a factor of 4. Table 3 shows the percentages for the 10 energy bins. These results demonstrate this methodology can reliably predict Fluence spectra from Peak spectra across the whole 5–200 MeV energy range. Additionally, the distributions shown in the bar plots can be treated as probability density functions (PDFs) of the Fluence predictions. Such PDFs can be calculated with a dense binning and implemented directly, or even fitted analytically, e.g., with Gaussians, in order to additionally derive confidence levels for the predictions.

Finally, the methodology and process described was also applied using the FORSPEF Event list with the RDSv2 dataset. This Event list contains 314 Events from 1984 to 2013 and it was defined from the cross-calibrated (with IMP-8) GOES proton data on the channel with effective energy of E _eff = 10 MeV. The 314 Peak and Fluence spectra were fitted, with the respective [a b c] triplets and residuals calculation, as well as the hypersurfaces derivation. We applied the methodology to an additional, and quite different, Event list in order to verify and demonstrate that no list specific characteristics affect the process or introduce some bias due to the way that Events are defined. Indeed, the results using the FORSPEF Event list are very similar to those using the SEPEM REL. This strongly indicates that our results are not list-specific and the Peak-Fluence relationships shown and used here are intrinsic and robust, meaning this approach can be employed regardless of the Event list used. All figures using the FORSPEF list with the same comparisons between predicted and measured Fluences with cross-plots and bar-plots, as in Figure 12, can be seen in the Supplementary materials. In Table 4 the percentages of predictions within a factor of 2 and 4 with the FORSPEF list are shown. We note that a key difference in the FORSPEF list is that, unlike in the SEPEM REL, multiple proton enhancements (sub-Events) are not compounded into a single Event, even if a new enhancement partially overlaps with a previous one. In this regard the FORSPEF list offers better statistics in terms of Peak-Fluence correlations at all energies (see Supplementary materials) as well as absolute number or Events with which to derive the predictive hypersurfaces. The deficiency of non-compounding Event definition is that the Fluence of a previous enhancement will be somewhat lower due to overlap. However, this is mitigated by the fact that in such cases the new enhancement usually occurs in the decay phase of the previous so less information is lost about the Fluence. These characteristics of the FORSPEF Event list are factors for the slightly higher prediction percentages yielded, especially in the interval.

As mentioned, an important aspect of these findings and this approach is that it can be used as an easy to implement additional tool or module in systems and models with methodologies that forecast SPE Peak fluxes. The only requirement is the calculation of the Peak-to-Fluence hypersurfaces from the SPE Event list used in each system. The hypersurfaces can then directly provide predictions and even confidence levels can be calculated from the resulting distributions such as those shown in the bar-plots in Figure 12.

3.2 Nowcasting of SPE Fluences

As discussed this approach can be applied in forecasting systems where a Peak spectrum is predicted ahead of the SPE. However, for nowcasting applications and systems the time window for a parameter prediction is of paramount importance, even excellent predictions are of little use if they come too late. In this regard, we investigate the applicability of our approach for SPE Fluence nowcasting, i.e., after an Event has been detected and near-real-time measurements are incoming, and we show that there is indeed practical potential. The process is similar to the validation process described above but it incorporates a per-measurement mode to account for the evolution of an SPE in time. The steps are analytically described below. For simplicity we consider a process where no peak forecast is available.

1. The fluence prediction process is initiated as soon as an SPE starts and proton flux intensity measurements become available.

2. At each time step the new flux intensity measurements are stored and the (so far) maximum values from all energy bins are considered as the (current) Peak spectrum.

3. At each time step when a new measurement becomes available the current Peak spectrum is fit with an exponential-cut-off power-law and its spectral parameters [ a _P , b _P , c _P ] are derived as well as the 10 residuals. The Peak spectrum is taken as the maximum flux at each energy channel (not necessarily coincident in time).

4. The Peak parameters and residuals are input in the pre-fitted hypersurfaces and the Fluence parameters and residuals are predicted. As in the previously described “leave-one-out” validation scheme, the hypersurfaces are pre-calculated without any of the variables of the “new” Event currently under consideration.

5. The Fluence spectrum is calculated and the current prediction at each time step is given.

As described the process is repeated with each new incoming measurement and this results in time-series of predictions for the Event Fluences. Figure 13 shows an example of the outputs of this process, Figure 13a shows the Event’s proton flux intensity time-series, Figure 13b shows the Fluence prediction time-series compared to the actual Fluences of the Event (dashed lines), and Figure 13c shows the calculated time-series of the ratio of Fluence over prediction. It can be seen that the prediction process starts with a big error as the initial measurements are quite far from the Peak values. However as the proton fluxes rise the prediction gets closer and closer to the actual Fluence values of the Event. After the SPE reaches peak values the prediction remains the same and is static, while the same is true for the ratio of Fluence over prediction. Finally, in Figure 14 the final predicted spectrum obtained at ~33 h after onset is shown for the sample Event in Figure 13.

Fig. 12 (a) Cross-plots of measured and predicted Event Fluences at 6 MeV, 26.3 MeV and 115 MeV, y = x equality lines are in blue and y = 4x and y = 0.25x lines are in red above and below respectively. (b) Distributions of the fraction of Fluence over predictions at the same energies.

Fig. 13 (a) Flux time-series of SPE, (b) Fluence prediction time-series are in solid lines and SPE data Fluences are shown in dashed lines for comparison and (c) time-series of Fluence over prediction ratio.

Naturally, for the application of this methodology in nowcasting the sooner an SPE reaches its maximum flux values the sooner a good prediction for the Fluence can be made. Therefore, we have performed a statistical analysis on the time occurrence of the Peaks in the SEPEM REL used here. Figure 15 shows the binned cumulative distributions of the time of Peak occurrence in hours, and as a fraction of the total Event duration at 6 MeV, 26.3 MeV and 115 MeV. More than 70% of the Events reach their Peak values at all energies within the first 30–50 h from the SPE onset and respectively more than 70% reach Peak values at all energies within the first 35%–50% of the total Event duration. However, we note here that this is a list-specific characteristic and it is dependent on the Event definition used. As previously discussed, the SEPEM REL contains multiple Events that are comprised of more than one enhancements, i.e., more than one overlapping sub-Events considered as one due to concurrency. They are not differentiated or separated, resulting in a number of Events having their Peak flux later in their duration, even though a sub-Event may have peaked earlier. In this sense, the SEPEM list is a worst-case scenario for this application. Contrary to this, since the FORSPEF list separates such multiple Events the statistics on the time of Fluence prediction are further improved. Figures with the distributions of Peak occurrence in hours, and as a fraction of the total Event duration after onset at all energies for the SEPEM REL and FORSPEF list can be seen in the Supplementary materials.

Fig. 14 Final Fluence prediction for the SPE shown in Figure 13 acquired at ~33 h after onset, compared to the actual SPE Fluence of the measured flux time-series.

Fig. 15 (a) Cumulative distributions of Peak occurrence in hours after onset at 6 MeV, 26.3 MeV and 115 MeV. (b) Cumulative distributions of Peak occurrence relative to total Event duration at the same energies. Red dashed lines denote the 70% limit.

These statistics show that our approach can indeed also be applied for nowcasting as it can provide accurate predictions of SPE Fluences within a reasonable timeframe, providing good predictions of Fluence for the majority of Events many hours or days before their end.

4 Summary

We have investigated in depth the relationships of Peak and Fluence spectra in SPEs using the SEPEM RDS v2.0 differential proton flux intensity dataset and two SPE Events lists, the SEPEM Reference Event List and the Event list used in the FORSPEF system. It is shown that within each Event the Peak spectrum exhibits a log-linear relationship with the respective Fluence spectrum across the range of proton energies of 5–200 MeV. By fitting the Peak and Fluence spectra with analytical exponential cut-off power-law functions it is found that this relationship is directly reflected in the close similarity of the spectral indexes and the scaling relationship between the spectral intensities. Also, it is shown that the spectral parameters are interdependent and exhibit high correlations among them. Furthermore the residuals from the analytical fitting are found to also exhibit high correlations/anticorrelations, meaning the Peak and Fluence spectra diverge from the strict analytical functions in similar ways. Based on these findings an approach is devised to predict a Fluence spectrum for a given Peak spectrum of an Event. This is done by predicting the Fluence spectral parameters as well as the residuals of the analytical fitting from the respective Peak spectral parameters and residuals. The approach is validated for all Events in a sound and consistent way and the prediction results are shown to agree well with the data across all energy bins. The predicted Fluences remain within a factor of 4 of the measured Fluences in more than 95% of all cases and at least 70% of the predictions are within a factor of 2 at all energies. Very similar results were also obtained using the FORSPEF Event list which strongly indicates that the results are not list-specific but based on intrinsic properties of the SPEs. This approach can be directly implemented in forecasting systems that already provide a prediction of the Peaks for an Event prior to its occurrence, or onset, in order to additionally derive a good Fluence prediction. Finally, we study the applicability of this approach for nowcasting systems which use real-time or near-real-time measurements of SPEs that are ongoing. It is shown that since the majority of the Events, more than 70%, reach their Peak values within 30–50 h from their onset our approach has application even in nowcasting systems as it can derive a reliable total Fluence spectrum prediction many hours or days before the end of the Event providing thus potentially useful and valuable information.

Acknowledgments

This work was supported by the European Space Agency (ESA) contract No. 4000120480/17/NL/LF/hh, “Solar Energetic Particle (SEP) Advanced Warning System (SAWS)”. The editor thanks two anonymous reviewers for their assistance in evaluating this paper.

Supplementary material

Figures S1 to S19 Access here

References

All Tables

All Figures

	Fig. 1 (a) Fluence-Peak pair showing strong log-linearity, (b) Fluence-Peak pair with larger residuals and (c) all Fluence-Peak pairs from the SEPEM list.
In the text

	Fig. 2 Slopes k versus constants m of the y = kx + m log-linear fits from all pairs of Peak and Fluence spectra.
In the text

	Fig. 3 (a) Peak spectrum derived from SEPEM RDS for an SPE Event and the analytical fit made with exponential cut-off power-law. (b) Same for respective Fluence spectrum. The fitted spectral parameters can be seen in the bottom left corners.
In the text

	Fig. 4 Cross-plots of the spectral parameters of Fluences and Peaks showing the relationships between them. (a) Spectrum intensity a , (b) power-law index b and (c) exponential index c .
In the text

	Fig. 5 Cross-plots of Peak spectral parameters for all 266 Events. (a) Spectrum intensity a P with power-law index b P and (b) power-law index b P with exponential index c P .
In the text

	Fig. 6 Cross-plots of Fluence spectral parameters derived from data and estimated from the estimator hypersurfaces, y = x equality lines are in blue. (a) Fluence spectrum intensity a F , (b) Fluence power-law index b F and (c) Fluence exponential index c F .
In the text

	Fig. 7 (a) and (b) cross-plots of fitting residuals of all Peak and Fluence spectra at 6 MeV and 115 MeV respectively, and (c) mapping of the correlation coefficients of all Peak residuals to all Fluence residuals at all energies.
In the text

	Fig. 8 Cross-plots of Fluence residuals derived from data and estimated from the hypersurfaces using the Peak residuals, y = x equality lines are in blue. (a) 6 MeV, (b) 12.6 MeV, (c) 38 MeV and (d) 115 MeV.
In the text

	Fig. 9 Comparison of Fluence spectra, derived from data in black, calculated from the estimated spectral parameters in blue (analytical spectrum) and in red is the final estimation which is the analytical spectrum adjusted with the estimated residuals.
In the text

	Fig. 10 The process used for the “leave-one-out” validation of the Fluence predictions.
In the text

	Fig. 11 Comparison of (a) mean spectra and (b) geometric mean spectra of Fluences from all Events in the SEPEM REL and all predicted Fluences in this work. Insets show percentage error at all energy bins.
In the text

	Fig. 12 (a) Cross-plots of measured and predicted Event Fluences at 6 MeV, 26.3 MeV and 115 MeV, y = x equality lines are in blue and y = 4x and y = 0.25x lines are in red above and below respectively. (b) Distributions of the fraction of Fluence over predictions at the same energies.
In the text

	Fig. 13 (a) Flux time-series of SPE, (b) Fluence prediction time-series are in solid lines and SPE data Fluences are shown in dashed lines for comparison and (c) time-series of Fluence over prediction ratio.
In the text

	Fig. 14 Final Fluence prediction for the SPE shown in Figure 13 acquired at ~33 h after onset, compared to the actual SPE Fluence of the measured flux time-series.
In the text

	Fig. 15 (a) Cumulative distributions of Peak occurrence in hours after onset at 6 MeV, 26.3 MeV and 115 MeV. (b) Cumulative distributions of Peak occurrence relative to total Event duration at the same energies. Red dashed lines denote the 70% limit.
In the text