Reduced accuracy in dendroglaciological mass balance reconstruction of Storglaciären since the 1980s

Abstract

Recent studies have raised concerns regarding the reconstruction of glacier mass balance using tree-ring data. The method relies on a stable relationship between both variables and summer (June–August) or melt season (May–September) temperature. However, with recent anthropogenic climate change the stability of this relationship is challenged by both, a divergence between tree-ring and temperature, as well as mass balance and temperature data. Establishing to what extent this divergence influences the reconstruction of mass balance using tree-ring data is important to assess the future viability and applicability of the method. In this paper we analyze the relationship between the Torneträsk tree-ring and Storglaciären mass balance records, their response to climate change, and investigate changes in the relationship. We provide evidence for a sensitivity loss in the Torneträsk record and quantify its impact on the reconstruction of summer mass balance of Storglaciären. We find that by including years post 1980, the amplitude of reconstructed variability is reduced, but it remains possible to explain the variance of the record significantly. Our results suggest that for glaciers without an extensive mass balance record the applicability of the method is challenged.

Introduction

Long-term mass balance measurements provide insight into past responses of glaciers to climatic changes and can thereby help to assess the future impact of increasing temperatures on glacier mass loss. To date, the longest continuous mass balance record of annual measurements is the one of Storglaciären in northern Sweden, dating back to 1946 (e.g. Holmlund and Jansson, 2002). The glacier’s net balance is comprised of in-situ measurements of accumulation and ablation, of which the latter is mainly determined by summer temperature (Holmlund et al., 2005), here defined as the mean June to August $(T_{\mathrm{JJA}})$ temperature. The linear relationship between the summer mass balance $(B_s)$ and mean summer air temperature allows for the reconstruction of $B_s$ using temperature records (Karlén, 1984).

While historical meteorological records provide important information on previous climatic conditions, they rarely date back further than 1850 CE and are sparsely distributed over the globe (Peterson and Vose, 1997). Attempts to reconstruct temperatures further back in time therefore rely on valuable incremental proxies such as tree rings, ice cores or varves, which reflect climatic conditions with high temporal resolution. In recent decades, especially tree rings have attracted scientific interest for the reconstruction of past climates (e.g. Briffa, 1999; Zaw et al., 2020). Different tree-growth indices, namely tree ring width (TRW) or maximum latewood density (MXD) have been used in multiple instances to reconstruct temperatures, for example for Scandinavia (Grudd, 2008; Gunnarson et al., 2011; Linderholm et al., 2015; Svarva et al., 2018) and the entire northern hemisphere (Anchukaitis et al., 2017; Esper et al., 2018; Wilson et al., 2016). Especially MXD records enable accurate reconstructions of temperature, as its correlation with temperature during the growing season (June to August) exceeds that of TRW (Esper et al., 2015; Grudd, 2008). In Fennoscandia, the extensive set of MXD measurements from the Torneträsk area provides the base for one of the main temperature reconstructions of the region (Grudd, 2008, Figure 1).

Figure 1.

Map of the study area in the northern Scandinavian Mountain Range. Displayed on the left are Storglaciären’s mass balance record and the Torneträsk chronology with temperature reconstruction by Grudd (2008). All data for the base map is provided by Lantmäteriet.

The reconstruction of $B_s$ from tree-ring records relies on the linear relationship of both parameters to summer air temperatures and has been applied in various regions including the Alps (Cerrato et al., 2019; Leonelli et al., 2011), Fennoscandia (Karlén, 1984; Linderholm et al., 2007; Linderholm and Jansson, 2007), and North America (Larocque and Smith, 2005; Luckman, 1996; Watson et al., 2006; Watson and Luckman, 2004; Wood and Smith, 2013). However, recent research suggests that with steeply increasing temperatures in the context of anthropogenic climate change, the northern hemispheric relationship between temperature and MXD is altered with sensitivity of MXD to temperature decreasing since approximately the 1980s (Cerrato et al., 2019).

This impairment in the relationship between MXD and temperature is known as the divergence problem and has been described in multiple instances for different tree species (Briffa, 2000; Briffa et al., 1998a, 1998b, 2004; Driscoll et al., 2005; D’Arrigo et al., 2004, 2006, 2008; Esper and Frank, 2009; George and Esper, 2019; Jacoby et al., 2000; Jacoby and D’Arrigo, 1995; Vaganov et al., 1999; Wilson and Luckman, 2003). First discovered by Jacoby and D’Arrigo (1995), it has since not been explicitly verified for the Torneträsk MXD record, since recent warming has not been as pronounced in Fennoscandia as elsewhere in the northern hemisphere.

