Varying contributions of greenhouse gases, aerosols and natural forcings to Arctic land surface air temperature changes

The Climatic Research Unit gridded time series (CRU) is a widely used climate dataset that covers most land domains worldwide. It is generated through the interpolation of monthly climate observations from more than 4000 weather stations (Harris et al 2020). In this study, the monthly mean temperatures from 1915 to 2014 are adopted from the CRU dataset (version 4.05) at a spatial resolution of 0.5°× 0.5° (https://crudata.uea.acuk/cru/data/hrg/cru_ts_4.05). To further evaluate the robustness of our findings, we repeat our analysis with two additional observational datasets: the Global land average temperature dataset from the University of Delaware (DGLT) and the Goddard Institute for Space Studies (GISS) surface temperature dataset. Detailed introductions of DGLT and GISS can be found in supplementary section 1.

To quantify the effects of natural and anthropogenic forcing agents, we use Arctic temperature simulations by the CMIP6 models driven with historical variations in GHGs, NATs and all forcing agents combined (ALL) (Eyring et al 2016). The OANT forcings are obtained by subtracting the simulations under GHGs and NATs from that under ALL forcing agents. After checking the CMIP6 data archive at https://esgf-node.llnl.gov/projects/cmip6/, only nine models (i.e. ACCESS-CM2, ACCESS-ESM, CanESM5, CESM2, GFDL-ESM4, GISS-E2-1-G, HadGEM3-GC31-LL, MRI-ESM2-0, and NorESM2-LM) provide the required simulations for the period 1915–2014 and are therefore selected for our detection and attribution analysis (supplementary table 1).

Previous studies usually assigned equal weights to each climate model, based on which the ensemble mean is used for analysis. However, different models often show different accuracies, thus pointing to the need for assigning different weights to each model. Here, based on Ridge Regression (RR), we obtain different weights for the nine climate models via building statistical relationship between simulated temperature under ALL forcing and observed temperature (supplementary sections 2, 3 and figure 1). The regression coefficients and weights of the nine models are presented in supplementary table 2. It should be noted that our weighting approach depends on the variable of interest, and it is possible that the CMIP6 models which are good in reproducing observed temperature are poor in terms of other climate variables. Similarly, we build the RR for each pixel and calculate the pixel-level weights of the nine climate models, based on which the multi-model weighted mean is obtained for each pixel. Before training the pixel-level RR, we resampled all model simulations to the resolution of 0.5°× 0.5° to be consistent with the observation dataset based on bilinear interpolation technique.

It is well known that changes in temperature are highly affected by naturally generated IV. In this study, we adopt the CSA to estimate the IVs of climate models, based on their preindustrial control simulations with constant external forcings (Resplandy et al 2015, Roberts et al 2015, Olonscheck and Notz 2017).

To obtain more reliable IV estimates, we first remove the least-squares linear trend from each model's preindustrial control simulation to minimize the effect of model drift. Here, the liner trend indicates the trend of each 100 year segment for preindustrial control simulation. We then evaluate the sensitivity of IV to the length of control simulations based on the power spectrum method (supplementary section 4). If the IV given by the control simulation is not equal to the total IV of this particular model, then such models are discarded. After performing the above two steps, a total of 31 032 years of preindustrial control simulations (i.e. piControl) by 48 CMIP6 models (supplementary table 3) are finally used to estimate the IV. We have also calculated the IV based on the nine models used for attribution analysis, though using all available models for calculating IV is well adopted by previous studies (Najafi et al 2015, Mueller et al 2018). Overall, similar estimate of IV is obtained with the relative difference ranging from 2.2% to 5.7% (supplementary figures 2 and 3), suggesting the robustness of our analysis.

In this study, we adopt the well-established approach of total least-squares optimal fingerprinting (TLSOF) technique for the detection and attribution of Arctic land surface air temperature changes to GHGs, NATs and OANT forcings (Hasselmann 1997, Allen and Stott 2003, Stott et al 2010, Hegerl and Zwiers 2011). TLSOF uses a generalized linear regression model, which assumes that the observed Arctic land temperature changes are determined by a linear combination of the covariates (GHGs, NATs and OANT agents) along with a random term. The temperature observations can be expressed as follows:

$$\begin{equation}{T_{{\text{obs}}}} = {\beta _1}{T_{{\text{ALL}}}} + {\beta _2}{T_{{\text{NAT}}}} + {\beta _3}{T_{{\text{GHG}}}} + \varepsilon \end{equation} \tag{ 1 }$$

where ${T_{{\text{obs}}}}$ is a vector of observed Arctic land temperature anomalies and ${T_{{\text{ALL}}}}$, ${T_{{\text{NAT}}}}$, and ${T_{{\text{GHG}}}}$ are estimates of the temperature anomalies in response to ALL, NAT and GHG forcings, respectively, with ${\beta _1},{\beta _2},{\beta _3}$ indicating the corresponding covariate coefficients. ε is the residual variability that is generated internally in the climate system. Since ${T_{{\text{ALL}}}}$ can be decomposed into ${T_{{\text{ALL}}}} = {T_{{\text{GHG}}}} + {T_{{\text{NAT}}}} + {T_{{\text{OANT}}}}$, equation (1) can be rewritten as follows:

$$\begin{align}{T_{{\text{obs}}}} & = \left( {{\beta _1} + {\beta _3}} \right){T_{{\text{GHG}}}} + \left( {{\beta _1} + {\beta _2}} \right){T_{{\text{NAT}}}} \nonumber \\ & \quad+ {\beta _1}{T_{{\text{OANT}}}} + \varepsilon .\end{align} \tag{ 2 }$$

The covariate coefficients associated with GHGs, NATs and OANT forcings, which are referred to as scaling factors, are defined as follows:

$$\begin{equation}\left\{ {\begin{array}{*{20}{l}} {{\beta _{{\text{GHG}}}} = {\beta _1} + {\beta _3}} \\ {{\beta _{{\text{NAT}}}} = {\beta _1} + {\beta _2}} \\ {{\beta _{{\text{ONAT}}}} = {\beta _1}} \end{array}.} \right.\end{equation} \tag{ 3 }$$

Scaling factors are statistically significant if their associated p values are below 10%. If the scaling factors are statistically significant and positive, this suggests that a positive signal is detected in the observations. A negative scaling factor suggests that no physical interpretation can be made or that some important physical processes and other external forcings are not well represented in the climate model (Bindoff et al 2013). If the scaling factors are close to unity, it suggests that the estimated model responses to ALL, NAT and GHG forcings are in accordance with the observed changes.

After estimating the signal amplitudes, we conduct a residual consistency test (RCT) to check whether the estimated regression residuals fit with the climate model noise represented by a second estimate of the IV covariance matrix. After removing all externally forced signals in the regression process, the regression residuals are expected to be consistent with the statistics. An F test is used to compare the variances of regression residuals with the IVs estimated from the control simulations.

In this study, a dynamic evolution analysis is used to explore the temporally changing features of the response signals and attributable trends across the Arctic. Specifically, the scaling factors and corresponding attribution trends of GHG, NAT and OANT forcings are calculated for different periods, based on a moving window approach. Here, the 100 year data (1915–2014) are divided into five different time slices (50 year, 60 year, 70 year, 80 year and 90 year), and one year is used as a sliding time window to analyze the dynamic changes in the response signals and attributable trends across the Arctic. A simple schematic diagram of the sliding window approach is presented in supplementary figure 4.

Figure 1 shows the observed and simulated annual temperature anomalies over the Arctic. The observed average temperature over the Arctic land areas shows a significant warming (p < 0.001) trend (1.20 °C/century) during the period 1915–2014, though a cooling trend from the 1940s to the 1970s and an accelerated warming trend after the 1970s are observed (figure 1). Generally, the climate model simulations under ALL forcings show similar change pattern with the observations, but overestimated the overall trend with a warming rate of 1.53 °C/century (figure 1). As expected, the ensemble weighted mean of simulated temperature under GHG forcings exhibits a monotonically warming trend at a rate of 3.25 °C/century during 1915–2014. The simulated temperatures under NAT forcing exhibit multi-decadal variations without significant change trends (p > 0.05) for the period 1915–2014, though a warming trend is found for the sub-period after 1980s. Under the OANT forcing, a significant (p < 0.001) cooling trend of −1.87 °C/century is detected during the study period 1915–2014.

Zoom In Zoom Out Reset image size

Spatially, both the observations and simulations with ALL forcing agents consistently show a warming trend from 1915 to 2014 over the Arctic, with the largest warming effects (above 2 °C/century) found in Chukotka, Ellesmere Island and Svalbard (figures 2(a)–(c)). The multimodel weighted mean responses of temperature to GHG forcings exhibit a larger warming trend than that of the observations (figure 2(c)). The GHG forcings alone contributed the most to the warming trend in the Arctic, especially at Bank Island, Severnaya Zemlya, Novaya Zemlya and Svalbard, with a GHG-induced warming trend above 4 °C/century (figure 2(c)). NAT forcings exerted no significant effects on Arctic land surface air temperature, except for northern part of Russia with a cooling trend ranging from −0.5 °C to −1 °C/century (figure 2(d)). The multimodel weighted mean responses of temperature to OANT exhibit a significant cooling effect over the whole Arctic, but with a distinct spatial heterogeneity (figure 2(e)).

Zoom In Zoom Out Reset image size

Based on the TLSOF, we quantify the contributions of GHGs, NATs and OANT forcing agents to the observed Arctic land surface air temperature changes. The scaling factors, which scale the simulated temperature responses to best reproduce the observed changes, are shown for the individual and multimodel-weighted mean responses with 90% CIs in figure 3(a). The wider CI means greater uncertainty level, while the narrower one indicates lower uncertainty level. It is found that the multimodel-weighted mean temperature responses to GHGs, NATs and OANT forcings are relatively well aligned with the observations, as indicated by the scaling factors of GHGs (β = 0.84), NATs (β = 0.80) and OANT forcings (β = 0.90) (figure 3(a)). However, large discrepancy is observed between the models, and notably, the CanESM5 model fail to detect the GHG, NAT and OANT signals (figure 3(a)). The RCT indicates that the residual variability that remains in the observations after removing the scaled responses is relatively consistent with the IV in the region as simulated in the climate model control simulations.

Zoom In Zoom Out Reset image size

The warming trends over 1915–2014 that are attributable to GHGs, NATs and OANT forcings are illustrated in figure 3(b). According to the multimodel-weighted mean, the GHG forcings alone have led to a warming of 2.72 °C (CI: 1.42 °C–4.03 °C) from 1915 to 2014, 1.68 °C (CI: 1.28 °C–2.18 °C) of which was offset by the cooling effects induced by OANT forcings. Compared to the findings of Najafi et al (2015), our calculations, which remove the influences of model drift when calculating internal climate variability, result in smaller GHG and OANT forcing effects. In addition, NAT forcings may have also contributed to 0.12 °C (CI: 0 °C–0.25 °C) of warming over the Arctic. By accounting for all forcing agents, a net warming effect (1.16 °C) close to the observed warming of 1.12 °C is obtained for the study period. The results derived from the best estimates indicate that the OANT forcing has offset 61.8% of Arctic land warming induced by GHG emissions between 1915 and 2014. To assess the robustness of our findings, we repeat our analysis with two additional observational datasets: DGLT and GISS. As shown in supplementary figure 5, the results obtained by using the DGLT and GISS datasets are very consistent with those obtained from the CRU dataset.

The flowchart of the detection and attribution of Arctic land surface air temperature changes at the spatial scale is presented in supplementary figure 6. We also examine the spatial patterns of IV based on the average of 48 individual models to reduce the influence of IV (supplementary figure 7). On the basis of the multimodel-weighted mean, the spatial patterns of the scaling factor and attributable trend are illustrated in figures 4 and 5, respectively. Figure 4 shows that the effects of GHGs, NATs and OANT forcings can be detected in most Arctic land areas. Overall, the multimodel mean temperature responses to GHG forcings are robustly detected (β close to 1) over the Arctic, except for Baffin Island and southern Greenland (figure 4(b)). In contrast, the NAT and OANT responses are overestimated in most regions of Arctic North America (e.g. Alaska, Northwest Territories and Nunavut), but underestimated in the Krasnoyarsk Krai of Russia (figures 4(e) and (h)). Above results indicate that the estimates of the responses to external forcing show large spatial heterogeneity, which is possibly related to the spatial difference in IV over the Arctic (supplementary figure 7)

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 5 shows the spatial distributions of the warming trends that are attributable to GHGs, NATs and OANT forcing agents. According to the best estimates, a large spatial difference is found regarding the GHG-induced warming trend across the Arctic. The results suggest that GHGs alone have warmed the Arctic by 1 °C–6 °C across the Arctic over the past century, especially in the regions of Ellesmere Island, Severnaya Zemlya, Svalbard and Krasnoyarsk Krai (figure 5(b)). The trend attributable to NAT forcings is close to zero over the whole Arctic, which is consistent with the weak long-term trends of solar and volcanic forcings over this period (figure 5(e)). The cooling effects of OANT forcings reach 1 °C–2 °C in most Arctic land areas except for Severnaya Zemlya (figure 5(h)). Similar spatial patterns are revealed for the temperature responses to GHGs, NATs and OANT forcings based on the DGLT and GISS observations (supplementary figure 8–11), suggesting the robustness of our findings.

To explore the temporal evolution dynamics of the response signals and attributable trends over the Arctic, the scaling factors and attributable trends are estimated by using moving windows with five different lengths (figure 6). The results indicate that the adoption of a longer time slice leads to smaller IVs (supplementary figure 12). In this section, the average IVs of 48 climate models are used to adjust the observed and simulated temperatures to eliminate the influence of IV on the attribution results. Like previous studies using a large ensemble of climate models to estimate the IV (Najafi et al 2015, Mueller et al 2018), we use all 48 models to account for the uncertainty of IV estimation instead of the nine models. After conducting a variance analysis, nonsignificant (p > 0.05) differences are observed among the IVs within different time slices, with IVs ranging from 0.0026 to 0.0037 °C/year (supplementary figure 13).

Zoom In Zoom Out Reset image size

Overall, smaller fluctuations in the GHG, NAT and OANT scaling factors are observed for longer time windows (figures 6(a)–(e)). The multimodel GHG, NAT and OANT responses are commonly underestimated (β> 1) between 1915 and the 1960s (figure 6(a)). For the period between the 1920s and 2014, the GHG and OANT signals are overestimated (β < 1) at all slice sizes (figures 6(a)–(e)). In addition, NAT scaling factors from above 1 to below 1 are observed for all window sizes except for 90 years (figures 6(a)–(e)). The strengths of the GHG and OANT signals have increased since the 1930s, with a scaling factor gradually approaching one (figures 6(a)–(e)).

The trends attributed to the GHGs, NATs and OANT forcings when using moving windows with different lengths are presented in figures 6(f)–(j). Due to the large uncertainties of the scaling factors, the best estimates of the attributable trends are used for analysis. For GHG-induced warming, the attributable trends significantly (p < 0.001) increase over time for all window sizes (figures 6(f)–(j)). The effect of OANT forcings is sensitive to the length of the utilized time slice, suggesting great differences in the amounts of aerosol emissions in different periods which could highly impact the relative contribution of each forcing component (figures 6(f)–(j)). The geographical distribution of anthropogenic aerosol emissions exhibited a shift from the western hemisphere to the eastern hemisphere around late 1970s. This evolution of aerosol forcing distribution is expected to have an impact on the Arctic air temperature trend and its spatial pattern. Notably, the sign of the NAT forcing effect has changed from cooling to warming, probably due to the reduction in volcanic activity toward the end of the study period (figures 6(f)–(j)). The revealed temporal evolution patterns based on the DGLT and GISS datasets are basically consistent with those obtained from the CRU dataset (supplementary figures 14 and 15).