Volcanic Impacts Dominate Bidecadal‐Multidecadal Temperature Variations During the Late Holocene in Northern Fennoscandia

A reconstruction of summer temperatures during two millennia in Northern Fennoscandia has been studied and 79 abrupt cold periods have been identified. Data regarding the volcanic eruptions are not a priori used. During 1750–2000 all of the abrupt, cold periods could posteriori be connected to candidate years of volcanic eruptions. It is concluded that the 79 abrupt cold periods are mainly, if not all, postvolcanic ones. Local impact of weaker, less global (Icelandic) or some unknown eruptions is evident in the results. In addition to the short cold periods, we found 42 longer bidecadal‐multidecadal oscillations during 50–1980 AD. The short abrupt cold periods and the longer oscillations were characterized as the atmospheric and oceanic responses to volcanic eruptions, respectively. The volcanic activity has clearly been a dominant contributor to the temperature variations in Northern Fennoscandia during the two last millennia.


Introduction
Aerosol dust loadings into the atmosphere caused by volcanic eruptions are often followed by abrupt cold peaks in sthe temperature reconstructions derived from tree rings (e.g. Briffa et al., 1998;D'Arrigo et al., 2006;Timmreck et al., 2009). To study the relation between volcanic eruptions and temperature changes, an estimate of the temperature without the volcanic impact is needed. The cold deviation is then estimated as the difference between this hypothetic, undisturbed state and the reconstructed temperature.
To get an estimate of the undisturbed temperature, properties of the temperature reconstruction are utilized. For instance, if the temperatures are given as anomalies, then the cold peaks can be estimated as the deviations from the zero line. Another alternative could be to use smoothed values of the temperature reconstruction. Generally, we will here define "baseline" to mean sets of temperature estimates used in computing the cold deviations.
The selection of a baseline determines the accuracy of the computed cold peaks and identification of the related volcanic events. The temperature reconstruction itself is prone to estimation errors, further decreasing the accuracy. In some cases, the short-range variation ("weather") can be so strong that the identification of the cold peak is prevented.
To avoid random variation, several eruption events or temperature reconstructions with variable spacing may be applied (e.g. Esper et al., 2015;Schneider et al., 2016;Sigl et al., 2015;Stoffel et al., 2015;Toohey et al., 2019). The significance of the results can be increased further by studying only very strong volcanic eruptions.
In local studies possibilities to damp random variations are limited because often only one temperature reconstruction is available. In the following, one extreme case is discussed as an example. Ice cap records from Arctic Canada suggest that intense ice growth took place between 1275 and 1300 AD, possibly linked to four preceding large volcanic events (Miller et al., 2012). Gennaretti et al. (2014) then suggested that the Samalas eruption 1257 could have been a precursor of the cold period. This is reasonable as the eruption was one of the strongest during the latest two millennia and took place just before the period studied in Miller et al. (2012). The abrupt cold change associated with the Samalas eruption is well seen in the local temperature reconstruction used. On the other hand, in their global study Timmreck et al. (2009) found that "the surface cooling was not substantially larger than for Pinatubo (≈0.4 K°)." Later, the analysis in Guillet et al. (2017) indicated warming and cooling over eastern Canada during 1258 and 1259, respectively. Indeed, Gennaretti et al. (2014) had to apply a sophisticated approach in order to show that the cooling in 1257 was significantly strong in their local temperature reconstruction.
Here the aim is to introduce a method for detecting abrupt cold periods based on one local reconstruction. Further, a simple method will be presented to test the significance of the cold peaks. Understanding the risks of random errors when applying locally one single temperature reconstruction, we concentrate only on the few cases of strong leading peaks; that is, we put emphasis on finding the first year of the significant temperature drop of an abrupt cold period. The method will be applied for a temperature reconstruction describing conditions in Northern Fennoscandia.

Data
The temperature reconstruction applied is adopted from Esper et al. (2012) for Northern Fennoscandia. The eruption data used are given in Sigl et al. (2015) and in the Online resource of Crowley and Unterman (2013). Sea-ice data are obtained from Kinnard et al. (2011), Supporting information S4. In addition to these, some Icelandic eruptions are taken from "A short history of earthquake activity on South Icelandic Seismic Zone (SISZ)" (Jón Frimann, 2011; http://www.jonfr.com/volcano/?p=745).

Method
In the following, we try to find significant cold peaks by using only one temperature reconstruction. More precisely, we concentrate on determining the first cold deviation (FCD) of abrupt cold periods. For that purpose, we estimate the baseline values and their random error. The intensity of the encountered cold drop is determined as the difference between the baseline and the temperature reconstruction. The significance of the peak is determined with the aid of the random variance estimated. In addition to that, the computations are extended to estimate the length of the cold period, simultaneously modifying the baseline.
To form the baseline needed to detect cold breaks, (summer) temperatures for year t 0 are linearly predicted from preceding values with a prediction error filter as follows: Here T(t) is the temperature series to be analyzed. T P (t) is the linearly predictable part of the time series T(t) and gives our estimate of the baseline. The filter length is m years.
The coefficients a i and the prediction error variance p m 2 (residual variance) are here determined by the Maximum Entropy Method often applied in spectrum analysis (Burg, 1967). For a Gaussian process, the resulting formula is the same as when optimally fitting an autoregressive model to a time series. The temperature reconstruction studied (T esp ) is given in Esper et al. (2012) and the deviations from the baseline are determined as T P (t)-T esp (t). The prediction method is sensitive to the selection of the length of the prediction error filter.
Imagine that the temperature reconstruction is smooth so that there are essentially no short-scale variations. Then any proper baseline should indicate no short-scale deviations. We selected the filter length in such a way that our baseline fulfills this condition. We smoothed T esp (t) (100-year cubic splines) to give T′ esp (t) and applied equation (1) by taking T(t) = T′ esp (t). With m = 15 the prediction T P (t) was nearly equal to T′ esp (t); that is, the method could identify the long-term variation as a predictable component and indicated no essential deviations from it. Accordingly, we will apply filter lengths of m = 15 years throughout this study.
The computations are made for temperatures during 0-2000 AD. For the leading 15 years, we have to take T P (t) = T esp (t); otherwise, T P (t) is determined from equation (1) by taking T(t) = T esp (t). The performance of the approach is illustrated in Figure 1a for the significant FCD found for 1068 AD.
A significant (97.5%, one-sided) cold year is termed as such, if Apart from 1068, no other year in Figure 1a fulfills this criterion. The cold period is here defined as to consist of successive years with T esp (t) ≤ T P (t). Counting from the FCD, the cold period 1068-1070 in Figure 1a lasts three years. It could have been longer if only T P (1071) had been slightly higher. In any case, the length seems to be an underestimation.
After 1068 T P (t) decreases rapidly in Figure 1a because it is, by nature of the linear prediction, contaminated by the prevailing cold temperatures, especially by the cold temperature of 1068. Such a decrease is far from the concept of an ideal baseline. Therefore, we repeat the computations modifying T(t) in equation (1) by replacing the value of T(1068) = −0.9°(thin red line in Figure 1a) with T P (1068) = 0.7°(thick blue line). Note that in this transformation the prediction for 1068 does not change whereas the following temperatures will be less contaminated by the coldness. We do similar changes for all significant cold years found during 2000 years within the cold periods. We then repeat the computations and get new values for the baseline and the prediction error variance. Note that these computations are only auxiliary to the main aim of identifying the first cold drop where no essential changes should be observable.
Computations are repeated several times until the seventh iteration step, each time replacing the significant cold years with the previous T P (t) and determining new coefficients a i plus the prediction error variance p m 2 .
After the seventh iteration step the results do not essentially change, this holds especially for the error variance and the testing results in equation (2). Note that in these computations, no eruption data are a priori used. Figure 1b illustrates the final result. The baseline is not changed much for the particular year 1068 AD ( Figure 1b) but is otherwise now much warmer and the length of the cold period is eight years, or close to 10 years.
The cold periods lasting only for one year will in the following be handled as if they were of a random origin, that is, then the series to be predicted is always T esp (t) and no replacement is made in the iterative computations.
In Figure 1 the smoothed reconstruction T′ esp (t) is compared with our baseline. Applying that smoothed version as the baseline, we would get 1066-1078 for the estimate of the period length, but on the other hand the intensity of the FCD would be reduced. The only aim of our baseline is to determine the intensity of the FCDs as well as possible. Thus, the success of our approach is to be measured by the number of the significant FCDs found.

Results
Hereafter, all of the results are based on equations (1) and (2) using m = 15 and applying seven iteration cycles.

Abrupt Cold Periods
We have identified significant cold peaks followed by at least one cold year. Altogether, 79 such cases were found (Figures 2-4). Eruption intensities and timings (Sigl et al., 2015, before 800 AD excluding the Southern Hemisphere cases; Crowley & Unterman, 2013, after 800 AD) within or one to two years before the cold period are only posteriori added into the illustrations; otherwise, the eruption data have not been applied. It is seen that often an eruption has taken place at the FCD; that is, the cold periods turn out to be caused by the volcanism. Exceptions are the well-known intense eruptions of 1257 (Samalas; Figure 3b), 1809 (unknown; Figure 4), and 1815 (Tambora; Figure 4) without any significant FCD in the temperature proxy used. However, weaker delayed peaks may be related to those eruptions, but they do not fulfill our criteria.   In 35 cases out of the 79 events, the length of the period is only two to three years. Thus, normally the periods are rather short, in agreement with Esper et al. (2013) (one to two years), Sigl et al. (2013) (one to two significant years), and the lifetime of stratospheric volcanic sulfur in Toohey et al. (2019) (one to two years in the north). Note that our period lengths are not directly comparable with the conventional ones where the lengths are often given by the number of significant cold postvolcanic years, whereas we include all successive years with negative deviations. On the other hand, the length is here counted only from the FCD forward ( Figure 1b). Further, any cold period may in our terms be coincident with several successive volcanic events. Esper et al. (2012) list the 30 coldest years using the very same data as we do. Their aim was to apply those peaks to compare different reconstructions. The case resembles ours as no eruption data were used. On the other hand, in their case no significance test was needed and they accepted any cold peaks whereas we count only significant FCDs. They found about 15 of our FCD years. Thus, the results and also the aims differ.

The Baseline
Our method (equation (1)) is designed to detect the FCDs exclusively and as well as possible. After a FCD, the linear prediction becomes contaminated by temperatures of the cold period. Therefore, the predicted values, i.e. the baseline, tend to show decreasing values (Figure 1) after the FCD. This is further illustrated in Figure 5, where we compare our analysis of the very strong Kuwae event (1453 AD) with the analysis in Esper et al. (2017). Note that in our case the years are defined to be cold whenever the nonpredictable temperature is negative.
The result in Esper et al. (2017) is based on an average over reconstructions from 20 sites and their computational technique differs from ours. Nevertheless, abrupt cold periods of both analyses agree well in the Kuwae case. In contradiction to that, the preceding (after 1435 AD) analogical strong cold period is not noticed in Esper et al. (2017), and, in fact, is not mentioned among the strong cases analyzed in the literature. Suitable eruption candidates are unknown or only very weak (1435 AD; Crowley & Unterman, 2013). That eruption is not visible in Figures 3b and 5 because the closest FCD follows only in 1438, that is, only after three years.
In the case of successive abrupt cold periods like in Figure 5, the baseline is continuously decreasing. Nevertheless, the information on the initial cold periods of the temperatures is not lost because the reconstruction is precisely the sum of the predictable (baseline) and nonpredictable parts. These parts seem to project the reconstructed temperatures during volcanic periods into two components resembling slow and fast modes. Such a property was not expected.
The oscillations of the baseline are clearly longer than the abrupt cold periods. For instance, the cooling baseline period starting with the FCD in 1095 AD (Figure 3a) includes three additional FCDs (1107, 1117, and 1127 AD). Altogether eight eruptions are listed for that period. The baseline decreases smoothly passing the abrupt cold periods. After several successive cold periods a minimum is reached and around 1175 AD the baseline slowly and smoothly recovers back to the initial temperature level. There are no FCDs during the warming phase and the baseline consistently reflects the course of the temperature reconstruction, although slightly delayed.
Similar quasi-periodic long variations were found in Rinne et al. (2014) by directly subtracting the abrupt cold periods caused by the known eruptions from a corresponding temperature reconstruction. No physical explanation of the result could be given. However, similarity to the oscillations of the late-summer seaice presented by Kinnard et al. (2011) was found. Accordingly, in the following we will compare our present baseline with their sea-ice oscillations. The baseline and sea-ice oscillations seem to be related to each other during ≈700-1550 AD. To study that more accurately, we scaled the sea-ice oscillations to fit with the temperatures of the cooling phase of the Kuwae event 1435-1470 AD ( Figure 5). The result is tested in an independent period of 700-1240 AD ( Figure 6) and we can observe the two curves following each other. The resolutions of the curves are different, being ≈40 years in the sea-ice construction, whereas the baseline temperatures are only implicitly smoothed by our method. However, because of the linear prediction method applied, the baseline tends to be delayed. Values in Kinnard et al. (2011) refer to a wider oceanic area and may therefore be subject to different variations compared to our local area in Fennoscandia and nearby seas. These features cause inaccuracies in the results so that a precise comparison between the curves cannot be made. For instance, the seaice seems to lead the volcanic oscillations of our baseline ( Figure 6). This is especially striking in the case of the Kuwae event with its hemispheric-wide impact area (Esper et al., 2017), as the volcanic forcing should be the leading factor. Hence, the sea-ice oscillations and the baseline show together a certain systematic, physically grounded behavior, although the accuracy of the data and computations do not allow for a more precise description.
The long oscillation after 1095 AD discussed above is intense (Figure 3a), the baseline decreasing by 1.5°and recovering then back to the initial level. The length of the oscillation is approximately 75 years. To study that symmetry in more detail, we illustrate in Figure

Abrupt Cold Periods and Eruption Data
We have found 79 abrupt cold periods by identifying their leading cold anomalies. Our tools, i.e. the baseline and prediction error variance, have made the high number possible, although only one local temperature reconstruction is studied. No eruption data have been utilized. Therefore, we can present only candidate volcanic events that would fit with the cold periods.
As far as we know, only volcanic eruptions can, in our Fennoscandian case, repeatedly cause abrupt strong cold peaks followed by several cold years. Therefore, we can in principle assume our 79 cold periods to be postvolcanic ones. Candidate eruptions are not found for 25 FCD cases during the 2000 years studied, of which 20 took place before 1000 AD.
To study that more closely, we picked all FCDs during 1750-2000 AD. It is seen that the data from recent centuries about the temperatures and eruption years are more precise than those from earlier times. Further, there are historical records available, too. We studied all cases in Figure 4 and found two FCDs without any indicated eruption, those of 1821 and 1949 (Figure 4). We find potential Icelandic candidates for these events, namely, Katla (1821) and Hekla (1947;  e.g., Watson et al., 2017). Ash fall of Hekla was observed in Fennoscandia (Salmi, 1948). In some cases the FCD is observable before the eruption, like in 1768 and 1962 (the latter preceding Mount Agung 1963). These cases could again be explained by Icelandic eruptions in 1766 (Hekla) and 1961 (Askja), respectively. The impact of Icelandic eruptions can become observable mainly in downstream areas, that is, in Northern Fennoscandia. Taking them into account, we can find candidate eruptions for all of the FCDs during 1750-2000. One specific case from earlier years is shown in Figure 5. As far as we know, the cold period with the FCD 1438 is not studied in the literature. It is, however, of the same structure and nearly as strong as the following cold period of the Kuwae case and can therefore be assumed to be an analogical event.
We conclude that the main part of our 79 FCDs, if not all, present postvolcanic cases considering the whole 2,000-year-long study period. In addition to them, in our approach many cases are not noted because only the leading FCDs are taken into account and some obvious cases had to be rejected like Pinatubo event (1991) or the very strong cold drop 536 AD (Figure 2b).

The Local Baseline Oscillations
Initially, the baseline values were constructed only to detect the years of the FCDs. Unexpectedly, the bidecadal-multidecadal oscillations of the final baseline turned out to include additional information about the impact of the volcanism.
In our approach, no forcing terms have been applied and so only the input data, a reconstruction of the Fennoscandian summer temperature, are used. Therefore, we cannot directly say anything about the actual contributions of the various applicable forcing effects. However, the importance of the volcanism can be deduced from the occurrence of the FCDs. We counted their numbers during different phases of the baseline. There were 67, 5, and 7 FCDs during the cooling, intermediate, and warming phases of the baseline, respectively. Thus, the abrupt cold periods concentrate in the cooling phases; that is, the volcanism is strictly connected to the baseline oscillations. When counting them, the baseline was shifted six years backward to take into account the lag of the baseline. In the five stronger oscillations of Figure 7, there are 14 and 0 FCDs during the cooling and warming phases of the baseline oscillation, respectively. Further, the cooling phases tend to start with a FCD and so the volcanism clearly is a principal contributor in the baseline oscillations leading to cold deviations.
In the following, the number and length of the baseline oscillations are studied case by case during two millennia. We measure the lengths of the baseline oscillations by determining the distance between the years of successive wave crests. In the following, only distances longer than 15 years are considered. The periods of 565-730 AD (Figures 2b and 2c) and 1275-1375 AD (Figure 3b) without any intense baseline oscillations will be omitted.
We find 42 such oscillations of the baseline during 50-1980 AD, the longest one starting in 413 AD ( Figure 2b) and being about 90 years in length. The total and average lengths of the oscillations are ≈1,665 and 40 years, respectively. In addition, the total length of the periods without any intense baseline oscillations becomes 265 years. Thus, 86% and 14% of the years 50-1980 AD are dominated by the presence and absence of the baseline oscillations, respectively. As the role of the volcanism is clear in the longbaseline oscillations and these dominate, it follows that the volcanism is a significant source of bidecadalmultidecadal variations during the two millennia.
The atmospheric and oceanic responses to the volcanic forcing are known to be short and long lasting, respectively (e.g., Swingedouw et al., 2015). Compare here also with Otterå et al. (2010), "… explosive volcanism has a strong influence on NH climate, not only for short term changes but also for multidecadal AMOtype changes." Our abrupt cold periods and long-baseline oscillations reflect thus the direct atmospheric and indirect oceanic volcanic impacts evident in the temperature reconstructions. Although the temperatures are continental, we nevertheless can see traces of the oceanic variations arising from adjacent arctic waters. As we are studying temperatures in the very north (≈70°N), the temperature variations in the ocean should have a further impact on the adjacent sea-ice extent. Thus, the order of the events in our scheme is initiating volcanic impacts, followed by temperature changes directly via atmosphere and indirectly via ocean and finally the sea ice variations. The outcome resembles that in modeling experiments in Slawinska and Robock (2018): "Volcanic eruptions are followed by a decadal-scale positive response of the Atlantic multidecadal overturning circulation, followed by a centennial-scale enhancement of the Northern Hemisphere sea ice extent" and further in PAGES 2k Consortium (2019): "Volcanic eruptions are found to coincide with strong multidecadal cooling trends in both reconstructions and model simulations, and this is then followed by an increased probability of strong warming trends, due to the recovery from the volcanic cooling." Further, Neukom et al. (2019) conclude, "A substantial portion of pre-industrial (1300-1800 CE) variability at multi-decadal scales is attributed to volcanic aerosol forcing." Haustein et al. (2019) continue to modern times and suggest that multidecadal ocean variability was unlikely to be the driver of observed changes in global mean surface temperature after 1850 A.D. Instead, virtually all (97-98%) of the global low-frequency variability (>30 years) can be explained by external forcing.
As the volcanic impacts can vary over different areas, we cannot derive conclusions valid over wider areas from our local case. For instance, the Tambora event (1815) is negligible in our case, although otherwise its impact was exceptionally strong on both sides of the Northern Atlantic. The Samalas eruption mentioned in the introduction serves as a similar example of varying impact areas. Further, Neukom et al. (2019) found no evidence for preindustrial globally coherent cold and warm epochs. Accordingly, the abrupt ice growth in Arctic Canada between 1275 and 1300 AD (Miller et al., 2012) may be linked to variations in our baseline, although we cannot study that because of remotely located areas.
The sea-ice oscillations in Kinnard et al. (2011) are determined for a wider area. Therefore, it is not necessary that their sea-ice oscillations are directly linked to our baseline variations. Like in the case of the Tambora eruption, impacts of the Kuwae event 1453 are observable on both sides of the Northern Atlantic (Esper et al., 2017). However, in the latter case the cooling was strong in our results, too, and so the impact can be expected to be visible in the near-by sea-ice variations. Thus, the sea-ice oscillations in Kinnard et al. (2011) can be expected to be at least temporarily connected to the Fennoscandian temperatures during certain periods, in our case during 700-1550 AD.
To summarize, the local temperature variations in Northern Fennoscandia are mainly caused by the volcanism and are seen as direct atmospheric and indirect oceanic impacts. We further find connections to the seaice variations in a multicentury time frame.

Pros and Cons of the Searching Method
Our method provides the means to detect the FCDs local to Northern Fennoscandia. They do not depend on the intensity of the potential eruption but on the strength of the impact observed. The high number of FCDs has helped to characterize bidecadal-multidecadal variations of the local baseline as the volcanic impact via arctic seas and to discuss those variations case by case.
The baseline generally gives a smoothed version of the temperature reconstruction but not precisely because the baseline oscillations tend to lag those of the temperature. Although the baseline is designed to be in the FCD years as accurate as possible, during the following cold years it may be less applicable. The baseline values can be unexpectedly high, that is, after 413 or 1902 AD, so the separation between the atmospheric and oceanic responses is only tentative during the abrupt cold periods. The length of the atmospheric response time tends anyway to be rather short.
A peculiarity of Northern Fennoscandia is that no FCDs could be found for the strong cases of 1257 (Samalas), 1809 (unknown), and 1815 (Tambora). In some cases a FCD is noticed, but nevertheless no suitable candidate eruption is found in the lists applied, that is, for oscillations after 830 ( Figure 2c) and 1175 AD (Figure 3a). Such deviations may obviously be caused by natural variations; otherwise, we can only speculate about the reasons. At 830 no strong cold peak is seen in Figure 2c; that is, the atmospheric response turns out to be weak. A possible explanation is that an oceanic response is nevertheless due to the advection of a cold temperature anomaly from a more remote area of the volcanic impact.
As only one temperature reconstruction has been studied, there are limitations in the results. FCDs may remain hidden: • if there is additional noise in the temperature reconstruction caused by strong weather variations or computational inaccuracies; • if a cold peak is observable but does not pass our significance test, that is, in the case of Pinatubo; • if the cold peak passes the test but is single and is thus not followed by a cold anomaly, that is, 536 AD; • if there are several FCDs within one continuous cold period, then only the first one is noted.

Conclusions
A reconstruction of summer temperatures during two millennia in Northern Fennoscandia has been studied and 79 abrupt cold periods have been identified. Data regarding the volcanic eruptions are not a priori used. Nevertheless, most of the cold periods found can be connected to eruption candidates, especially during the most recent 250 years. Exceptionally, the strong eruptions of 1257 (Samalas), 1809 (unknown), and 1815 (Tambora) seem to have no significant direct impact on the Fennoscandian temperature climate. On the other hand, local impact of weaker, less global (Icelandic), or unknown eruptions is evident in the results.
In addition to the short cold periods, there are longer bidecadal-multidecadal temperature oscillations. We found 42 such oscillations during 50-1980 AD, the longest one lasting about 90 years. Their total length is 1,665 years. Two periods without any clear longer temperature oscillations were understood to be periods of a less intense volcanic impact. Their total length is 265 years. Hence, the volcanism has clearly been a dominant contributor to local temperature variations in Northern Fennoscandia during the two last millennia. Multicentury similarity to the sea-ice oscillations was observed. The short abrupt cold periods and the longer oscillations were characterized as the atmospheric and oceanic responses to volcanic impacts, respectively.