The reason for the sensitivity loss of MXD is still debated. Schweingruber (1988) noted that both MXD and TRW are subject to physiological limitations at certain temperature thresholds that vary between tree species and may play a role in the divergence. Additionally, along with climate in general, insect and moisture stress induced by recent warming trends, as well as increases of carbon dioxide, nitrogenous compounds or ultra violet radiation, and the timing of snow melt (Briffa et al., 1998b; Taubes, 1995; Vaganov et al., 1999) may all influence the pattern of divergence. Unfortunately, this divergence might impair the reconstruction of previous climate variations (Loehle, 2009) and accordingly also that of $B_s.$ It is therefore important to reassess the existence and magnitude of a similar loss in MXD temperature-sensitivity in Fennoscandia and to quantify its possible impact on mass balance reconstructions of Storglaciären.

In this study we update the existing reconstructions of Storglaciären $B_s$ (Linderholm et al., 2007; Linderholm and Jansson, 2007) by introducing additional measurements of MXD and $B_s$ for the period 1980–2010 into the analysis. We further reassess the applicability of the method in the face of a decrease in the strength of the MXD– $B_s$ relationship in the Torneträsk chronology and quantify its impact on reconstructions of $B_s$ of Storglaciären.

Study Area

Located approximately 150 km north of the Arctic Cirlce in the Kebnekaise massif, northern Sweden, Storglaciären terminates in the Tarfala Valley ( $67.{91}^o\mathrm{N,}$ $18.{61}^o\mathrm{E,}$ Figure 1). The regional climate of the area is influenced by the North Atlantic Oscillation and due to the rainshadow of the Scandinavian Mountain Range in the West is more continental than the climate on the western side of the mountain range (Chen and Hellström, 1999). However, the relative proximity to the coast is still reflected in a milder climate than shown for other locations at this latitude (Linderholm et al., 2007).

With the long continuous mass balance record of Storglaciären dating back to 1946 and the meteorogical measurements from the Tarfala Research Station, the Tarfala valley provides an ideal study area. Storglaciären’s mass balance (Figure 1) is known to be well correlated with other glaciers from the Tarfala Mass Balance Program and can therefore be seen as representative of the wider area (Linderholm et al., 2007), even though this relationship might have decreased in recent years (Brugger, 2007).

Data and methods

Mass balance data

With over 70 years of data, the record of mass balance measurements at Storglaciären (Figure 1) is to date the longest continuous record of distributed mass balance measurements in the world (measurements are described in detail in Holmlund and Jansson, 1999). Its net balance is determined by separate consistent measurements of summer and winter balance at the beginning and end of each melt season (May to September). However, with MXD being determined mainly by summer temperatures, we focus on $B_s,$ expressed in meters of water equivalent (m w.e.). The record is compiled from measurements at multiple locations throughout the glacier. Values in the observed $B_s$ (1946–2010, obtained from the Bolin Centre database) range from 0.67 to 2.90 m w.e. (mean: 1.67 m w.e.).

Tree-ring and meteorological data

The Torneträsk MXD chronology is one of the longest continuous MXD records in the world, covering the period from 441 CE to 2010 CE (Grudd, 2008; Schweingruber et al., 1988). The sampling site is located at Lake Torneträsk, approximately 60 km northeast of Storglaciären (Figure 1). We utilize the time period from CE 1500–2010 for the $B_s$ reconstruction. Until 1980, the chronology is comprised of MXD series from 65 living and fossil specimens of scots pine (Pinus sylvestris), of which fossil trees represent the period of 441–1400 CE supported by living trees from 1400 onward. For the period 1980–2010 more than 35 additional series from living trees extend the record to its full length (Briffa et al., 1990; Grudd, 2008).

