Is there a direct solar proton impact on lower-stratospheric ozone?

We investigate Arctic polar atmospheric ozone responses to solar proton events (SPEs) using MLS (Microwave Limb Sounder) satellite measurements (2004–now) and WACCM-D (Whole Atmosphere Community Climate Model) simulations (1989–2012). Special focus is on lowerstratospheric (10–30 km) ozone depletion that has been proposed earlier based on superposed epoch analysis (SEA) of ozonesonde anomalies (up to 10 % ozone decrease at ∼ 20 km). SEA of the satellite dataset provides no solid evidence of any average SPE impact on the lower-stratospheric ozone, although at the mesospheric altitudes a statistically significant ozone depletion is present. In the individual case studies, we find only one potential case (January 2005) in which the lower-stratospheric ozone level was significantly decreased after the SPE onset (in both model simulation and MLS observation data). However, similar decreases could not be identified in other SPEs of similar or larger magnitude. Due to the input proton energy threshold of > 300 MeV, the WACCM-D model can only detect direct proton effects above 25 km, and simulation results before the Aura MLS era indicate no significant effect on the lower-stratospheric ozone. However, we find a very good overall consistency between WACCM-D simulations and MLS observations of SPE-driven ozone anomalies both on average and for the individual cases including January 2005.


Introduction
In the near-Earth space, solar wind charged particles are guided by the Earth's magnetic field and are able to precipitate into the middle and upper atmosphere in the polar re-gions. Such a kind of precipitation creates the spectacular aurora but also produces considerable amounts of HO x (H, OH, HO 2 ) and NO x (N, NO, NO 2 ) through ion-neutral chemistry (e.g., Verronen and Lehmann, 2013). HO x and NO x increases lead to ozone loss through catalytic reactions in the mesosphere and upper stratosphere, respectively (Sinnhuber et al., 2012). Moreover, in polar winter, NO x has a long chemical lifetime due to limited photodissociation by solar radiation. NO x produced by energetic particle precipitation (EPP) in the mesosphere to lower thermosphere is transported down to the stratosphere by the Brewer-Dobson circulation inside the polar vortex , causing depletion of upper-stratospheric ozone (Damiani et al., 2016). A number of studies have confirmed EPP's remarkable role in ozone depletion directly during large EPP events (e.g., Funke et al., 2011) and indirectly due to descending NO x (e.g., Randall et al., 2007). Thus, many advanced chemistryclimate models are now including EPP forcing, in order to correctly represent the ozone distribution in the polar stratosphere and mesosphere (Matthes et al., 2017;Stone et al., 2018). Solar proton events (SPEs) are one of the main types of EPP. During SPE, particles (mainly protons) with energies from tens to hundreds of megaelectronvolt (MeV) precipitate into the atmosphere at geomagnetic latitudes larger than 60 • for days. Such high-energy particles mainly affect the atmosphere at altitudes of 35-90 km, providing direct ionization forcing on the polar middle atmosphere. Large SPEs have been studied since the 1960s up to today using satellite observations and model simulation. In addition to tens of percent of ozone loss observed at altitudes above 35 km (Jackman et al., 2001;Seppälä et al., 2004;Verronen et al., 2006), a strong SPE can reduce total ozone by 1 %-3 % for months after the event Jackman et al., 2014).
Recently, Denton et al. (2018a, b) presented statistical studies of average ozone changes from 191 SPEs between 1989 and 2016 using ozonesonde measurements. Superposed epoch analysis of ozone anomalies at polar stations (Sodankylä, Ny-Ålesund and Lerwick) indicated that SPEs occurring during winter are causing ozone decrease by 5 %-10 %, on average, at 20 km altitude. This effect is not produced in the current models because SPE-induced ionization rates are insignificant at this altitude even during the largest events with high proton energies from 300 to 20 000 MeV . Denton et al. (2018a, b) included also a large number of very small SPEs in their analysis. Such ozone decreases have not been observed in the case studies of very extreme (particles with energies > 10 MeV are greater than 10 000 pfu, particle flux units) SPEs, e.g., the 2003 "Halloween" event, from either simulation or satellite observation (Funke et al., 2011, and references therein). Recently, statistical analysis based on simulations has found no evidence of such low-altitude ozone impact (Kalakoski et al., 2020). Moreover, from the chemical aspect, we also rather expect ozone increase in the lower stratosphere due to the enhanced NO x interfering with chlorine-driven catalytic ozone loss .
Here we investigate the proposed SPE-induced direct depletion on lower-stratospheric (10-30 km) ozone using ozone data from the Microwave Limb Sounder (MLS) instrument aboard the Aura satellite and the Whole Atmosphere Community Climate Model (WACCM-D) simulations. We proceed to evaluate ozone changes at altitudes of 10-70 km caused by SPEs both statistically (superposed epoch analysis) and individually (case by case). The MLS ozone data, WACCM-D atmospheric simulation and SPE datasets are presented in Sect. 2. In order to cross-check ozone depletion at 20 km reported based on the ozonesonde data, statistical ozone responses from MLS satellite measurements are firstly provided in Sect. 3. Following that, MLS and WACCM-D ozone changes after individual SPEs are given in Sect. 4. Finally, we summarize our results and conclusions in Sect. 5.

