Large amounts of snowfall in winter over northern continents have important implications for human society and natural ecosystems, i.e., causing disruptions to transport, energy supply, and power transmission, as well as affecting runoff and water supply for coming seasons. Winter snowfall is primarily determined by air temperature and availability of atmospheric moisture. It is generally accepted that winter snowfall would be expected to decrease over northern continents due to greenhouse warming (Solomon et al., 2007). Observations and climate model simulations show that winter has been warming greater than any other seasons (not only with respect to the mean temperature, but also the daily minimum and maximum temperature), particularly in northern high-latitudes (Rangwala et al., 2013; Serreze and Barry, 2011; Knutti and Sedláček, 2013). This acts to reduce the fraction of precipitation that falls as snow and increase snowmelt. Using historical simulations and RCP4.5 scenario projections of CMIP5, Krasting et al. (2013) pointed out that projected snowfall changes may be apparent over most regions in the Northern Hemisphere in the 21st century. O’Gorman (2014) contrasted the response of mean snowfall and snow extreme to the climate change by comparing the model simulated snowfall change in the warm stage of 2081–2010 in the RCP8.5 scenario and control stage of 1981–2000 in the historical simulations. In the quite cold regions, the fraction of both mean snowfall and snow extreme is increasing in a warmer environment, while in regions that are not so cold, the fraction of mean snowfall decreases dramatically.

In contrast, for recent winters, anomalously large snowfall has been experienced in large parts of North America, Europe, Siberia, and East Asia. An increasing body of recent research indicates significant changes in the Arctic environment, and that an accelerated decline of Arctic sea ice might be associated with the recent snowy winters (Liu et al., 2012; Francis and Vavrus, 2012; Ghatak et al., 2012; Cohen et al., 2010). Recent studies argued that the dramatic Arctic sea ice loss adds tremendous extra heat to the atmosphere, weakening the north-south temperature gradient and hence weakening westerly winds. The weakened westerly winds tend to enhance blocking circulations, favoring more frequent incursions of cold air mass from the Arctic into northern continents (Liu et al., 2012; Francis and Vavrus, 2012; Mori et al., 2014; Kug et al., 2015). In addition, the dramatic Arctic sea ice loss greatly increases the moisture flux from the ocean to the atmosphere, providing an enhanced local moisture source for some regions surrounding the Arctic Ocean (Liu et al., 2012). Together, these changes could lead to increased winter snowfall.

Here we present results that suggest such an emerging climate change will be a persistent feature in a warming climate from a new phase of the internationally coordinated Coupled Model Intercomparison Project (CMIP5, Taylor et al., 2012). Some of the key improvements over the previous phase (CMIP3) include a more diverse set of model types, a number of improvements in physics, numerical algorithms, and configurations, and a new set of scenarios called Representative Concentration Pathways (RCPs). Sixteen CMIP5 models on the Program for Climate Model Diagnosis and Intercomparison (PCMDI) data portal that provide pre-industrial simulations, “all-forcings” historical simulations and projection simulations under the RCP4.5 scenario are analyzed. The RCP 4.5 is a medium-mitigation emission scenario that tends to stabilize direct radiative forcing at ~4.5 W/m² (~650 ppm CO₂ equivalent) near the end of the 21st century (Moss et al., 2010). The RCP4.5 is more conservative than the Special Report on Emission Scenarios (SRES) A1B (end-of-century CO₂ of 720 ppm).

In this study, we focus on the snowfall flux (precipitation of all forms of water in the solid phase), which reflects both snow cover and depth, and is more closely tied to climate forcings than snow cover. For instance, in winter, snow cover may have no variability when terrestrial Arctic is all snow covered, whereas snowfall flux may still have considerable variability. Figure 1a shows changes in the winter snowfall flux (averaged for December-January-February) in the extra-tropical Northern Hemisphere during 2000–2019 relative to 1980–1999 for the multi-model ensemble mean of the 16 CMIP5 models. Compared to the last two decades of the 20th century, the winter snowfall during the first two decades of the 21st century is clearly predicted to increase in the northern Alaska, central and eastern Canada, Greenland, northeastern Europe, Siberia, and northern and western China, and decrease in most regions of the United States, Mexico, western and central Europe and southeastern China. Moreover, each 20-year segment change in the winter snowfall flux extending into the 21st century (Figs 1b-d) exhibits a steady increase of the winter snowfall in these areas, and both the pattern and magnitude of the winter snowfall change in Figs 1b-d resemble those in Fig. 1a. We also note that the magnitude of the predicted winter snowfall increase in these areas is obviously reduced during the last two decades of the 21st century relative to 2060–2079, which coincides with the timing of the approximate stabilization of direct radiative forcings under the RCP4.5 (Moss et al., 2010).