To use summer temperature as the main denominator between tree growth and glacier mass balance, we obtained temperature measurements from Tarfala Research Station for 1946–2010 from the Bolin Centre database. Mean monthly values are available for the summer season (JJA) only since 1964, and for all months since 1965, when automated data collection was first established at the station. With the Torneträsk record being significantly correlated with climate data from the entire Fennoscandian North, there is no conflict in connecting Torneträsk tree growth to summer temperatures in the Tarfala Valley (Briffa et al., 1992; Gouirand et al., 2008; Grudd, 2008). We thus utilize the MXD record, standardized using signal-free regional curve standardization to preserve low frequency variability in the data, to attempt an updated reconstruction.

Regression analysis and mass balance reconstruction

To assess the impact of the impaired relationship between MXD and $T_{\mathrm{JJA}}$ on mass balance reconstructions of Storglaciären, we examined temporal changes in the relationship of both MXD and $B_s$ with temperature. In this context we correlated the Torneträsk MXD chronology and Storglaciären’s $B_s$ separately with monthly mean temperatures from Tarfala Research Station for the period 1946–2010. Additionally, seasonal means for the summer (June-August) and melt (May-September) season ( $T_{\mathrm{JJA}}$ and $T_{\mathrm{MJJAS}},$ respectively) were included in the analysis. A Spearman’s rank correlation coefficient $(\rho )$ was used to identify the month or season when both MXD and $B_s$ show the strongest correlation with temperature. We accounted for the possibility of a lagged response of MXD and $B_s$ to temperature changes by introducing time lags of up to 20 years into the correlations analysis. Subsequently, the variable with the highest correlation with MXD was used as the predictor in separate least-square linear regression models with $B_s$ and MXD as response variables in order to quantify the explanatory power of temperature over both response variables, respectively. To account for temporal variations, the linear regression models were set up using a 30-year moving-window approach with annual steps. In doing so, changes in the relationship of temperature with both variables over time can be investigated by deriving the coefficient of determination $(r^2)$ for each 30-year window. We used bootstrapping, a commonly used resampling method to derive estimates of confidence intervals by random sampling with replacement from a dataset to obtain confidence intervals for each 30-year window. Because of seasonal variability in the time series data, a block bootstrapping approach with a random sample size of 2000 was favored over normal bootstrapping. In this approach, samples are drawn with replacement from a set of independent blocks with a length generated from a geometric distribution (Berkowitz and Kilian, 2000). Based on the fact that each window’s model shows significance $(\alpha =0.05)$ and to stay consistent with previous reconstructions attempts of glacier mass balance of Storglaciären (Linderholm et al., 2007; Linderholm and Jansson, 2007), we used a simple linear regression model to predict $B_s$ based on MXD. The reconstruction was validated by testing whether the reconstructed distribution is significantly different from the measured $B_s$ using a non-parametric Mann-Whitney $U$ test (e.g. Birnbaum, 1956). Furthermore, the accuracy of the reconstruction was evaluated by calculating the root mean square error (RMSE). The time series were split into the periods 1946–1980, 1981–2010 and the entire record from 1946 to 2010 to investigate the change in sensitivity of the relationship between the variables due to changes in climate since the 2007 reconstructions (Linderholm et al., 2007; Linderholm and Jansson, 2007).

Results

Monthly mean temperatures for July and August, as well as the seasonal mean for the summer and melt season correlate significantly with both MXD and $B_s$ for all investigated time periods (1946–1980, 1981–2010, and 1946–2010). No significant correlation with either MXD or $B_s$ was found for the individual mean temperatures of May and September. The highest correlation with both parameters for all seasons (except for the 1981–2010 $B_s$ record where it is $T_{\mathrm{MJJAS}}$ ) was found with $T_{\mathrm{JJA}}$ . While $B_s$ and MXD show very strong correlations for 1946–1980 ( $\rho =0.85$ and $\rho =0.9,$ respectively), the coefficients over the later period (1981–2010) are notably lower ( $\rho =0.71$ and $\rho =0.7,$ respectively). Note the decrease of 0.2 for MXD.

The decrease in the strength of the relationship of the two parameters is evident in the $r^2$ of the linear regression models (Figure 2). The necessity of introducing a time lag into the regression analysis was discarded, since we found no significant lagged response of the variables to temperature changes. While summer temperature is able to explain $B_s$ in a stationary fashion until the 1988 window, the temperature sensitivity of the MXD follows a downward trend over the entire time period. However, even with the identified decrease in the correlation over time, $T_{\mathrm{JJA}}$ is able to explain MXD with an $r^2$ of 0.64.