O 3 profile measurements by MLS
MLS onboard the Earth Observing System (EOS) Aura satellite measures ozone emission at 240 GHz, providing ozone volume mixing ratios at 55 pressure levels since 15 July 2004 (Waters et al., 2006). Vertical profiles are retrieved from the MLS observations every 165 km along the polar orbit at altitudes between 8 and 90 km, with a vertical resolution of ∼ 3.2 km. In this work, we use version 4.2 ozone data measured at 261-0.02 hPa (∼ 10-70 km) to calculate the daily averaged ozone density profile at northern high latitudes (60-90 • N). Readers who are interested in the MLS data quality are referred to Livesey et al. (2018).

O 3 from WACCM-D simulations
WACCM is a global circulation model, including fully coupled dynamics and chemistry. Here, we use version 4 of WACCM with resolution of 1.9 • latitude by 2.5 • longitude, with 88 vertical levels reaching from surface to 6 × 10 −6 hPa (≈ 140 km). Overview of the model and the description of climate and variability in long-term simulation was presented by Marsh et al. (2013), with details of model physics in MLT (mesosphere to lower thermosphere region) and the response of the model to radiative and geomagnetic forcing during solar maximum and minimum described by Marsh et al. (2007). The simulation results presented here are from WACCM-D, a variant of WACCM with a more detailed set of lower ionospheric chemical reactions, aiming at better reproduction of observed effects of EPP on MLT neutral composition Andersson et al., 2016).
We use SD-WACCM-D specified dynamics configuration, with Modern-Era Retrospective analysis for Research and Applications (MERRA) (Rienecker et al., 2011) meteorological fields to force dynamics at every time step up to about 50 km. Simulation covers years 1989-2012 and uses forcings from auroral electrons (E < 10 keV), solar protons (E < 300 MeV) and galactic cosmic rays for energetic particle precipitation. The SPE ionization rates are based on proton flux measurements from the Geostationary Operational Environmental Satellite (GOES) (see e.g., Jackman et al., 2011, for the calculation method). The WACCM-D SPE effects on neutral species are compared to satellite observations in Andersson et al. (2016). Note that WACCM-D has not been validated below 20 km. Nevertheless, in Andersson et al. (2016) the HNO 3 response above 15 km to single SPE onset was reasonable compared to MLS data. We also stress that protons with energies over 300 MeV are not included in the simulation. The 300 MeV protons mostly affect the atmosphere at around 25 km (Turunen et al., 2009). As 300 MeV is the upper limit of the proton energies considered in our model simulation, the WACCM-D simulation presented here can therefore only investigate the impact of direct proton forcing at altitudes above 25 km. For more details of the simulation setup, see Kalakoski et al. (2020).