A linear least-squares fit regression is then applied to the winter snowfall flux anomalies for the multi-model ensemble mean in the extra-tropical Northern Hemisphere for the period 2000–2099 to capture the trend. Consistent with every 20-year segment change, significant positive trends in the winter snowfall are found in the northeastern Europe and central and northern Asia and northern North America (>99% significance, Fig. 2a). This pattern of winter snowfall trend resembles that in Fig. 2c of Krasting et al. (2013). The winter snowfall averaged over northern Eurasia and northern North America (the positive territories within thick black contours in Fig. 2a), respectively, increases about 24.0% and 21.8% by the end of the 21st century under the RCP4.5. We also noted that the increase trend of the averaged snowfall over northern Eurasia is flattering from ~2060 to the end of 21st century, which is corresponding to the stabilization of direct radiative forcings in RCP4.5 since 2060 (Fig. 2b). The increase is also gets flattering for the averaged snowfall over the northern North America from ~2050 (ten years earlier than in the northern Eurasia), as shown in Fig. 2c. To confirm the robustness of the trend pattern derived from the multi-model ensemble mean, we repeat the trend analysis for the winter snowfall flux for each individual CMIP5 model. Encouragingly, all the 16 CMIP5 models show spatial distributions of the winter snowfall trend that are very similar to that of the multi-model ensemble mean (not shown here) although the regional details differ somewhat among different models (i.e., the Midwestern United States). Using the previous phase (CMIP3) model simulations, one study also showed an increase of snow water equivalent in extended winter over the coldest areas of the Northern Hemisphere continents under the SRES A1B scenario (Räisänen, 2008).

To identify the relationship between the continental snowfall in boreal winter and the Arctic sea ice in autumn during the 21st century, the Arctic sea ice extent are gained by the simulated sea ice concentration (100% is regarded when the value is greater than 15%) multiplied by corresponding grid area. Figure 3a shows time series of autumn Arctic sea ice extent from the multi-model ensemble mean for the period 2000–2099 (thick black line). The projected autumn sea ice extent exhibits a persistent and linear decreasing trend until the late 2070s, and then tends to level off, resulting in about a 60.8% reduction in sea ice cover by the end of the 21st century. To what extent the increase of winter snowfall in the northern high-latitude continents during the 21st century (Fig. 2) can be attributed to the diminishing autumn Arctic sea ice in a warming climate is a key question. To address this, we calculate the fraction of the winter snowfall flux trend that is linearly congruent with the autumn Arctic sea ice extent by (1) regressing the winter snowfall flux anomalies onto the index of the standardized autumn Arctic sea ice extent anomalies, and (2) multiplying the resulting regression coefficients by the trend in the index of the standardized autumn Arctic sea ice extent. To be consistent with Fig. 2, the time series we used in Fig. 3 is still from 2000 to 2099. As shown in Fig. 3b, the decline of autumn Arctic sea ice does favor an increase of the winter snowfall in the northeastern Europe, central and northern Asia and northern North America (>99% significance). The signature of the autumn Arctic sea ice strongly reflects both the pattern and amplitude of the trends in the snowfall flux during winter season (Fig. 3b vs Fig. 2b). Roughly all of the winter snowfall flux increases during the 21st century averaged over the northern Eurasia and North America, respectively, are linearly congruent with the trend in the autumn Arctic sea ice. To confirm the robustness of the regression pattern derived from the multi-model ensemble mean, we repeat the regression analysis for each individual model. Interestingly, all the 16 CMIP5 models show the winter snowfall changes that are associated with the reduction of sea ice are in very good agreement with that of the multi-model ensemble mean.

