Research Article
What can the annual ^{10}Be solar activity reconstructions tell us about historic space weather?
^{1}
Department of Meteorology, University of Reading,
Reading, UK
^{2}
100 Mt. Jellore Lane,
Woodlands,
2575, Australia
^{*} Corresponding author: l.a.barnard@reading.ac.uk
Received:
11
August
2017
Accepted:
20
February
2018
Context: Cosmogenic isotopes provide useful estimates of past solar magnetic activity, constraining past space climate with reasonable uncertainty. Much less is known about past space weather conditions. Recent advances in the analysis of ^{10}Be by McCracken & Beer (2015, Sol Phys 290: 305–3069) (MB15) suggest that annually resolved ^{10}Be can be significantly affected by solar energetic particle (SEP) fluxes. This poses a problem, and presents an opportunity, as the accurate quantification of past solar magnetic activity requires the SEP effects to be determined and isolated, whilst doing so might provide a valuable record of past SEP fluxes.
Aims: We compare the MB15 reconstruction of the heliospheric magnetic field (HMF), with two independent estimates of the HMF derived from sunspot records and geomagnetic variability. We aim to quantify the differences between the HMF reconstructions, and speculate on the origin of these differences. We test whether the differences between the reconstructions appear to depend on known significant space weather events.
Methods: We analyse the distributions of the differences between the HMF reconstructions. We consider how the differences vary as a function of solar cycle phase, and, using a KolmogorovSmirnov test, we compare the distributions under the two conditions of whether or not large space weather events were known to have occurred.
Results: We find that the MB15 reconstructions are generally marginally smaller in magnitude than the sunspot and geomagnetic HMF reconstructions. This bias varies as a function of solar cycle phase, and is largest in the declining phase of the solar cycle. We find that MB15's excision of the years with very large ground level enhancement (GLE) improves the agreement of the ^{10}Be HMF estimate with the sunspot and geomagnetic reconstructions. We find no statistical evidence that GLEs, in general, affect the MB15 reconstruction, but this analysis is limited by having too few samples. We do find evidence that the MB15 reconstructions appear statistically different in years with great geomagnetic storms.
Key words: solar activity reconstructions / cosmogenic isotopes / solar energetic particles
© L. Barnard et al., Published by EDP Sciences 2018
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Records of solar magnetic activity are critical for further developing our understanding of the Sun's magnetic field, as well as of the resulting space weather and space climate phenomena driven by solar magnetic activity. Past solar magnetic activity can be estimated from a range of different proxies, with varying reliability, temporal resolution, and extent into the past. Three well established proxies are found in sunspot observations, geomagnetic activity records, and the concentrations of cosmogenic radionuclides stored in natural archives (Owens et al., 2016a,b). These estimates of historical solar variability are invaluable for placing the comprehensive spaceage solar observations in a proper context (Barnard et al., 2017), and for investigating plausible scenarios of future space climate change (Abreu et al., 2008; Barnard et al., 2011; Steinhilber et al., 2012; Matthes et al., 2016).
Both the sunspot observations and geomagnetic activity records rely on a network of human observers and so can provide information on the Sun as far back as reliable measurements have been recorded and preserved (Lockwood, 2013; Usoskin, 2013). In practice these sunspot and geomagnetic activity records provide information on solar variability extending back several centuries, albeit with increasing uncertainty into the past, owing to changes in instrumentation, and observing procedures, and problems of intercalibrating the results of different observers.
A different kind of proxy is available in the form of the concentrations of cosmogenic radionuclides stored in natural archives, such as ^{14}C fixed in tree trunks, and ^{10}Be stored in polar ice sheets. These records do not rely on the past efforts of an observer network, but instead require present day retrievals of environmental samples which require significant levels of expertise and expense to process and interpret. For example, converting the deposited cosomogenic iosotope concentrations into the cosmogenic isotope production rate throughout Earth's atmosphere requires a detailed understanding of the transport and deposition processes of that cosmogenic isotope through the Earth system (McCracken, 2004; Heikkilä et al., 2007). Additionally, the inversion procedure needed to understand the implications for solar activity requires a model of galactic cosmic ray (GCR) transport through the heliosphere, and estimates of the local interstellar spectrum of GCRs. However, using suitable inversion techniques with cosmogenic radionuclides it is possible to infer solar variability over several millenia. Although, due to the timescales of transport and deposition of cosmogenic radionuclides through the Earth system, this extended view into the past is at lower temporal resolution and higher uncertainty than is achieved with the sunspot and geomagnetic reconstructions (Lockwood, 2013; Usoskin, 2013). Recently annual records of solar activity extending back several centuries have been developed from cosmogenic isotope records, for example, from ^{10}Be by McCracken & Beer (2015).
Given the complexity and independence of these different proxy records, they each illustrate a remarkably consistent picture of longterm space climate variability over the past few centuries (Owens et al., 2016a,b). However, there is much more uncertainty about the historical space weather environment. Geomagnetic observations provide a good quantitative record of geomagnetic storms back to 1868 (Kilpua et al., 2015; Lefèvre et al., 2016; Vennerstrom et al., 2016) and increasingly sparse and qualitative information on geomagnetic storms back to approximately 1845 (Lockwood, 2013).
Much less is known about the historic occurrence of solar energetic particle (SEP) events and ground level enhancements (GLEs), with the first observation of GLEs in ionisation chambers in 1942 (Forbush, 1946). It was previously thought that enhanced atmospheric nitrate production due to SEP fluxes may leave a measurable increase in polar ice nitrate concentrations, and that consequently impulsive enhancements in polar ice nitrates could serve as a proxy for historic SEP events (Dreschhoff & Zeller, 1990; Shea et al., 1999; McCracken et al., 2001; Kepko et al., 2009). Subsequent investigations questioned the reliability of this potential proxy record (Wolff et al., 2012; Duderstadt et al., 2014), and currently it is widely held that polar ice nitrate concentrations are not a valid proxy for SEP fluxes.
It has also been argued, both through observations and physical modelling, that extreme SEP events may measurably perturb the ^{14}C and ^{10}Be cosmogenic isotope records (Usoskin et al., 2006; Miyake et al., 2012, 2013; Usoskin et al., 2013; Herbst et al., 2015; McCracken & Beer, 2015; Dee et al., 2016). Therefore, cosmogenic isotopes may provide the only viable proxy record for historic SEP fluxes. Such a proxy would be of significant value, as the effective mitigation of space weather hazards requires the development of relevant extreme event scenarios, and it may be a poor approximation to extrapolate our experience of current space weather into the future, when solar activity could be significantly different. In addition, the detection, accurate quantification, and removal of SEPgenerated cosmogenic isotopes would generate a better record of past variations in GCR fluxes.
We aim to assess if there is evidence for an SEP induced bias in the annual ^{10}Be reconstruction of McCracken & Beer (2015) (hereafter MB15). MB15 provided two annual reconstructions of the average heliospheric magnetic field (HMF), based on the analysis of ^{10}Be. The first assumed no ^{10}Be production from SEP events, while the second suppressed impulsive enhancements in ^{10}Be thought to be due to GLEs and atmospheric nuclear bomb tests. The motivation for this is two fold; firstly, to develop a better understanding of uncertainties in the MB15 reconstruction and assess how effective the SEP excision from ^{10}Be record was; secondly, if a bias is present, can it be used to infer anything about historical SEP fluxes, and does it motivate the pursuit of a next generation of reconstructions with improved temporal resolution and/or sensitivity?
2 Data
2.1 Cosmogenic isotope HMF reconstruction
McCracken & Beer (2015) produced annually resolved estimates of the nearEarth HMF magnitude from 1391–1983, derived from the analysis of ^{10}Be concentrations provided by the Dye3 icecore and North Greenland Ice Core project (NGRIP). Furthermore, MB15 presented evidence of impulsive enhancements in ^{10}Be concentrations following known GLE events (in 1942, 1949, and 1956), as well as after highaltitude nuclear bomb tests (in 1962). By analysing the magnitude of the impulsive ^{10}Be enhancements associated with the solar events, MB15 conceived an algorithm to identify other years in the ^{10}Be record which may have been affected by extreme space weather, over the period 1800–1983. We note that MB15 also analysed historic great geomagnetic storm (GGMS) occurrence, which provided additional circumstantial evidence that the years of impulsive ^{10}Be enhancements were associated with extreme space weather events. Consequently, MB15 considered it prudent to provide two estimates of the HMF; the first, B_{c}_{1}, was based on the original ^{10}Be record, without further consideration of the impulsive enhancements assumed to be associated with extreme space weather; the second, B_{c}_{2}, aimed to suppress the effect of the impulsive ^{10}Be enhancements by replacing the data in these years with interpolated estimates based on surrounding years. Here we use both the B_{c1} and B_{c2} records, over the period 1845–1983.
2.2 Geomagnetic HMF reconstruction
Over the last two decades, much progress has been made in estimating the nearEarth solar wind and HMF conditions from geomagnetic variability (Lockwood et al., 1999; Svalgaard & Cliver, 2005, 2010; Lockwood, 2013; Lockwood et al., 2013a,b, Lockwood et al., 2014a,b; Holappa et al., 2014a,b; Svalgaard, 2014; Mursula et al., 2016). Typically these techniques work by establishing a statistical relationship between the observed geomagnetic variability indices and nearEarth solar wind conditions. These relationships are calibrated over the space age, for which we generally have insitu spacecraft observations of the nearEarth solar wind. Then, to estimate historic solar wind conditions, the statistical relationships are extrapolated, applying them to geomagnetic observations from before when insitu observations are available. For technical reasons, related mainly to factors such as the shortterm variability in the HMF orientation and the variability of the spatial distribution of ionospheric conductivities, it has so far only been possible to estimate annual means of the HMF magnitude from geomagnetic indices (Lockwood, 2013; Lockwood et al., 2016b). Focusing on annual reconstructions of the HMF intensity, Owens et al. (2016a) recently reviewed the progress in this area, and published a consensus estimate of annual mean HMF computed from geomagnetic indices. The consensus was computed from a weighted mean of the estimates provided by Lockwood et al. (2013b) and Svalgaard (2014), and spans the period 1845–2013. We use the Owens et al. (2016a) consensus geomagnetic reconstruction, B_{g}, in our comparison with the B_{c1} and B_{c2} records.
2.3 Sunspot HMF reconstruction
The annual mean nearEarth HMF magnitude can also be estimated from sunspot observations. One approach to this is similar in nature to the geomagnetic reconstructions, which is to establish a statistical relationship between the observed HMF and sunspot number. For example, it has been shown that an approximately linear relationship exists between the annual mean HMF magnitude and the square root of the annual mean sunspot number (Svalgaard & Cliver, 2005). Another approach is to model the open solar magnetic flux (OSF), the source of HMF, as a continuity equation (Vieira & Solanki, 2010; Owens & Lockwood, 2012). For example, in the Owens & Lockwood (2012) model, the source of the OSF is parameterised in terms of the sunspot number, while the loss is parameterised in terms of the heliospheric current sheet tilt. With an estimate of the OSF, an empirical relationship between the OSF and HMF is used to compute the HMF. An additional consideration with these reconstructions is which sunspot record to use? Both the international and group sunspot numbers have been significantly revised recently, and multiple different methodologies and corrections have been proposed, with varying degrees of acceptance (Clette & Lefèvre, 2016; Cliver, 2016; Lockwood et al., 2016a). Owens et al. (2016a) also produced a consensus HMF reconstruction from sunspot observations, spanning the period 1750–2013. Owens et al. (2016a) used both the Svalgaard & Cliver (2005) and Owens & Lockwood (2012) methodologies with four records of the sunspot number (2 international and 2 group numbers) to compute the consensus HMF reconstruction, B_{r}, which we use here. The aim of using both reconstruction methods with multiple sunspot records was to help compute the uncertainty in HMF estimate from both the sunspot records and HMF reconstruction methods.
2.4 Ground level enhancement catalogue
We use a catalogue of GLEs provided by the National Geophysical Data Center, available at https://goo.gl/ztktlZ. This list is compiled from both the world wide neutron monitor network, whose spatial coverage varies in time, and also a sparser network of neutron monitors and ionisation chambers from before the development of the world wide network. Here we use this record to determine years with and without GLEs, over the period 1942–1983. We note that the sparsity of neutron monitors both before and during the development of the world wide network means that there is a chance, particularly earlier in the catalogue, that GLEs were missed. Note also that there is a >500fold dynamic range between the smallest (<10%) and largest (5000%) GLE in the catalogue, and we may expect that the contributions of the smaller (<200%) GLE will be proportionally smaller than those considered by MB15.
2.5 AA geomagnetic index
We use the aa geomagnetic index to analyse the occurrence of geomagnetic storms. We choose to use the aaindex as it is the longest continuously available index that is also sensitive to geomagnetic storms. The aaindex is a 3 h range index, which, at any given time, is computed from geomagnetic variability observed at two approximately antipodal observatories, one in the UK and one in Australia. However, both the UK and Australian observatories have moved twice over the duration of the aaindex, and so in fact the aaindex is a composite index formed from 6 pairings of crosscalibrated observatories. We obtained the 3hourly aaindex data for the period 1868present from http://isgi.unistra.fr/indices_aa.php.
3 Methods
3.1 Fractional differences
To compare the different HMF reconstructions we compute the fractional differences of the form (1) where B_{x} refers to either the geomagnetic or sunspot reconstructions, B_{g} or B_{r}, and B_{cn} refers to either B_{c1} or B_{c2}, the MB15 reconstructions that have not or have been modified to account for impulsive ^{10}Be enhancements.
3.2 Solar cycle phase and polarity cycle phase calculations
We will compare how the differences between the HMF reconstruction methods vary as a function of both the Schwabe solar cycle phase, and also the polarity cycle phase. By polarity cycle, we mean the cycle in the modulus of the dipole component of the solar magnetic moment. The solar cycle phase is computed according to (2) where t is time in years, and t_{start} and t_{stop} are the start and end dates of the Schwabe cycle, taken from https://en.wikipedia.org/wiki/List_of_solar_cycles. The polarity cycle phase is computed according to (3)where t_{n, max} and t_{n+1, max} are the times of the sunpot maxima in the Schwabe cycle, for the nth and the nth + 1 cycles respectively. The −1 in the numerator accounts for the fact that, on average, the solar polar fields reverse approximately 1 year after sunspot maximum (Lockwood et al., 2012). Of course, for the more recent cycles we could use magnetogram data to estimate polarity cycle timings more accurately. But as this is not possible for the entire record, we choose to use this more approximate method, which we can apply to the whole data sequence uniformly.
3.3 Geomagnetic storm definition
To define geomagnetic storms in the aaindex, we follow a procedure similar to Kilpua et al. (2015). We define storms in the aaindex as periods where aa exceeds 100 nT. Storm start and end times are set to the times of the nearest 3hourly value above the 100 nT threshold (approximately the 99th percentile of the full aa distribution). We also compute the peak aaindex value in each storm, and integrated intensity. Secondary storm maxima with no intermediate drop below the storm threshold are not split into separate events, and there is also no minimum time threshold between consecutive storms. This procedure generates 1919 storms over the period 1868–1983. Finally, we are interested particularly in a record of the largest storms, so we select from this list only storms with maximum intensities greater than the 90th percentile of storm maximum intensities. This results in a list of 186 large geomagnetic storms, which we hereafter refer to as GGMS. The 90th percentile threshold is arbitrary, but the results we later discuss based on this list of GGMS are robust to moderate changes to this threshold (between the 85 and 95th percentiles).
4 Results
Panels A–D of Figure 1 show the time series of the four HMF reconstructions, with B_{c1}, B_{c2}, B_{g}, and B_{r} in panels A, B, C and D, respectively. This window is restricted to the period of overlap between the four reconstructions, 1845–1983. The red squares mark the years in B_{c1} that are affected by impulsive ^{10}Be enhancements (McCracken & Beer, 2015). Some of these points clearly appear to be outliers, being widely separated from the local trend, more than uncertainty in the values would appear to explain.
Each of the HMF reconstructions shows similar Schwabe cycle variations, although these are better resolved, with smaller uncertainties for B_{g} and B_{r}. Note that this is not a criticism of the B_{c1} and B_{c2} series; considering the very different nature of the B_{c1} and B_{c2} series, we consider it remarkable that they agree with B_{g} and B_{r} so closely.
Panel E of Figure 1 shows kernel density estimates of the distribution of each of the HMF reconstructions, using the same color scheme as panels A–D. Kernel density estimation provides an empirical estimate of a parameter's probability distribution function by assigning a specified analytical distribution function (e.g. Gaussian, tophat, cosine etc.) to each observation, and assuming that the normalised sum of these individual distributions provides a fair representation of the parameter's true distribution function (Wilks, 1995).
Overall, the distributions of these four reconstructions appear quite similar, with similar modal values. Both the B_{c1} and B_{c2} series appear to be located at slightly smaller HMF magnitudes (≈0.1 nT); the median values of B_{c1} and B_{c2} are 6.08 nT and 6.16 nT, while the median values of B_{r} and B_{g} are 6.27 nT and 6.30 nT. The largest differences between B_{c1} and B_{c2} and B_{g} and B_{r} are found in the lower tail of the distributions, where there is more probability at lower values in the B_{c1} and B_{c2} series. Conversely, the upper tails of the four distributions appear quite similar. The impact of the modification applied to B_{c2} to account for the ^{10}Be impulsive enhancements is apparent, with the low tail having less probability at lower HMF values than for B_{c1}.
Figure 2 presents the time series of the fractional differences between the HMF reconstructions, F_{g1}, F_{g2}, F_{r1}, and F_{r2} in panels A–D, respectively. Each time series appears to show cyclic variations which correspond approximately with Schwabe solar cycles, which will be analysed in more detail in Section 4.1. Furthermore, we note the close correspondence between F_{g1} and F_{r1} series, as well as the F_{g2} and F_{r2} series, which arises due to the remarkable agreement between B_{g} and B_{r}. Finally, we note that, visually, there does not appear to be any secular trend in the fractional difference time series.
Kernel density estimates of the distributions of the fractional differences are shown in Figure 3. The distributions for F_{g1} and F_{g2}, and F_{r1} and F_{r2} are shown in panels A and B, respectively. The correspondingly colored dots mark the values of fractional differences in the years of impulsive ^{10}Be enhancements, where the vertical coordinates are random numbers, to make clear the spread of the values. For F_{g2} and F_{r2}, these points are clustered closer to the zerodifference line than for F_{g1} and F_{r1}, showing that in these years B_{c2} is in improved agreement with B_{r} and B_{g} relative to B_{c1}. This demonstrates the impact of the McCracken & Beer (2015) procedure to suppress the effects of the impulsive ^{10}Be enhancements. The four distributions are approximately centered on zero, with similar modal values of 0.01, 0.08, 0.00, and 0.07 for F_{g1}, F_{r1}, F_{g2} and F_{r2}, respectively. The distributions are moderately asymmetric, with each showing more weight below the modal value than above. It is not clear what causes this asymmetry, but it could possibly result from extra ^{10}Be production due to SEP fluxes. Were extra ^{10}Be to be produced by SEPs, that was assumed to be due to GCRs, it would lead to a underestimate of the HMF magnitude in the MB15 reconstructions. Such a negative bias would increase the occurrence of negative fractional differences, and so might explain the asymmetry shown in Figure 3.
Fig. 1 (A–D) Time series of the annual HMF reconstructions from 1845–1983, derived from; ^{10}Be concentrations, B_{c1} (with impulsive enhancements), and B_{c2} (with suppressed impulsive enhancements); geomagnetic observations, B_{g}; and sunspot observations, B_{r}. Error bars are 1σ errors of the mean. The red squares in panel (A) mark years where B_{c1} ≠ B_{c2}. (E) Kernel density estimates of the distributions of four different HMF reconstructions. 
Fig 2 Time series of the fractional differences between the B_{c1} and B_{c2} reconstructions with both B_{g} and B_{r}. Computation of the fractional differences is explained in Section 3.1: F_{g1} is the fractional difference between B_{c1} and B_{g}, F_{g2} between B_{c2} and B_{g}, F_{r1} between B_{c1} and B_{r}, and F_{r2} between B_{c2} and B_{r}. Error bars are 1σ errors of the mean. 
Fig. 3 (A) Kernel density estimates of the distributions of fractional differences between B_{g} with B_{c1} and B_{c2}, F_{g1} and F_{g2}, respectively. The purple and orange dots show the fractional differences in years (McCracken & Beer, 2015) identified as being affected by SEP events (the vertical coordinates of these points is randomised, which otherwise overlap). (B) The same as panel (A), except showing the fractional differences of B_{r} with B_{c1} and B_{c2}. 
4.1 Differences with Schwabe cycle phase and polarity cycle phase
Figure 2 revealed that the fractional differences between ^{10}Be reconstructions with the geomagnetic and sunspot reconstructions appear to show regular variations with the solar cycle. Here we investigate this further by looking at composites of the HMF and fractional differences over the solar cycle. Each year of the reconstructions is assigned a solar cycle phase, as described in Section 3.2. We divide the solar cycle into 10 equally spaced and contiguous phase bins, and for each bin, we compute the mean of each HMF reconstruction, as well as the mean of the fractional differences F_{g1}, F_{r1}, F_{g2} and F_{r2}. The overlapping extent of the 4 reconstructions spans solar cycles 9–21, although neither the beginning of cycle 9 nor the end of 21 are fully sampled. Because of this, and because of the variable length of the Schwabe cycle, each phase bin is not equally sampled, and minimum number of samples is 11 annual values, the maximum is 16 values, while 6 of the 10 bins have 13 or 14 samples. However, although the phase sampling is not equal, the differences are modest, and each bin has a reasonable number of samples, so we assume that this doesn't influence our subsequent analysis. Figure 4 presents the results of the solar cycle composites, with the HMF, F_{g1} and F_{r1}, and F_{g2} and F_{r2}, shown in panels A, B, and C, respectively. In each panel, the error bars are 1σ errors of the mean computed in each phase bin. This confirms that there is cyclic variability in the differences between the ^{10}Be reconstruction and the geomagnetic and sunspot reconstructions, which peaks in the declining phase of the solar cycle, around φ_{s} = 4 rads. Furthermore, the close agreement between the independently derived B_{g} and B_{r} series (panel A), suggests to us that these differences arise from the ^{10}Be reconstruction. Additionally, as they are present for both B_{c1} and B_{c2}, with nearly identical magnitudes, we conclude that these differences and are not caused or resolved by the procedure to suppress the impulsive ^{10}Be enhancements. This result is robust to changes in the number of solar cycle phase bins, with qualitatively similar results obtained with 8–12 phase bins.
The question then, is which processes could cause the systematic disagreement between the ^{10}Be reconstructions with the geomagnetic and sunspot reconstructions in the declining phase of the Schwabe cycle? One plausible candidate is the modelling of the GCR heliospheric transport, results from which are employed by McCracken & Beer (2015) to estimate the HMF from the ^{10}Be observations. GCRs propagate diffusively through the turbulent HMF to nearEarth space, but they also undergo systematic drift patterns due to large scale structure and polarity of the HMF. Full details of the MB15 inversion procedure and GCR transport model are provided in McCracken & Beer (2015). However, to summarise, in their description of the ^{10}Be inversion procedure, McCracken & Beer (2015) detail how the systematic effects of the Hale cycle modulations remain in the HMF series, as the inversion procedure necessarily averaged over the positive and negative phases of the Hale cycle. These polarity cycle effects are largest in the declining phase of the Schwabe cycle, as the new polarity cycle commences with the polarity reversal of the solar polar magnetic field.
Compositing over Schwabe cycles aliases the effects of the differing heliospheric GCR modulations with polarity cycle phase (Thomas et al., 2015). This aliasing can be avoided by instead compositing over both positive and negative polarity cycles, which we do here following a similar procedure to Owens et al. (2015). Following the standard nomenclature, we here on refer to the positive and negative polarity cycles as qA > 0 and qA < 0, respectively (where q is the charge polarity of the GCR and A is the solar field polarity with A > 0 and A < 0 describing field that is towards the solar poles in, respectively, the southern and northern hemisphere). As with the composite analysis presented in Figure 4, we split the positive and negative polarity cycles into 10 phase bins, and compute the means of the HMF and the fractional differences in each phase bin. Subdividing the data further into positive and negative polarity cycles means that the sampling in each phase bin is poorer than for the Schwabe cycle analysis, but each phase bin still has between 6 and 9 samples. The results of this compositing are shown in Figure 5, with the composites of HMF, F_{g1} and F_{r1}, and F_{g2} and F_{r2}, shown in panels A, B, and C, respectively. Here the largest systematic difference between ^{10}Be reconstructions and the geomagnetic and sunspot reconstructions appears in qA > 0 cycles.
Our interpretation of the systematic difference between the MB15 reconstructions with the geomagnetic and sunspot reconstructions is that this is primarily an artefact of the GCR heliospheric transport modelling used in the MB15 ^{10}Be inversion. Furthermore, this analysis shows that the systematic difference doesn't persist over the whole solar cycle, but is apparent in the declining phase of the Schwabe cycle, and can lead to differences of approximately 10–20% relative to the geomagnetic and sunspot reconstructions. Recent observations of the Hale cycle modulation of the GCR intensity should make it possible to use a GCR heliospheric transport model that includes particle drift effects (e.g. Strauss & Effenberger 2017) to examine long term trends, possibly reducing the differences seen in Figures 4 and 5.
Fig. 4 (A) Composite of the HMF reconstructions with Schwabe solar cycle phase. (B) Composite of the fractional differences of B_{g} and B_{r} with B_{C1}. (C) Composite of the fractional differences of B_{g} and B_{r} with B_{C2}. For each panel, the error bars are 1σ errors of the mean. In panels B and C, the black dashed line shows the zero fraction difference line. 
Fig. 5 Panels A–F show composites of the HMF reconstructions and the fractional differences F_{g1}, F_{r1}, F_{g2}, F_{r2}, with polarity cycle phase. The composites are split according to solar magnetic polarity, with qA > 0 cycles in panels A–C, and qA < 0 cycles in panels D–F. For each panel, the error bars are 1σ errors of the mean. 
4.2 HMF reconstructions differences conditional on large space weather event occurrence
In this section we will analyse the distributions of the reconstruction differences F_{g1} and F_{r1}, conditional on whether or not large space weather events were know to occur. As we wish to assess the possible sensitivity of the ^{10}Be record to significant space weather events, we exclude F_{g2} and F_{r2} from further analysis, as the B_{c2} series was processed to minimise the effects of known GLE events.
4.2.1 Historical space weather
As we look further into the past, particularly prior to the space age, reliable records of space weather events become harder to obtain (Barnard et al., 2017). Furthermore, the past development of observatories and instrumentation often means that it is not possible to obtain homogeneous quantitative records of space weather occurrence over long time periods. Here we will assess the historic occurrence of large space weather events in terms of the occurrence of GLEs (Sect. 2.4) and GGMS defined in the aa geomagnetic index (Sects. 2.5 and 3.3). Given the hypothesis that strong SEP events and GLEs may bias the cosmogenic isotope reconstructions, we would, ideally, only consider records of SEPs and GLEs. However, the GLE record only extends back to 1942, having a short 41 year overlap with the MB15 reconstructions, which end in 1983. Therefore we will also consider GGMS, the record of which extends back to 1868, having 115 years of overlap. In doing so we make an assumption that in years where GGMS occur, it is more likely that strong SEP events and GLEs will also occur. MB15 provided some circumstantial evidence that supports this assumption, showing that the majority of the 3σ increases in ^{10}Be before 1942 were associated with GGMS.
Figure 6 summarises these two metrics of historic space weather. Panel A shows a time series of the 3hourly aa index (black), with the maxima of GGMS marked by red circles. Furthermore, the timings of known GLEs are marked by the vertical cyan lines. Panels B and C shows these same data, but as a function of Schwabe cycle phase and polarity cycle phase. Additionally, Figure 7 presents an alternative view of these data. Panel A shows the annual count of GGMS and GLEs, while panels B and C show the count of GGMS and GLES with both Schwabe cycle phase and polarity cycle phase. For the later two, the Schwabe and polarity cycles were split in to 10 equal width phase bins, and the number of events in each bin was counted. This figure is intended to show the general trend of the occurrence of GGMS and GLE in time and with solar cycle phase. Note that we do not take account of the slightly uneven sampling of φ_{s} and φ_{p}, nor do we scale the counts to account for the 115 year duration of the GGMS record relative to the 41 year duration of the GLE record; the GGMS and GLE counts should therefore be compared in waveform, but not in amplitude. Taken together, Figures 6 and 7 reveal that there is cycletocycle variability in the occurrence of GGMS and GLE, but that typically these events are more likely to occur near the middle of the Schwabe cycle. A similar result was found for the occurrence of gradual SEP events (Barnard & Lockwood, 2011).
Fig. 6 (A) Time series of the aa geomagnetic index (black), from 1867–2014, with the maxima of great geomagnetic storms marked by red circles. Also shown are the timings of known GLE events (cyan lines). Panel (B) Shows the same data as (A), except as a function of Schwabe cycle phase. Panel (C) Shows the same data as (A), except as a function of polarity cycle phase. We note that, due to the availability of ^{10}Be data, we only use the geomagnetic and GLE observations up until 1983. 
Fig. 7 (A) Time series of the annual count of great geomagnetic storms (red), and GLE events (cyan). (B) The number of great geomagnetic storms (red) and GLE events (cyan) as a function of Schwabe cycle phase. (C) The number of geomagnetic storms (red) and GLE events (cyan) as a function of polarity cycle phase. In panels (B) and (C), phase bins 2π/10 wide were used. No account was taken of the uneven sampling of solar cycle phase, or the fact the GLE record only spans 1937present; these plots only show coarse trends with solar/polarity cycle phase, and the absolute number of geomagnetic storms and GLE events should not be directly compared. We note that, due to the availability of ^{10}Be data, we only use the geomagnetic and GLE observations up until 1983. 
4.2.2 HMF reconstruction differences with GLEs
To discern whether GLEs, in general, measurably affect the MB15 reconstruction we will use a statistical approach that compares the distributions of F_{g1} and F_{r1} in years with and without GLE events. We are restricted to the 41 year overlap period between the MB15 reconstructions and the GLE record, for which 19 years contained a GLE and 22 years did not contain a GLE. For both the GLE and noGLE conditions, we compute the empirical cumulative distribution function (ECDF) of F_{g1} and F_{r1}. Computation of ECDFs is discussed in Appendix A.1. Furthermore, we compute 100 bootstrap estimates of the ECDFs by drawing random samples from the F_{g1} and F_{r1} series, with same number of years as the GLE and noGLE conditions, but irrespective of whether or not GLEs occurred (more information on the bootstrap procedure, and an example, are provided in Appendix A.2). If the ECDFs of F_{g1} and F_{r1} for the GLE and noGLE conditions are markedly different from the ECDFs of the 100 bootstrap samples, then this is some evidence that years with and without GLEs affect the MB15 reconstruction differently. Taken in conjunction with a viable mechanism for GLEs affecting ^{10}Be, such a difference might suggest that GLEs do, in general, affect the MB15 reconstruction. These ECDF data are presented in Figure 8, with the F_{g1} and F_{r1} data analysed in the left and right columns, and the noGLE and GLE conditions given in the top and bottom rows. The observed ECDFs for the two conditions are drawn in red and blue for F_{g1} and F_{r1}, while the bootstrap ECDFs are drawn in grey. In each case, the ECDFs of the GLE and noGLE conditions are located within the spread of the bootstrap ECDFs, and this test provides no evidence that GLEs, in general, affect the MB15 reconstruction. This may be due to the limited sample of years. Additionally it may be that only the largest GLE can significantly perturb the ^{10}Be record. For example, MB15 estimate that an individual GLE needs exhibit a >100% increase in Neutron Monitor count rate(s) before it is likely to generate a quantity of ^{10}Be large enough to stand out from other sources of variability at annual resolution. Only approximately 15% of recorded GLE have been estimated to meet this criteria. However, we also note that at annual resolution it is possible that an accumulation of smaller GLE, and SEP events that do not register as GLE, could plausibly deposit the same net fluence of ^{10}Be producing energetic charged particles. Ideally we would further subdivide the data to assess the impact of GLE magnitude, but unfortunately there are too few samples for this to be reasonable.
Fig. 8 (A) The ECDF of F_{g1} is given in red, computed for only the 22 years without GLE events. The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 22 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 19 years with GLE events. 
4.2.3 HMF reconstruction differences with GGMS
We will now assess whether the distribution of F_{g1} and F_{r1} values depends on the occurrence of GGMS. We do not hypothesise that GGMS are directly affecting the production of cosmogenic isotopes, but do assume that the GGMS record serves as a proxy for the occurrence of large space weather events which may in turn affect the ^{10}Be record. There is an overlap of 115 years, from 1868–1983, between the GGMS record and the HMF reconstructions, of which 67 contained a GGMS and 48 did not contain a GGMS. Figure 9 presents the ECDF data with the F_{g1} and F_{r1} data analysed in the left and right columns, and the noGGMS and GGMS conditions given in the top and bottom rows. In this instance, for both F_{g1} and F_{r1}, there is a clear shift in the location of the ECDFs for the GGMS and noGGMS conditions relative to the bootstrap ECDFs. The ECDFS of F_{g1} and F_{r1} are shifted towards more negative values in years with GGMS, as expected by the increased ^{10}Be production idea.
A twosample twotailed KolmogorovSmirnov (KS) test was used to assess the hypothesis that the distributions of F_{g1} under the GGMS and noGGMS conditions were drawn from the same underlying distribution (Wilks, 1995) (see Appendix A.3). Essentially, this test computes the probability (pvalue) that the differences between the observed distributions would be this large due to drawing two random samples from the same underlying distribution. The test cannot confirm that the distributions are different, but does quantify how unusual it would be to obtain these samples if they came from the same underlying distribution.
In this instance, the KS test returns a pvalue of 0.004, meaning a 0.4% chance that differences this large would be observed under the noGGMS and GGMS samples if they were in fact drawn from the same distribution. Therefore, this test provides evidence that there are significant differences between the distribution of F_{g1} under the noGGMS and GGMS conditions. This analysis was repeated for F_{r1}, which returned a pvalue of 0.057, again showing that the differences between the F_{r1} distribution for the noGGMS and GGMS conditions are larger than we would expect due to random sampling. These and further statistics are summarised in Table 1. We cautiously interpret these results as evidence that space weather events do bias the MB15 reconstruction towards lower values of the HMF.
However, this analysis could be influenced by aliasing between the Schwabe/polarity cycle variations in the fractional differences F_{g1} and F_{r1} (Figs. 4 and 5) and the occurrence of GGMS (Fig. 7); GGMS are more likely to occur in years where the polarity phase composites show F_{g1} and F_{r1} are, on average, more negative. Therefore, we investigate this further by splitting the polarity cycle into active and quiet periods, where the occurrence frequency of GGMS varies much less within these periods. We define the active phases as φ_{p} ≤ π/2 and φ_{p} ≥ 3π/2 and quiet phases as π/2 < φ_{p} < 3π/2. There are 59 active years (46 with GGMS and 13 with noGGMS), and 56 quiet years (21 with GGMS and 35 with noGGMS). Figure 10 repeats the ECDF analysis for active years, and Figure 11 shows the results for quiet years. Subsampling by polarity cycle phase means that there are fewer samples in each category, and the results are less clear, particularly for the noGGMS condition in active years, which has only 13 samples. However, for both the active years and the quiet years the same trend is apparent; the GGMS distribution tends to be located at more negative values, to the left hand side of the bootstraps elements, than the noGGMS condition, which tends to be located to the right hand side of the bootstrap elements. We also repeated the KS test procedure on these active and quiet phase subsets, with the statistics included in Table 1. This reveals that for the active and quiet samples there is a much higher probability that differences between the GGMS and noGGMS conditions could be due to random sampling, and so we have low confidence that these distributions are actually different.
Fig. 9 (A) The ECDF of F_{g1} is given in red, computed for only the 48 years without GGMS events. The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 48 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 67 years with GGMS events. 
KolmogorovSmirnov test statistics comparing GGMS and NoGGMS conditions.
Fig. 10 (A) The ECDF of F_{g1} is given in red, computed for only the 13 years without GGMS events, from polarity phases corresponding to high solar activity (φ_{p} ≤ π/2 and φ_{p} ≥ 3π/2). The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 13 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 46 years with GGMS events. 
Fig. 11 (A) The ECDF of F_{g1} is given in red, computed for only the 35 years without GGMS events, from polarity phases corresponding to low solar activity (π/2 ≤ φ_{p} ≤ 3π/2). The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 35 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 21 years with GGMS events. 
5 Discussion and conclusions
We have compared the McCracken & Beer (2015) (MB15) annual ^{10}Be HMF reconstructions with two independent HMF reconstructions derived from sunspot records and geomagnetic activity. The four series span different periods of time, and the comparison is made only over the period of overlap, spanning 1845–1983. More restricted windows are also considered, in our analysis of GGMS (1868–1983), and also GLEs (1940–1983).
Comparing the distributions of yearly values of each HMF reconstruction revealed that the MB15 reconstructions display a modest negative bias relative to both the sunspot and geomagnetic reconstructions (Fig. 1). However, this bias is smaller than other sources of variability between the reconstructions, as demonstrated by the distributions of the fractional differences between the sunspot and geomagnetic reconstructions with the MB15 series (Fig. 3); these distributions were approximately centered on zero, and displayed a small asymmetry, such that there was more probability below the modal values than above. This comparison also demonstrated that the MB15 procedure to mitigate the impact of impulsive ^{10}Be enhancements on the HMF reconstruction does improve the agreement with the sunspot and geomagnetic HMF reconstructions; the fractional differences between B_{c2} with B_{g} and B_{r} are smaller than the differences computed with B_{c1}.
The time series of the fractional differences (Fig. 2) are suggestive of a solar cycle trend. Compositing the fractional differences as a function of solar cycle phase revealed a systematic variation in the mean fractional differences over the solar cycle. These variations are similar for each of the F_{g1}, F_{g2}, F_{r1}, and F_{r2} series. Given the independence of the sunspot and geomagnetic reconstructions, this implies that the solar cycle variation arises primarily from aspects of the MB15 reconstructions. We note that as the variations are similar for both B_{c1} and B_{c2}, they appear to be independent of the MB15 procedure to suppress impulsive ^{10}Be enhancements. The most likely explanation is that the solar cycle variation in the reconstruction differences is due to limitations in the ^{10}Be inversion procedure, specifically with the GCR transport model. MB15 describe how it was necessary to employ GCR transport coefficients which were an average over both solar polarity cycle, while it is known that GCR transport to nearEarth space varies significantly and systematically with polarity cycle phase.
We considered whether there was statistical evidence that GLEs and GGMS do, in general, affect the ^{10}Be record, and subsequent estimates of the HMF. The hypothesis for GLEs affecting the ^{10}Be record is simple, with the GLE acting as new temporary source of energetic charged particles, increasing the total population of particles capable of producing ^{10}Be in the upper atmosphere. MB15 presented evidence of an association between some impulsive ^{10}Be enhancements and known GLE events, as well as upper atmospheric nuclear bomb tests. We considered a larger set of GLE events, and investigated if the differences between the reconstructions were significantly different in years with and without GLE events. Using a twosample KS test, we could not resolve statistically significant differences between the distributions of reconstruction differences conditional on whether or not GLE occurred. This null result could arise from the limited sample size available to us, with only 41 years of overlap between the HMF record and GLE catalogue. Alternatively it could also indicate that GLEs do not, in general, significantly perturb the ^{10}Be record, which is consistent with modelling work performed by (Herbst et al., 2015). For example, as previously discussed, MB15 estimate that an individual GLE must exhibit a >100% increase in neutron monitor count rate before it can generate enough ^{10}Be to stand out above background variability and statistical uncertainty. Unfortunately we were unable to further subdivide the GLE data to test for the effect of GLE magnitude as there are too few events to reasonably do so. Given the limited extent of overlap between these records, and the generally low occurrence of GLE events, it is likely that that the question of if, and by how much, GLE affect the ^{10}Be record is best answered by physical modelling of ^{10}Be production and transport for different GLE scenarios. Such modelling work has begun, for example Usoskin & Kovaltsov (2012); Herbst et al. (2015) but will hopefully be extended to include more detailed GCR transport modelling, for more GLE scenarios and a wider range of solar activity states.
We also tested whether the differences between the reconstructions were significantly different in years with and without GGMS. In this instance, it was found that there was a statistically significant difference between the distributions of fractional differences in years with and without GGMS; the fractional residuals tended to be more positive in years with GGMS (consistent with an underestimate of the MB15 HMF by additional production of ^{10}Be by SEPs/GLEs). However, the correlation between the occurrence frequency of GGMS and polarity cycle phase means that this result may arise from aliasing between the observed solar cycle trend in the reconstruction differences and any effects due to space weather events. To try and isolate this possibility, we further divided the data into “active” and “quiet” phases, and repeated the procedure to assess if the reconstruction differences were significantly different in years with and without GGMS. The same trend persists, with the distribution of reconstruction differences typically being more positive in years with GGMS, but this is less clear, being heavily impacted by the sample size reduction from subsampling. Applying the same KS testing procedure, we could not resolve statistically significant differences between the reconstruction differences in these categories. Therefore we conclude there is some evidence that the MB15 reconstructions do vary in response to GGMS, but that this variation is small compared to other sources of variability. There are several plausible hypothesis about why we may expect the MB15 HMF reconstruction to depend on the occurrence of GGMS. Most simply, it is reasonable to assume that the occurrence of GGMS is a weak proxy for the occurrence of SEP and GLE events, as the energetic CMEs that drive GGMS may also possibly generate significant SEP fluxes and GLEs. Correspondingly, Figure 7 shows a weak correlation between the annual count of GLE and GGMS. However, there are also other possible mechanisms, which may have a nonnegligible impact. For example, magnetospheric dynamics during GGMS modifies the access of energetic charged particles into the magnetosphere, and consequently the upper atmosphere. Finally, we do not rule out that changes in the global structure of the heliosphere over the solar cycle may have some role in explaining these differences.
Recent research has shown that analysis of cosmogenic isotopes may be the only available proxy for historic extreme SEP events (Usoskin et al., 2006; Miyake et al., 2012, 2013; Usoskin et al., 2013; Herbst et al., 2015; McCracken & Beer, 2015; Dee et al., 2016). This would be a hugely valuable data source in constraining the limits of extreme SEP events, which have quite poor statistics over the space age (Barnard & Lockwood, 2011), and are critical for assessing the hazards and risks to our continuing exploration and development of the solar system (Cannon et al., 2013), especially for manned missions. Furthermore, the accuracy of reconstructions of the historic GCR environment and large scale solar magnetic activity will improve if the impacts of SEP fluxes on these records are resolved and removed.
As future technological advances facilitate cosmogenic isotope records with higher temporal resolution and higher precision it should be expected that limitations in the cosmogenic isotope inversion procedures become increasingly apparent. We consider the results of this study strong motivation for investigating the role of more advanced GCR transport models in the inversion of high resolution cosmogenic isotope records, as well as developing robust procedures to identify the possible impacts of space weather events on cosmogenic isotope records.
Acknowledgements
L.B, M.J.O and M.L thank the Science and Technologies Facilities Council for support under grant ST/M000885/1. K.G.McC thanks the International Space Science Institute (ISSI) for ongoing support. The editor thanks Olga Malandraki and Konstantin Herbst for their assistance in evaluating this paper. The authors also thanks these two reviewers for providing helpful evaluations of this article. The analysis code supporting this article is available at https://github.com/LukeBarnard/ReconCompare.
Appendix A: Statistical methods
Below we provide further details regarding some statistical methods employed throughout the analysis in this article. Further details on all of these techniques are described in detail in Wilks (1995).
A.1 Empirical cumulative distribution functions
The ECDF provides a nonparametric estimate of a variables cumulative distribution. For a variable, x, the ECDF, F(x), is calculated as (4) where i is the order statistic of x, and n is the total number of samples.
A.2 Bootstrap resampling
In Section 4.2 we use bootstrap resampling to help assess whether the ECDFs of the reconstruction differences show a dependence on whether or not GLEs of GGMs occurred. Here we illustrate this process by way of an artificial example. The left panel of Figure A.1 shows a scatter plot of the daily mean aa index against the daily sunspot number, R, from 1868 to 2015. All these data are shown with black dots, except those values with R > 400. In this example we are interested in whether or not the distribution of aa is different for R > 400, than for other times. The right panel of Figure A.1 shows the ECDF for all daily mean aa values in black, and the ECDF for aa values corresponding to R > 400 in red, for which there are 66 values. Comparing only the black and red lines, it is difficult to assess whether there is evidence that the R > 400 condition significantly affects the distribution of aa. However, the grey lines are 100 bootstrap estimates of the aa ECDF, calculated by randomly sampling 66 values without replacement from the full aa distribution. Unsurprisingly these grey lines of the bootstrap estimates all cluster around the full aa distribution. However, the red line lies predominantly outside of the cluster of bootstrap lines, suggesting that the R > 400 condition does yield a sample which is unusual compared to a random sample of the same size from the full distribution.
Fig. A.1 The left panel shows a scatter plot of the the daily mean aa index agaisnt the daily sunspot number, R. The black dots show all values, while those in red mark values with R > 400. The right hand panel shows the ECDF of all aa values in black, of the aa values corresponding to R > 400 in red, and for 100 bootstrap estimates of the aa distribution in grey. 
A.3 Twosample KolmogorovSmirnov tests
The twosample KS test is a statistical hypothesis test which can be used to compare two observed distributions, and assess the plausibility of the nullhypothesis that both observed distributions are in fact random samples from the same underlying distribution. As previously mentioned, the KS test cannot confirm that the distributions are different, but does quantify how unusual it would be to obtain two sets of samples if they came from the same underlying distribution. This is a nonparametric test, with the test statistic D being given by (5) where F_{1}(x_{1}) and F_{2}(x_{2}) are the ECDF's of the two observed distributions, such that D is the maximum difference between them. From the D statistic, a pvalue can be calculated which gives the probability of obtaining this maximum difference between the two distributions, had they been randomly sampled from the same underlying distribution. Consequently, small pvalues imply that the null hypothesis is not very plausible, and the differences between F_{1}(x_{1}) and F_{2}(x_{2}) can probably not be explained by random sampling. Of course, this does not mean any alternative hypotheses are necessarily true, just that there is a difference that likely merits further investigation.
References
 Abreu JA, Beer J, Steinhilber F, Tobias SM, Weiss NO. 2008. For how long will the current grand maximum of solar activity persist ? Geophys Res Lett 35: 2–5. DOI:10.1029/2008GL035442, http://www.agu.org/pubs/crossref/2008/2008GL035442.shtml. [NASA ADS] [CrossRef] [Google Scholar]
 Barnard L, Lockwood M. 2011. A survey of gradual solar energetic particle events. J Geophys Res 116: 1–13. DOI:10.1029/2010JA016133, http://www.agu.org/pubs/crossref/2011/2010JA016133.shtml. [CrossRef] [Google Scholar]
 Barnard L, Lockwood M, Hapgood MA, Owens MJ, Davis CJ, Steinhilber F. 2011. Predicting space climate change. Geophys Res Lett 38: 7–12. DOI:10.1029/2011GL048489. http://www.agu.org/pubs/crossref/2011/2011GL048489.shtml. [Google Scholar]
 Barnard L, Owens MJ, Scott CJ. 2017. The space environment before the space age. Astron Geophys 58: 2.12–2.16. DOI:10.1093/astrogeo/atx056. [CrossRef] [Google Scholar]
 Cannon P, Angling M, Barclay L, Curry C, Dyer C, et al. 2013. Extreme space weather : impacts on engineered systems and infrastructure. Tech Rep R Acad Eng www.raeng.org.uk/spaceweather, ISBN:1903496950. [Google Scholar]
 Clette F, Lefèvre L. 2016. The new sunspot number: assembling all corrections. Sol Phys 291: 2629–2651. DOI:10.1007/s112070161014y, http://arxiv.org/abs/1510.06928. [CrossRef] [Google Scholar]
 Cliver EW. 2016. comparison of new and old sunspot number time series. Sol Phys (1994) 291, DOI:10.1007/s1120701609297. [Google Scholar]
 Dee M, Pope B, Miles D, Manning S, Miyake F. 2016. Supernovae and singleyear anomalies in the atmospheric radiocarbon record. Radiocarbon 1006: 1–10. DOI:10.1017/RDC.2016.50, http://www.journals.cambridge.org/abstract_S0033822216000503. [Google Scholar]
 Dreschhoff GAM, Zeller EJ. 1990. Evidence of individual solar proton events in Antarctic snow. Sol Phys 127: 333–346. DOI:10.1007/BF00152172. [CrossRef] [Google Scholar]
 Duderstadt KA, Dibb JE, Jackman CH, Randall CE, Solomon SC, Mills MJ, Schwadron NA, Spence HE. 2014. Nitrate deposition to surface snow at Summit, Greenland, following the 9 November 2000 solar proton event. J Geophys Res Atmos 119: 6938–6957. DOI:10.1002/2013JD021389. [CrossRef] [Google Scholar]
 Forbush SE. 1946. Three unusual cosmicray increases possibly due to charged particles from the Sun. Phys Rev 70: 771–772. DOI:10.1103/PhysRev.70.771. [NASA ADS] [CrossRef] [Google Scholar]
 Heikkilä U, Beer J, Feichter J. 2007. Modeling cosmogenic radionuclides ^{10}Be and ^{7}Be during the Maunder Minimum using the ECHAM5HAM General Circulation Model. Atmos Chem Phys Discuss 7: 15341–15372. DOI:10.5194/acpd7153412007, http://www.atmoschemphysdiscuss.net/7/15341/2007/. [CrossRef] [Google Scholar]
 Holappa L, Mursula K, Asikainen T. 2014a. A new method to estimate annual solar wind parameters and contributions of different solar wind structures to geomagnetic activity. J Geophys Res Space Phys 119: 9407–9418. DOI:10.1002/2014JA020599. [NASA ADS] [CrossRef] [Google Scholar]
 Holappa L, Mursula K, Asikainen T, Richardson IG. 2014b. Annual fractions of highspeed streams from principal component analysis of local geomagnetic activity. J Geophys Res Sp Phys 119: 4544–4555. DOI:10.1002/2014JA019958. [CrossRef] [Google Scholar]
 Herbst K, Heber B, Beer J, Tylka A. 2015. Modelling the production of cosmogenic radionuclides due to galactic and solar cosmic rays. In: Proceedings of the 34th Int. Cosm. Ray Conf. of Science. vol. 236, https://pos.sissa.it/236/537/pdf. [Google Scholar]
 Kepko L, Spence HE, Smart DF, Shea MA. 2009. Interhemispheric observations of impulsive nitrate enhancements associated with the four large groundlevel solar cosmic ray events (1940–1950). J Atmos SolarTerr Phys 71: 1840–1845. DOI:10.1016/j.jastp.2009.07.002. http://adsabs.harvard.edu/abs/2009JASTP.71.1840K. [CrossRef] [Google Scholar]
 Kilpua EKJ, Olspert N, Grigorievskiy A, Käpylä MJ, Tanskanen EI, Miyahara H, Kataoka R, Pelt J, Liu YD. 2015. Statistical study of strong and extreme geomagnetic disturbances and solar cycle characteristics. ApJ 806: 272. DOI:10.1088/0004637X/806/2/272, http://stacks.iop.org/0004637X/806/i=2/a=272?key=crossref.874e87f7ac81b1a48c1541d391970e36. [CrossRef] [Google Scholar]
 Lefèvre L, Vennerstrøm S, Dumbović M, Vršnak B, Sudar D, Arlt R, Clette F, Crosby N. 2016. Detailed analysis of solar data related to historical extreme geomagnetic storms: 18682010. Sol Phys 291, DOI:10.1007/s1120701608923. [Google Scholar]
 Lockwood M, Stamper R, Wild MN. 1999. A doubling of the Sun's coronal magnetic field during the past 100 years. Nature 399: 437–439. DOI:10.1038/20867. [NASA ADS] [CrossRef] [Google Scholar]
 Lockwood M, Owens MJ, Barnard L, Davis CJ, Thomas S. 2012. Solar cycle 24: what is the Sun up to ? Astron Geophys 53: 3.09–3.15. DOI:10.1111/j.14684004.2012.53309.x. [Google Scholar]
 Lockwood M. 2013. Reconstruction and prediction of variations in the open solar magnetic flux and interplanetary conditions. Living Rev Sol Phys 10: DOI:10.12942/lrsp20134. [CrossRef] [Google Scholar]
 Lockwood, M, Barnard L, Nevanlinna H, Owens MJ, Harrison RG, Rouillard AP, Davis CJ. 2013a. Reconstruction of geomagnetic activity and nearEarth interplanetary conditions over the past 167 year − part 1: a new geomagnetic data composite. Ann Geophys 31: 1957–1977. DOI:10.5194/angeo3119572013, http://www.anngeophys.net/31/1957/2013/. [CrossRef] [Google Scholar]
 Lockwood M, Barnard L, Nevanlinna H, Owens MJ, Harrison RG, Rouillard AP, Davis CJ. 2013b. Reconstruction of geomagnetic activity and nearEarth interplanetary conditions over the past 167 year − part 2: a new reconstruction of the interplanetary magnetic field. Ann Geophys 31: 1979–1992. DOI:10.5194/angeo3119792013, http://www.anngeophys.net/31/1979/2013/. [NASA ADS] [CrossRef] [Google Scholar]
 Lockwood M, Nevanlinna H, Barnard L, Owens MJ, Harrison RG, Rouillard AP, Scott CJ. 2014a. Reconstruction of geomagnetic activity and nearEarth interplanetary conditions over the past 167 year − part 4: nearEarth solar wind speed, IMF, and open solar flux. Ann Geophys 32: 383–399. DOI:10.5194/angeo323832014, http://www.anngeophys.net/32/383/2014/. [NASA ADS] [CrossRef] [Google Scholar]
 Lockwood M, Nevanlinna H, Vokhmyanin M, Ponyavin DI, Sokolov S, Barnard L, Owens MJ, Harrison RG, Rouillard AP, Scott CJ. 2014b. Reconstruction of geomagnetic activity and nearEarth interplanetary conditions over the past 167 year − Part 3: Improved representation of solar cycle 11. Ann Geophys 32: 367–381. DOI:10.5194/angeo323672014, http://www.anngeophys.net/32/367/2014/. [CrossRef] [Google Scholar]
 Lockwood M, Owens MJ, Barnard L, Usoskin IG. 2016a. An assessment of sunspot number data composites over 1845–2014. Astrophys J 824: 54. DOI:10.3847/0004637X/824/1/54, http://stacks.iop.org/0004637X/824/i=1/a=54?key=crossref.5961a0acc8d7b507dd06a81a140f33ff. [NASA ADS] [CrossRef] [Google Scholar]
 Lockwood M, Owens MJ, Barnard LA, Bentley S, Scott CJ, Watt CE. 2016b. On the origins and timescales of geoeffective IMF. Space Weather 14: 406–432. DOI:10.1002/2016SW001375. [CrossRef] [Google Scholar]
 Matthes K, Funke B, Anderson ME, Barnard L, Beer J, et al. 2016. Solar forcing for CMIP6 (v3.1). Geosci Model Dev Discuss 6: 1–82. DOI:10.5194/gmd201691, http://www.geoscimodeldevdiscuss.net/gmd201691/. [CrossRef] [Google Scholar]
 McCracken KG. 2004. Geomagnetic and atmospheric effects upon the cosmogenic ^{10}Be observed in polar ice. J Geophys Res 109: A04101, DOI:10.1029/2003JA010060. [Google Scholar]
 McCracken KG, Beer J. 2015. The annual cosmicradiation intensities 1391–2014; the annual heliospheric magnetic field strengths 1391–1983, and identification of solar cosmicray events in the cosmogenic record 1800–1983. Sol Phys 290: 3051–3069. DOI:10.1007/s112070150777x. [NASA ADS] [CrossRef] [Google Scholar]
 McCracken KG, Dreschhoff GAM, Zeller EJ, Smart DF, Shea MA. 2001. Solar cosmic ray events for the period 1561–1994 1. Identification in polar ice, 1561–1950. J Geophys Res 106: 21,585–21,598. DOI:10.1029/2000JA000237, http://www.agu.org/pubs/crossref/2001/2000JA000237.shtml. [NASA ADS] [CrossRef] [Google Scholar]
 Miyake F, Nagaya K, Masuda K, Nakamura T. 2012. A signature of cosmicray increase in ad 774775 from tree rings in Japan. Nature 486: 774–776. DOI:10.1038/nature11123. [CrossRef] [Google Scholar]
 Miyake F, Masuda K, Nakamura T. 2013. Another rapid event in the carbon14 content of tree rings. Nat Commun 4: 1748. DOI:10.1038/ncomms2783, http://www.ncbi.nlm.nih.gov/pubmed/23612289. [NASA ADS] [CrossRef] [Google Scholar]
 Mursula K, Holappa L, Lukianova R. 2016. Seasonal solar wind speeds for the last 100 years: unique coronal hole structures during the peak and demise of the Grand Modern Maximum. Geophys Res Lett 44, DOI:10.1002/2016GL071573. [Google Scholar]
 Owens MJ, Lockwood M. 2012. Cyclic loss of open solar flux since 1868: the link to heliospheric current sheet tilt and implications for the Maunder Minimum. J Geophys Res Space Phys 117: 1–9. DOI:10.1029/2011JA017193. [CrossRef] [Google Scholar]
 Owens MJ, McCracken KG, Lockwood M, Barnard L. 2015. The heliospheric Hale cycle over the last 300 years and its implications for a “lost” late 18th century solar cycle. J Space Weather Space Clim 5: A30. DOI:10.1051/swsc/2015032. [CrossRef] [Google Scholar]
 Owens MJ, Cliver EW, McCracken KG, Beer J, Barnard L, et al. 2016a. NearEarth heliospheric magnetic field intensity since 1750. part 1: sunspot and geomagnetic reconstructions. J Geophys Res Space Phys 121, DOI:10.1002/2016JA022529. [Google Scholar]
 Owens MJ, Cliver EW, McCracken KG, Beer J, Barnard L, et al. 2016b. NearEarth heliospheric magnetic field intensity since 1750. part 2: cosmogenic radionuclide reconstructions. J Geophys Res Space Phys 121, DOI:10.1002/2016JA022550. [Google Scholar]
 Shea MA, Smart DF, Dreschhoff GAM. 1999. Identification of major proton fluence events from nitrates in polar ice cores. Radiat Meas 30: 309–316. DOI:10.1016/S13504487(99)000578, http://linkinghub.elsevier.com/retrieve/pii/S1350448799000578. [CrossRef] [Google Scholar]
 Steinhilber F, Abreu JA, Beer J, Brunner I, Christl M, et al. 2012. 9400 years of cosmic radiation and solar activity from ice cores and tree rings. Proc Natl Acad Sci U.S.A. 109: 5967–5971. DOI:10.1073/pnas.1118965109, http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=3341045&tool=pmcentrez&rendertype=abstract. [NASA ADS] [CrossRef] [Google Scholar]
 Strauss RDT, Effenberger F. 2017. A hitchhiker's guide to stochastic differential equations. Space Sci Rev 212, DOI:10.1007/s112140170351y. [Google Scholar]
 Svalgaard L. 2014. Errors in scale values for magnetic elements for Helsinki. (March). Ann Geophys 32, 633–641, DOI:10.5194/angeo326332014, http://arxiv.org/abs/1403.2707. [NASA ADS] [CrossRef] [Google Scholar]
 Svalgaard L, Cliver EW. 2005. The IDV index: its derivation and use in inferring longterm variations of the interplanetary magnetic field strength. J Geophys Res 110: 103. DOI:10.1029/2005JA011203. [NASA ADS] [CrossRef] [Google Scholar]
 Svalgaard L, Cliver EW. 2010. Heliospheric magnetic field 1835–2009. J Geophys Res 115: 111. DOI:10.1029/2009JA015069, http://doi.wiley.com/10.1029/2009JA015069. [NASA ADS] [CrossRef] [Google Scholar]
 Thomas SR, Owens MJ, Lockwood M, Barnard L, Scott CJ. 2015. NearEarth cosmic ray decreases associated with remote coronal mass ejections. ApJ 801: 5. DOI:10.1088/0004637X/801/1/5, http://stacks.iop.org/0004637X/801/i=1/a=5?key=crossref.2044e2d221606382f66ca6a2df3966e5. [CrossRef] [Google Scholar]
 Usoskin IG. 2013. A history of solar activity over millennia. Living Rev Sol Phys 10: DOI:10.12942/lrsp20131, http://www.livingreviews.org/lrsp20131. [NASA ADS] [CrossRef] [Google Scholar]
 Usoskin IG, Kovaltsov GA. 2012. Occurrence of extreme solar particle events: assessment from historical proxy data. ApJ 757: 92. DOI:10.1088/0004637X/757/1/92, http://stacks.iop.org/0004637X/757/i=1/a=92?key=crossref.8b9f2d4f1c17e9be6434333d4857e488. [NASA ADS] [CrossRef] [Google Scholar]
 Usoskin IG, Solanki SK, Kovaltsov GA, Beer J, Kromer B. 2006. Solar proton events in cosmogenic isotope data. Geophys Res Lett 33: 2–5. DOI:10.1029/2006GL026059, http://www.agu.org/pubs/crossref/2006/2006GL026059.shtml. [Google Scholar]
 Usoskin IG, Kromer B, Ludlow F, Beer J, Friedrich M, Kovaltsov GA, Solanki SK, Wacker L. 2013. The AD775 cosmic event revisited: the Sun is to blame. A&A 552: L3. DOI:10.1051/00046361/201321080, http://www.aanda.org/10.1051/00046361/201321080. [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Vennerstrom S, Lefèvre L, Dumbović M, Crosby N, Malandraki OE, et al. 2016. Extreme geomagnetic storms − 18682010. Sol Phys 291, DOI:10.1007/s112070160897y. [CrossRef] [Google Scholar]
 Vieira LEA, Solanki SK. 2010. Evolution of the solar magnetic flux on time scales of years to millenia. A&A 509: A100. DOI:10.1051/00046361/200913276. [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Wilks DS. 1995. Statistical methods in the atmospheric sciences. Academic Press, ISBN 978012751965–4. [Google Scholar]
 Wolff EW, Bigler M, Curran MAJ, Dibb JE, Frey MM, Legrand M, McConnell JR. 2012. The Carrington event not observed in most ice core nitrate records. Geophys Res Lett 39: 1–5. DOI:10.1029/2012GL051603, http://www.agu.org/pubs/crossref/2012/2012GL051603.shtml. [NASA ADS] [CrossRef] [Google Scholar]
Cite this article as: Barnard L, McCracken KG, Owens MJ, Lockwood M. 2018. What can the annual ^{10}Be solar activity reconstructions tell us about historic space weather? J. Space Weather Space Clim. 8: A23
All Tables
All Figures
Fig. 1 (A–D) Time series of the annual HMF reconstructions from 1845–1983, derived from; ^{10}Be concentrations, B_{c1} (with impulsive enhancements), and B_{c2} (with suppressed impulsive enhancements); geomagnetic observations, B_{g}; and sunspot observations, B_{r}. Error bars are 1σ errors of the mean. The red squares in panel (A) mark years where B_{c1} ≠ B_{c2}. (E) Kernel density estimates of the distributions of four different HMF reconstructions. 

In the text 
Fig 2 Time series of the fractional differences between the B_{c1} and B_{c2} reconstructions with both B_{g} and B_{r}. Computation of the fractional differences is explained in Section 3.1: F_{g1} is the fractional difference between B_{c1} and B_{g}, F_{g2} between B_{c2} and B_{g}, F_{r1} between B_{c1} and B_{r}, and F_{r2} between B_{c2} and B_{r}. Error bars are 1σ errors of the mean. 

In the text 
Fig. 3 (A) Kernel density estimates of the distributions of fractional differences between B_{g} with B_{c1} and B_{c2}, F_{g1} and F_{g2}, respectively. The purple and orange dots show the fractional differences in years (McCracken & Beer, 2015) identified as being affected by SEP events (the vertical coordinates of these points is randomised, which otherwise overlap). (B) The same as panel (A), except showing the fractional differences of B_{r} with B_{c1} and B_{c2}. 

In the text 
Fig. 4 (A) Composite of the HMF reconstructions with Schwabe solar cycle phase. (B) Composite of the fractional differences of B_{g} and B_{r} with B_{C1}. (C) Composite of the fractional differences of B_{g} and B_{r} with B_{C2}. For each panel, the error bars are 1σ errors of the mean. In panels B and C, the black dashed line shows the zero fraction difference line. 

In the text 
Fig. 5 Panels A–F show composites of the HMF reconstructions and the fractional differences F_{g1}, F_{r1}, F_{g2}, F_{r2}, with polarity cycle phase. The composites are split according to solar magnetic polarity, with qA > 0 cycles in panels A–C, and qA < 0 cycles in panels D–F. For each panel, the error bars are 1σ errors of the mean. 

In the text 
Fig. 6 (A) Time series of the aa geomagnetic index (black), from 1867–2014, with the maxima of great geomagnetic storms marked by red circles. Also shown are the timings of known GLE events (cyan lines). Panel (B) Shows the same data as (A), except as a function of Schwabe cycle phase. Panel (C) Shows the same data as (A), except as a function of polarity cycle phase. We note that, due to the availability of ^{10}Be data, we only use the geomagnetic and GLE observations up until 1983. 

In the text 
Fig. 7 (A) Time series of the annual count of great geomagnetic storms (red), and GLE events (cyan). (B) The number of great geomagnetic storms (red) and GLE events (cyan) as a function of Schwabe cycle phase. (C) The number of geomagnetic storms (red) and GLE events (cyan) as a function of polarity cycle phase. In panels (B) and (C), phase bins 2π/10 wide were used. No account was taken of the uneven sampling of solar cycle phase, or the fact the GLE record only spans 1937present; these plots only show coarse trends with solar/polarity cycle phase, and the absolute number of geomagnetic storms and GLE events should not be directly compared. We note that, due to the availability of ^{10}Be data, we only use the geomagnetic and GLE observations up until 1983. 

In the text 
Fig. 8 (A) The ECDF of F_{g1} is given in red, computed for only the 22 years without GLE events. The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 22 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 19 years with GLE events. 

In the text 
Fig. 9 (A) The ECDF of F_{g1} is given in red, computed for only the 48 years without GGMS events. The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 48 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 67 years with GGMS events. 

In the text 
Fig. 10 (A) The ECDF of F_{g1} is given in red, computed for only the 13 years without GGMS events, from polarity phases corresponding to high solar activity (φ_{p} ≤ π/2 and φ_{p} ≥ 3π/2). The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 13 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 46 years with GGMS events. 

In the text 
Fig. 11 (A) The ECDF of F_{g1} is given in red, computed for only the 35 years without GGMS events, from polarity phases corresponding to low solar activity (π/2 ≤ φ_{p} ≤ 3π/2). The grey lines show 100 bootstrap estimates of the F_{g1} ECDF, computed by randomly sampling 35 years from the F_{g1} series. Panel B has the same structure as panel A, but instead shows the ECDF of F_{r1} in blue. Panels (C) and (D) have the same structure as (A) and (B), but instead show the ECDFs of F_{g1} and F_{r1} from only those 21 years with GGMS events. 

In the text 
Fig. A.1 The left panel shows a scatter plot of the the daily mean aa index agaisnt the daily sunspot number, R. The black dots show all values, while those in red mark values with R > 400. The right hand panel shows the ECDF of all aa values in black, of the aa values corresponding to R > 400 in red, and for 100 bootstrap estimates of the aa distribution in grey. 

In the text 