Figure 2.

(a) Centered 30-year moving-window for each time step (step $=1\mathrm{yr}$ ) for the linear regression between $T_{JJA}$ and MXD, $T_{JJA}$ and $B_s$ and MXD and $B_s.$ Shaded areas indicate the bootstrapped confidence intervals of each regression. Years on the y axis indicate the central year for each 30 year time window. (b) Spearman rank correlation coefficients for the moving windows.

The observed drop in explanatory capacity of $T_{\mathrm{JJA}}$ on MXD reaches its lowest point in the 1986 window ( $r^2=0.49$ ), strongly impacting the correlation between $B_s$ and MXD (lowest $r^2=0.28$ ). Nevertheless, in all three regressions the explained variance remains significant over the entire period (Figure 3). Testing the difference between reconstruction and observation using a Mann-Whitney $U$ test (1946–1980: $p=0.96;$ 1946–2010: $p=0.87$ ) revealed no significant stochastic difference between the distributions for the overlapping period of 1946–2010. However, the accuracy of the prediction decreased with the introduction of post 1980 years, showing RMSEs of 0.29 m w.e. (1946–1980) and of 0.33 m w.e. (1946–2010).

Figure 3.

(a) Reconstructed $B_s$ with 10-year rolling mean (saturated lines) of both reconstructions back to CE 1500 (desaturated lines) and (b) both reconstructions compared to the 1946–2010 observations.

Discussion

Correlation analysis and divergence

Judging from our linear regression, a decrease in temperature sensitivity of MXD over the entire correlation period becomes especially clear with the introduction of more recent years into the moving windows. While other studies also mention this period as the starting point of the divergence problem (e.g. Briffa et al., 1998b), none of these studies indicates the existence of a pronounced divergence in the Torneträsk MXD record, especially after the introduction of younger trees in the most recent tree-ring record to lower the mean cambial age (Grudd, 2008). The question is thus raised, whether the divergence found in our data provides first evidence of the divergence problem reported in other studies or if some other phenomenon causes the relationship between MXD and summer temperature to decline. Having additionally tested the regression for other temperature and MXD records of the area, however, we do not believe the decline is a phenomenon caused by the choice of record, and thus argue that it does show some overall resemblance to the divergence problem observed in North America (D’Arrigo et al., 2004; Wilson and Luckman, 2003) and the European Alps (Büntgen et al., 2006, 2008; Frank and Esper, 2005a, 2005b). With the inclusion of additional young trees for the extension to 2010, there should further be no concerns regarding an apparent reduction in MXD temperature sensitivity that can be connected to a high mean cambial age toward the more recent part of the record (Grudd, 2008).

Another possible cause for the divergence in the second half of the twentieth century could be the varying annual sample size of the Torneträsk chronology. Toward the end of the chronology, fewer trees are included, possibly lowering the signal-to-noise ratio and influencing the correlation. As the sample size decreases, the risk for erroneous results in the tree-ring data is likely to increase. However, since a sensitivity loss is detectable as early as the 1960s, where a representative number of samples is still given, this might not be the determining factor for the relationship’s decline.

A loss in sensitivity to summer temperature during recent decades can also be seen in the relationship of $B_s$ and $T_{\mathrm{JJA}},$ especially after the 1977–2006 window. A similar decrease in correlation was reported by Cerrato et al. (2020) for the Careser glacier in the European Alps when the outlier of 2003 was removed from the dataset. This speaks for a more widespread change in the response of glaciers to climate which Cerrato et al. (2020) attribute to the decay of the Careser glacier in their case. Thus, the reduction in the temperature sensitivity of Storglaciären could potentially too stem from the recent mass loss of the glacier (Koblet et al., 2010), resulting in a non-linear response to climate.

Influence of divergence on the mass balance reconstruction

Our reconstruction shows a similar synchronicity between reconstruction and observation as the reconstruction by Linderholm et al. (2007). Since their study only includes values until 1980, their reconstruction may more closely resemble the observed data than our extended version. Our $B_s$ reconstruction further shows a somewhat lower mean (1.55 m w.e.) than the reconstruction by Linderholm et al. (2007) (1.67 m w.e.) as well as a higher standard deviation (0.46 compared to 0.34 m w.e.) and higher variance (0.22 as opposed to 0.12 m w.e.). The difference in the statistics may arise from the calibration of the linear regression equation parameters or tree-ring standardization method applied by Linderholm et al. (2007).