Solar proton events
The data of solar proton events (SPEs) used in this study are based on NOAA GOES proton flux observations. Figure 1 presents 261 SPEs recorded from 1975 to date, including their onset time, fluxes detected in space, approximated time of duration and average ionization rates to the atmosphere at two altitudes. Here, the onset of a SPE is defined as the time when 5 min average proton fluxes with energies > 10 MeV are greater than 10 particle flux units (1 pfu = 1 particle/cm 2 /s/sr) at the geosynchronous orbit. For the estimation of SPE duration and its impact on the atmosphere, we use the daily average ion pair production rates at ∼ 1 hPa (∼ 46 km, upper panel of Fig. 1) and ∼ 12 hPa (∼ 29 km, lower panel of Fig. 1). These ionization rates are calculated from GOES proton flux observations using the energy deposition methodology described in, for example, Jackman et al. (2011). The SPE durations presented here were calculated as the period when the ionization rates at ∼ 1 and 12 hPa are larger than 2 ion pair/cm 3 /s before the next event starts. The average ionization rates in Fig. 1 were then derived by averaging the ionization rates at 1 and 12 hPa during this period. Our study used 49 events that occurred after the launch of Aura MLS (July 2004-now) and 177 events that occurred in the complete WACCM-D simulation period (January 1989-December 2012) to evaluate the ozone changes following SPEs. It is clearly demonstrated in Fig. 1 that these SPEs are more frequent near solar maximum years. A majority of the events are with flux less than 400 pfu, and their impacts to the atmosphere below 1 hPa are small. It is worth mentioning that although these SPEs seem to have no preference in occurring season, their seasonal distribution varies by months and should be considered during the interpretation (Fig. A1 in the Appendix).