To further test the strong sea ice-snowfall relationship, a singular value decomposition analysis (SVD, a multivariate statistical analysis, used to capture a more compact representation between the multivariate datasets by reducing the large number values to significant fewer values, and provide insight exhibitions into spatial and temporal variations, thus, it is a method for detecting coupled variability between different components of climate system; Bretherton et al., 1992), is performed on the autumn Arctic sea ice concentration anomalies and the winter snowfall flux anomalies north of 20°N for the period 2000–2099 to identify the dominant coupled pattern. Figure 4 shows spatial patterns and temporal coefficients of the first SVD mode that accounts for 99.1% of the total covariance. The first SVD mode of the autumn of Arctic sea ice has negative values almost everywhere in the Arctic, and a poleward increase of variability (Fig. 4a). The corresponding time coefficients show a substantial upward trend (>99% significance), indicating the reduction of autumn Arctic sea ice during the 21st century. The largest reduction is in the central Arctic Ocean. Associated with the basin-wide sea ice reduction, the first SVD mode of the winter snowfall flux shows a clear feature of increasing snowfall in much of northern high-latitude continents, which resembles the regression pattern (Fig. 4b). The correlation between the time series of the first SVD mode of the autumn Arctic sea ice and winter snowfall flux is 0.98 (>99% significance, 0.76 with the trend removed). The consistency of the regression and SVD analyses attests the central role of the diminishing autumn Arctic sea ice in the increased winter snowfall. Moreover, the sea ice-snowfall relationship appears to be enhanced during the 21st century (Fig. 4c).

Comparisons between the projection and pre-industrial (considered as unforced variability since external forcings are fixed) simulations suggest that the linear decreasing trend in the Arctic sea ice extent during the 21st century is primarily attributed to increased greenhouse gases in the atmosphere, even though the role of natural variability in the recent dramatic decline of Arctic sea ice remains a critical research question (Liu et al., 2013; Notz and Marotzke, 2012; Stroeve et al., 2012). The Arctic Oscillation (AO) is considered as the dominant mode of natural variability over northern high latitudes (Thompson and Wallace, 2001). It has been a strong descriptor of the atmospheric circulation in the extra-tropical Northern Hemisphere (Thompson et al., 2000). In winter, negative phases of the AO favor anomalously cold temperature over the eastern North America and northern Eurasia (Cohen et al., 2012), and increased storminess in the mid-latitudes, which could lead to increased snowfall there. Here we repeat the aforementioned regression analysis for the winter AO for the period 2000–2099. The winter AO index is defined as the leading principle component of the simulated sea level pressure anomalies poleward of 20°N. All the 16 CMIP5 models show the opposing sea level pressure pattern between the Arctic and the northern mid-latitudes, a distinct feature of the winter AO (not shown here). On average (the average of the winter AO calculated by multi-model simulations), the models show the winter AO remains approximately constant throughout the 21st century (Fig. 5a). There is no evidence of a tendency towards negative territory of the AO during the winter season in response to increased greenhouse gases in the atmosphere, even though Arctic sea ice continues to decrease. Moreover, none of the individual models show a significant negative AO response in winter in a warming climate. The lack of AO trends in the future is consistent with other studies using the CMIP5 models (Cattiaux and Cassou, 2013). As shown in Fig. 5b, changes in the winter snowfall flux that are linearly congruent with the winter AO resemble the patterns as in Figs 2a and 3b, but are extremely small as compared to those associated with the autumn Arctic sea ice. This is also supported by the regression analysis for each individual model (Appendix Fig. A1). In short, we find changes in the winter AO contributing much smaller than the Arctic sea ice to the projected increased winter snowfall.

The evidence demonstrated in this study, based on a new generation of more complex climate model simulations, reveals a robust tendency towards increasing snowfall during the winter season over northern high-latitude continents in a warming environment. In fact, the emerging signal of anomalously large snowfall over North America, Europe, Siberia, and East Asia in recent winters exhibits some similarities with the snowfall flux change documented in this study. Our results suggest that the majority of the projected increased winter snowfall in northern high-latitude continents during the 21st century can be attributed to the diminishing Arctic sea ice. The emergence of a strong sea ice-snowfall feedback will be persistent and enhanced in coming decades. This would enhance the likelihood of heavy snowfall events in winter. It also remains to be seen whether other aspects of the simulated response become detectable as Arctic sea ice continues to decline.

The authors acknowledge the climate modeling groups all across the world, the WCRP Working Group on Coupled Modeling, and the Program for Climate Model Diagnosis and Intercomparison for making CMIP5 model outputs available.