The underestimation of the 1946–2010 based reconstruction can likely be traced back to the weaker correlation of MXD and $B_s$ toward the end of the time series and exemplifies the influence of the sensitivity loss of MXD during the late 20th century on the reconstruction. Note that in a least-square linear regression the amplitude of the reconstruction decreases by the square root of the unexplained variance (Esper et al., 2005). Since the reduction in explained variance since the 1980s can be attributed to the divergence between summer temperature and both $B_s$ and MXD, we consider the reduced amplitude in the reconstruction a consequence of the divergence. Although more sophisticated scaling models have been shown to better capture the amplitude (Cerrato et al., 2020), one aim of this study was to evaluate whether the methods previously used to reconstruct $B_s$ of Storglaciären are still applicable. Despite the reduced amplitude being a consequence of the linear regression, we find our results to be in agreement with the findings by Cerrato et al. (2020). In their reconstruction, the authors present a higher $B_s$ when using data until the 1980s, compared to subsequent decades. They further emphasize the leverage of outliers in the relationship between $B_s$ and temperature. However, we did not find similar notable changes to the correlation after removing outliers, likely due to our time series being significantly longer.

Apart from the model uncertainty, differences and uncertainties in the standardization of the tree-ring data and chronology development, as well as calibration methods could potentially also result in reconstruction differences (Cook and Peters, 1997; D’Arrigo et al., 2006).

No matter any slight disagreements between previous reconstructions of Storglaciären’s summer mass balance using MXD and our present attempt, it becomes clear that an introduction of more recent years into the reconstruction will alter the ability to detect extremes and reflect less of the variability of the record. By the use of more sophisticated scaling models this reduction in amplitude may be circumvented at the expense of variance precision (Esper et al., 2005).

Glacial conditions captured by the reconstruction

The reconstruction of $B_s$ back to 1500 CE allows for some comparison with local climatic events to assess whether such events are captured by our reconstruction. Notably, cold episodes during the Little Ice Age (Fuentes et al., 2018) in the early and mid 1600s, as well as during the early 18th century are captured by our reconstructions in the form of low $B_s$ values. They are much more pronounced, however, in the 1946–1980 based time series as anomalously low or high values are smoothed out by including the years of lower sensitivity. Generally lower $B_s$ values like those expected for the early 20th century, which was characterized by glacial advances (Linderholm and Jansson, 2007), are underestimated in the 1946–2010 reconstruction. Conversely, they are displayed in the 1946–1980 based reconstruction. Still, the subsequent retreat of glaciers after 1910 (Holmlund and Jansson, 1999) is depicted well by both reconstructions as an increase in $B_s.$ Accelerated mass loss in the 1940s as described by Huss et al. (2015) for the Swiss Alps is further represented in the reconstructions by elevated $B_s$ values. Thus, although the 1946–2010 reconstruction shows lower amplitudes, it is still able to properly represent temperature variability in the reconstructed $B_s.$

Conclusions

In this study we reassess the reconstruction of $B_s$ of Storglaciären in the face of an impaired relationship between mean summer air temperatures and MXD, as well as $B_s$ in connection with current climate change.

We are able to provide evidence of a divergence for the Torneträsk chronology, identifying a loss of sensitivity in the $T_{\mathrm{JJA}}-\mathrm{MXD}$ relationship since approximately the 1980s. This divergence has a clear adverse effect on the MXD based reconstruction of $B_s,$ reducing the explained variance from $r^2=0.63$ (1946–1980 based reconstruction) to $r^2=0.48$ by including the years 1981–2010 into the analysis. The reconstruction of $B_s$ is further influenced by a reduction in temperature sensitivity of $B_s$ itself. Thus, future reconstruction efforts are hampered by both divergences and should be attempted with caution. Utilization of more sophisticated models than linear regression may help to reduce the loss in reconstructed $B_s$ amplitude introduced by the use of less sensitive years for introduction, especially for areas where long term observations are not available for comparison. A simple linear regression model would likely not be able to capture the observed and therefore past $B_s$ extremes satisfactorily and generally underestimate $B_s$ observations.