Statistical O 3 response from MLS
Similar to the method used by Denton et al. (2018b), we applied a superposed epoch analysis to the MLS daily ozone anomalies. The superposed epoch analysis, also referred to as composite analysis in geophysics, is used to acquire vari-ation of a time series before and after an event or a chain of certain kinds of events. The point of time when the event begins is the epoch time. In this case, the epoch times are the onset times of individual SPEs during the MLS operating period. All available ozone data were binned as a function of epoch time and altitude, with temporal resolution of 1 d. The pre-and post-epoch spans used here are 30 and 60 d, respectively. For the selected sets of SPEs, all the binned ozone datasets were averaged to represent the effect of the SPEs. This method excludes natural ozone variations that are larger than the span scale. Since SPE-driven effects are expected to take place on daily to monthly timescales, variations caused by, for example, QBO (quasi-biennial oscillation) can be excluded. However, seasonal variations must be excluded before using superposed epochs. Thus, the daily profile climatology calculated from the ozone data was subtracted from the daily ozone data. Different from Denton et al. (2018a, b), to make sure SPEs are "isolated" from the previous events, events that happened within 10 d of the previous SPE were excluded.
In order to test the statistical significance of the obtained results, a Monte Carlo test was implemented. Instead of using SPE onset as epoch times, the analysis was rerun using 2000 random sets of epoch times. SPE-epoch averaged variations larger than 95 % of the 2000 randomized results are considered significant and robust (reported as > 95 % confidence), suggesting that these extracted signatures are likely not random but related to SPE or driven by some other external forcing.   Fig. 2a and c) and for the ones occurring in winter (November-April) (13 out of 19 events, Fig. 2b and d) within the instrument's operational period. Robust averaged anomalies (> 95 % confidence) are presented within the white thick lines. Spatial distribution of statistically robust anomalies is similar in all-SPE epochs and winter-SPE epochs. The depletion is more pronounced for winter epochs. This, of course, could be a statistical effect due to the much lower number of events used in the study but is also expected due to two facts: (1) ozone recovery is slower due to less production from O 2 photodissociation and (2) the largest SPEs with a flux > 1000 pfu that cause more ozone depletion to happen to occur in NH winter. Among all the SPEs during MLS measurement period, ∼ 3/4 of large SPEs are in NH wintertime (see Fig. A1). In both upper panels, closely following the SPE onset, very pronounced ozone depletion appears above 50 km for over 5 d. This is the direct ozone loss caused by the SPE-induced HO x enhancement. The number of extreme SPEs is relatively small, which explains the absence of the long-lasting ozone depletion that would be expected between 40 and 50 km from enhanced amounts of NO x . While the upper-stratospheric ozone depletion signature is not seen in the statistical average, a 5 %-10 % decrease of ozone is present below 30 km, including ozone loss around 20 km similar to that reported by Denton et al. (2018a, b). However, since this variation starts already several days before the epoch, we cannot exclude the possibility that the whole robust variation in the stratosphere is more related to other phenomena in the northern polar cap, e.g., to changes in the strength of polar vortex or related chemical effects. We will discuss this in more detail in Sect. 4.
A superposed epoch analysis of WACCM-D ozone anomalies from SPEs during 1989-2012 has been reported by Kalakoski et al. (2020); thus, we will not repeat it here. In their results, the epoch-averaged WACCM-D ozone anomalies showed the same robust depletion at above 50 km. Since their analysis included also the very large SPEs that occurred in 1989-2004 (see Fig. 1 in this study or the list of 15 largest SPEs in Table 1 in , long-term ozone depletion in the upper stratosphere was clearly detected as well. However, there was no robust ozone loss below 30 km found in the WACCM-D simulations.

O 3 response to individual SPEs
Considering the limited number of SPE events during MLS era and the high variability of stratospheric ozone influenced by, for example, sudden stratospheric warmings (SSWs) or heterogeneous chemistry on polar stratospheric cloud (PSC) surfaces, particularly during winter, in this section we analyze ozone responses to individual SPEs.
Similar to the analysis presented in Sect. 3, ozone anomalies presented here were calculated by subtracting daily climatology from daily averaged ozone data from MLS and WACCM-D. For Figs. 3 and 4, to make the results from MLS and WACCM-D simulation comparable, WACCM-D daily ozone was calculated using simulation profiles at MLS observation time and location. The climatology from MLS and WACCM-D were derived from their overlapping time period to guarantee a comparable background. For Figs. 6, A2 and A3, the subtracted daily mean climatology from MLS and WACCM-D were derived from the MLS data period and the WACCM-D simulation period, respectively. Then, instead of applying superposed epoch analysis on multiple SPEs, ozone anomalies are presented 30 d before and 60 d after onset of individual SPE. For estimating the statistical significance of the ozone anomalies found in the individual SPEs, we applied a similar Monte Carlo approach as in the case of SEA, i.e., the variance of 6000 random "onset" times was used as a measure for a significant anomaly. It is worth noting, however, that this method recognizes all "statistically significant" anomalies larger than the random background variation, whether the anomaly is due to SPE or, for instance, due to exceptional dynamical and chemical anomalies, which have a similar occurrence probability as SPEs.
Anomalies following all individual SPEs can be found in Figs. A2 and A3. In general, SPEs with proton fluxes < 400 pfu do not cause visible daily ozone depletion in the mesosphere (below 75 km) nor in other altitudes. Ozone changes following individual SPEs are more pronounced during winter. Figures 3 and 4  Ozone depletion is most pronounced following the January 2005 event. For this event, we also observe a robust lowerstratospheric ozone loss from MLS following SPE for the first time: ozone is depleted by ∼ 1 ppmv (∼ 15 %) at 20-35 km and by ∼ 0.15 ppmv (> 20 %) below 15 km 5 d after SPE onset.
Overall, wintertime ozone variation below 35 km is rather complicated. Year-to-year variability of stratospheric polar ozone is mostly controlled by dynamical and chemical processes; both are essentially coupled to temperature changes. Factors that modify polar temperature, e.g., sudden strato-spheric warming (SSW) and El Niño-Southern Oscillation (ENSO), are essentially planetary wave perturbations that modulate the strength of polar vortex. The probabilities of major SSWs and, on the other hand, springs with extremely strong polar vortex are at similar levels as the one of SPEs. Thus, ozone variations by these events will be seen as robust signatures in our study as well, yet they do not necessarily coincide with onsets of SPEs with proton fluxes > 400 pfu and > 10 000 pfu, respectively. The large SPE in January 2012 (Figs. 3 and 4) is severe enough to destroy stratospheric ozone. However, the stratospheric ozone anomalies at that time were dominated by dynamical ozone enhancement from SSW in 17 January 2012 . One of the most pronounced examples of extreme strong polar vortex impact is the well-reported ozone depletion during spring of 2011, which can be observed in ozone anomalies around the two small SPEs that occurred in March 2011 (see Fig. A2). The lower-stratospheric polar vortex was the strongest (in either hemisphere) in the previous 32 years (Manney et al., 2011). A large volume of polar stratospheric clouds (PSCs) converted chlorine reservoirs to ozone-destroying species, leading to extraordinary low ozone levels in the stratosphere (Pommereau et al., 2018). Similarly, robust anomaly seen after January 2016 SPE can be explained by cold 2015-2016 winter anomaly. We are confident to exclude SPE's influence on the anomaly in both cases because (firstly) the signal is not following SPE onset and (secondly) these SPEs are such small events that ozone loss was not observed -not even in the mesosphere. These robust non-SPE signals are included in the superposed epoch analysis performed in Sect. 3, contributing to the robust anomalies below 30 km in Fig. 2.
Identifying sources of the robust ozone anomaly below 35 km following the SPE beginning on 16 January 2005 is difficult. With a moderate cold winter temperature causing more ozone loss, coincident anomalies of robust dynamical ozone changes following the SPE exist. Meanwhile, an extremely large (over 270 %) ground-level enhancement (GLE) of neutrons occurred during the SPE period on 20 January 2005 . Ionization rate reached 500 cm −3 s −1 at 30 km for 1 d due to the very high energy protons (300-20 000 MeV) that caused the GLE (Usoskin et al., 2011 carried out a detailed study of January 2005 SPE's influence on the northern polar atmosphere using WACCM3 simulation, and they reported an ozone column decrease of less than 0.01 % by GLE protons, while the ozone changes below 50 km observed in MLS data were attributed to seasonal changes. The MLS ozone anomalies we observe are, on the contrary to the analysis in Jackman et al. (2011), not due to seasonal changes. To identify whether the anomalies are due to direct SPE effect or not, relative ozone response from MLS and WACCM-D simulation in MLS observation time and location to 16 January 2005 SPE are compared in Fig. 5. As WACCM-D simulations are carried out in the specified dynamics mode, any dynamical variations of ozone, including ozone chemistry, are expected  Fig. S4 were removed). The black and red thick line area demonstrates direct ozone anomalies with > 95 % confidence after the Monte Carlo test from WACCM-D and MLS data, respectively. to be reproduced by the model well. But any direct proton impacts below 25 km would not be reproduced at all since protons with energies > 300 MeV are not included in the model input. So significant differences between the model and MLS response might indicate a direct proton impact. As shown in Fig. 5, ozone responses below 20 km are very similar between results derived from these two data sources (5 d after SPE onset, close to the GLE event), indicating no significant proton effect. We do see some differences between 20 and 30 km, which might demonstrate a possible direct proton effect. However, we would like to point out that compared to MLS, WACCM-D holds a > 20 % overestimation of northern polar cap ozone below 30 km in January-April (see right panel of Figs. 5, S4 in the Supplement and Fig. 1 in Froidevaux et al., 2019). Such differences may implicate a transport-related issue in the model (Froidevaux et al., 2019) and therefore weaken our confidence to confirm the robust signal difference between MLS and WACCM-D at 20-30 km as the evidence of direct SPE impact. Readers who are interested in the ionization rate of this case are referred to Fig. A4.
Nonetheless, in our study the robust MLS ozone destruction signature in the lower stratosphere following the January 2005 SPE is unique not only when compared to other SPEs cases after 2004 but also when large and extreme SPEs before 2004 are included (see the WACCM-D simulation result presented in Fig. 6). Further research needs to be done to confirm the dynamical and chemical factors that led to ozone destruction below 35 km in January 2005.

Conclusions
Recent studies have reported observations of up to 10 % average decrease of lower-stratospheric ozone at 20 km altitude following solar proton events (SPEs). However, mechanisms which could cause such a large low-altitude impact are not clear. We used the Aura MLS satellite ozone datasets from 2004 to date and WACCM-D model simulations from 1989 to 2012 to analyze SPE-driven ozone changes. In our approach, stratospheric and mesospheric daily ozone anomalies (10-70 km) were examined over the epochs of SPEs by applying (1) a superposed epoch analysis (SEA) for all the cases and (2) a case-by-case analysis for individual events. Statistical significance of the anomalies found in the ozone levels was estimated by employing a Monte Carlo approach.
Arctic polar ozone destruction in the mesosphere and upper stratosphere can be directly observed from satellite measurement anomaly, when following SPEs in September-April with proton fluxes > 400 pfu and > 1000 pfu, respectively. We observe 5 %-10 % ozone destruction below 30 km altitude in MLS SEA results. However, the depletion appears before the epoch time, i.e., SPE onset. We argue that such lower-stratospheric ozone losses are rather caused by an unusually stable and strong polar vortex, together with sufficient ozone-depleting reservoirs of chlorine as confirmed by the case-by-case study. In the case-by-case study, we find a very good overall consistency between SPE-driven ozone anomalies derived from the WACCM-D model simulations and the Aura MLS data. Despite the fact that the model can only detect direct proton effects above 25 km due to the input proton energy threshold of 300 MeV, the good consistency enables us to generalize the study also to the SPEs before the Aura MLS era. From 1989 to date, robust lowerstratospheric ozone decrease after SPEs was observed only once in ozone anomaly, i.e., following the January 2005 SPE. Ozone was depleted by ∼ 1 ppmv (∼ 15 %) at 20-35 km and by ∼ 0.15 ppmv (> 20 %) below 15 km 5 d after SPE onset. We further investigated this case by comparing WACCM-D and MLS data. Since WACCM-D is not expected to observe direct SPE impact below 25 km, a consistent ozone depletion below 15 km demonstrated that direct SPE impact is less likely to be the reason for this robust ozone loss. The source of ozone loss above 20 km, however, is not fully confirmed. We stress that the January 2005 event was followed by a GLE, but the SPE was not the strongest on record by far. The exact mechanisms of the suggested lower stratosphere impact following this event are currently unclear. The simulation results indicate that even for the strongest SPEs in our record, there is no significant effect on the lower-stratospheric ozone as such.
Although it remains unclear to what degree the lower ozone decrease in January 2005 was caused by the SPE, and how much was due to other natural variability, we suspect that the observed statistically significant lower-stratospheric ozone impact is most likely by chance coincident with the SPE onset. We note that further research on the January 2005 SPE case is necessary to solidly confirm the EPP dynamical and chemical factors that led to ozone destruction below 35 km, but this is outside the scope of this paper. Data availability. MLS ozone data used in this study are available at https://mls.jpl.nasa.gov/products/o3_product.php (last access: 19 November 2020, Livesey et al., 2018). Proton fluxes and solar proton events are available from https://www.ngdc.noaa.gov/stp/ satellite/goes/index.html (last access: 19 November 2020, NOAA, 2019). Daily geomagnetic activity Ap-index used in Fig. 1 can be found at https://www.ngdc.noaa.gov/stp/GEOMAG/kp_ap.html (last access: 19 November 2020, NOAA, 2019). The SPE-induced ionization rate dataset is available at https://solarisheppa.geomar.de/ solarprotonfluxes (last access: 19 November 2020, Jackman, 2011).
Author contributions. JJ and AK formed the idea of the work. JJ performed the analysis and wrote the paper with contributions from NK, PTV and AK. MES provided the WACCM-D model data. All the authors had intensive discussions about the method and results during the research.
Competing interests. The authors declare that they have no conflict of interest.
Acknowledgements. The work of Antti Kero is funded by the Tenure Track Project in Radio Science at Sodankylä Geophysical Observatory/University of Oulu. We would like to thank the MLS ozone teams for providing the ozone data. This work was carried out as a part of the International Space Science Institute (ISSI) project "Space Weather Induced Direct Ionisation Effects On The Ozone Layer". We appreciate the fruitful discussion with the group members in this project. We would like to thank the editor and the two anonymous referees for dedicating their time to this work.
Financial support. This research has been supported by the Kvantum Institute for the Mesospheric Monitoring of Ozone (MeMO) project.
Review statement. This paper was edited by Gabriele Stiller and reviewed by two anonymous referees.