Changes in Winter Temperature Extremes From Future Arctic Sea‐Ice Loss and Ocean Warming

Observed rapid Arctic warming and sea‐ice loss are likely to continue in the future, unless and after greenhouse gas emissions are reduced to net‐zero. Here, we examine the possible effects of future sea‐ice loss at 2°C global warming above pre‐industrial levels on winter temperature extremes across the Northern Hemisphere, using coordinated experiments from the Polar Amplification Model Intercomparison Project. 1‐in‐20‐year cold extremes are simulated to become less severe at high‐ and mid‐latitudes in response to Arctic sea‐ice loss. 1‐in‐20‐year winter warm extremes become warmer at northern high latitudes due to sea‐ice loss, but warm by less than cold extremes. We compare the response to sea‐ice loss to that from global sea surface temperature (SST) change also at 2°C global warming. SST change causes less severe cold extremes and more severe warm extremes globally. Except northern high latitudes, the response to SST change is of larger magnitude than that to Arctic sea‐ice loss.

• Less severe winter cold extremes in northern mid-and high-latitudes in response to future Arctic sea-ice loss • Winter hot extremes increase in severity over high latitudes due to future Arctic sea-ice loss, but warm less than cold extremes • In a majority of the latitudes, both cold and hot extremes warm more in response to future global sea surface temperature change than due to sea-ice loss

Supporting Information:
Supporting Information may be found in the online version of this article. species, and human livelihoods, and they are projected to become more frequent in future (Serreze et al., 2021). Work is needed to investigate changes in winter warm extremes due to future Arctic sea-ice loss.
Furthermore, there is uncertainty about the influence of Arctic amplification on atmospheric circulation and mid-latitude severe weather (Cohen et al., 2020;Overland et al., 2021). For example, coupled atmosphere-ocean models suggest that Arctic sea-ice loss intensifies the wintertime Siberian High, but the temperature response is not robustly simulated (Labe et al., 2020;Screen & Blackport, 2019;Screen et al., 2018). Uncertainty comes from the different climate models, different forcings and methodologies, and in some cases, relatively small ensembles used (Cohen et al., 2020;Overland et al., 2016;Smith et al., 2019). This provides a strong rationale for using coordinated experiments in a large multi-model ensemble.
The Polar Amplification Model Intercomparison Project (PAMIP) provides a set of coordinated experiments designed to understand the causes as well as the consequences of polar amplification (Smith et al., 2019). It is a contribution to the Coupled Model Intercomparison Project Phase 6 (Eyring et al., 2016). By running standardized experiments in different climate models and generating large ensembles from each model, PAMIP helps to provide a better estimate of the forced response and to quantify model uncertainty (Screen et al., 2018). PAMIP simulations have been used to study, for instance, the effects of Arctic sea-ice loss and/or warming on the North Pacific jet stream (Ronalds et al., 2020), poleward heat transport (Audette et al., 2021), the wintertime Siberian High (Labe et al., 2020), and mid-latitude westerly winds (Smith et al., 2022).
Here, we utilize the atmosphere-only PAMIP experiments for the first time to assess the respective responses of boreal winter cold and warm extremes to future Arctic sea-ice loss and sea surface temperature (SST) change associated with 2°C global mean warming above pre-industrial levels. We focus on land regions in the Northern Hemisphere, where extreme temperatures have direct impacts on their communities. Using daily temperature output from 10 PAMIP models, each of which having up to 200 ensemble members, we examine the changes in 1-in-20-year cold and warm events. Expanding on previous studies, we study both cold and hot extremes and also examine the respective responses to sea-ice loss and SST change. Consideration of the response to SST change is important, as the local cooling in response to sea-ice loss proposed in earlier studies may be overwhelmed by warming due to global SST change.

PAMIP Experiments
We compare model-simulated temperatures between three PAMIP atmosphere-only time slice experiments. First, we use an experiment forced by present-day (i.e., 1979-2008 climatological) SSTs and sea-ice concentration (Smith et al., 2019), denoted as "pd" hereafter. Second, we use an experiment forced by present-day SSTs but future Arctic sea-ice concentration representative of 2°C global average warming above pre-industrial levels. This experiment is denoted as "futArcSIC." Third, we make use of an experiment in which climate models are forced by future SSTs representative of 2°C global warming but sea-ice concentration at the present-day level. This experiment is referred to as "futSST." We note that 2°C global average warming above pre-industrial levels is equivalent to 15.7°C in absolute global mean temperature, and that sea ice thickness changes are not included this these experiments (Smith et al., 2019). All of these experiments are 1-year time slices with radiative forcing from the year 2000. As such, comparing futArcSIC with pd provides an estimate of extreme temperature changes due to future Arctic sea-ice concentration loss, whereas comparing futSST with pd provides an estimate of changes due to future SST change.
These experiments are run by climate models with a minimum of ∼100 winters to generate large ensembles that are suitable for studying climate extremes (Smith et al., 2019). We make use of daily minimum (tasmin) and maximum (tasmax) near-surface air temperature outputs from 10 climate models, as listed in Table S1 in Supporting Information S1. Specifically, we focus on the respective changes in minimum tasmin and maximum tasmax in boreal winter (December-January-February or DJF) due to future Arctic sea-ice loss and SST change. All included models have daily tasmin and tasmax outputs for pd and futArcSIC. A subset of six models also have outputs for futSST. More than half of the models have at least 200 ensemble members. We use a maximum of 200 members from each model to compute the differences in 1-in-20-year winter minimum and maximum temperatures at each model grid point due to Arctic sea-ice loss and SST change. A 1-in-20-year event has a 5% chance of occurring in any given year, and we use it to represent extremes. A maximum of 200 members is a large enough sample size for this return period, but more members could have been used from some models ( Table S1 in Supporting Information S1) (Thompson et al., 2017). By focusing on DJF minimum and maximum temperatures, we avoid averaging out the extremes in seasonal means (Francis, 2021). We conduct an additional return period analysis at the regional scale (Sections 2.2 and 2.3).
The included models have different atmospheric horizontal resolutions, ranging from 0.83° × 0.56° in HadGEM3-GC31-MM (Andrews et al., 2020) to ∼2.8° in CanESM5 (Swart et al., 2019). For all grid cells in the Northern Hemisphere, we calculate the difference in 1-in-20-year minimum and maximum temperatures between the PAMIP experiments in individual models, as well as the multi-model mean difference (giving each model equal weight). When considering the individual models, we compute the temperature difference in the models' native grids. When considering the multi-model mean, we regrid (through nearest-neighbor regridding) all model results to CanESM5's grid because it is the coarsest among the studied models, before computing the multi-model mean difference.

Regions
We perform analyses in 14 selected regions in the northern mid to high latitudes. These regions are selected from a pre-defined set of regions that are ∼2 Mm 2 in size and designed for examining climate extremes and their impacts (Stone, 2019). The regions are shown in Figure S1 in Supporting Information S1.

Return Period Analyses
We compute return periods by sorting each temperature ensemble of DJF minimum tasmin (and maximum tasmax) in ascending (descending) order and dividing the size of the ensemble by the ranks of the temperature values within the sorted series. We find the difference in 1-in-20-year temperature between experiments at each model grid point. We test whether the two samples of temperatures (i.e., not just the 1-in-20-year values) from different experiments are significantly different using a Kolmogorov-Smirnov (KS) test (Daniel, 1990).
For the regional analysis, we produce and compare return period curves from the pd simulations and the ERA5 reanalysis (Hersbach et al., 2020). We find the regional mean DJF minimum tasmin and maximum tasmax by area-weighted averaging values across native grid cells whose grid point values are within the boundary of each region. Since the present-day conditions in pd are based on 1979-2008 climatology, we extract ERA5 data from the same time period for comparison. This comparison is not completely like-for-like because inter-annual variability exists in ERA5 but not in pd, which has constant boundary forcing. To remove the climate change signal from the regional ERA5 time series and approximately isolate internal variability, we fit a linear trend to the corresponding DJF mean tasmin (and tasmax) time series and remove this trend from the 1979-2008 DJF minimum tasmin (maximum tasmax) time series. This ensures that the trends in the winter season, not just in the extremes, are removed. We then add the regional 1995-2005 average DJF minimum tasmin (maximum tasmax) value in ERA5 to the detrended data, to obtain absolute temperatures for comparison with model output. We choose the 1995-2005 decade because it is centered on year 2000, the year from which radiative forcing is used in the PAMIP time slice experiments.
The modeled pd data do not need detrending because they come from large ensembles of time slice simulations and use a constant radiative forcing. To bias-correct data from each model, we remove from each ensemble member the bias between ensemble-mean regional-mean DJF minimum tasmin (maximum tasmax) and the corresponding 1979-2008 mean regional-mean ERA5 value. We then find the return period curves based on bias-corrected pd data and detrended ERA5 data.
We estimate the uncertainty associated with the ERA5 return period curve by resampling with replacement the ERA5 distribution 1,000 times, though acknowledging that uncertainty sampling in the extremes is limited by the observations. The comparison between individual model return period curves and the ERA5 90% confidence interval enables us to identify models that simulate present-day winter temperature extremes reasonably well in the selected regions. Figures 1a and 1c show this comparison for the North European Economic Area (EEA), for which four and two models (indicated by dashed lines) are excluded in model selection for cold and warm extremes, respectively, because their return period curves are outside the ERA5 envelope at a majority of return periods.
To assess the effects of future Arctic sea-ice concentration loss and SST change on regional winter extremes, we find the return period curves using the futArcSIC and futSST ensembles, respectively. Example return period curves from futArcSIC, futSST, and pd simulated by IPSL-CM6A-LR for the North EEA region are shown Figures 1b and 1d. For each model and region, we find the temperature difference between futArcSIC and pd, and between futSST and pd, at the 20-year return period. For analyses involving futArcSIC, we report the temperature  differences from the individual models, as well as the multi-model mean across all 10 models and the mean across a subset of models that simulate the present day well (according to ERA5). This subset varies from region to region and between cold and warm extremes ( Figure S2 in Supporting Information S1). For analyses involving futSST, we mainly report the multi-model mean temperature difference across the six models for which there is output for this experiment (Table S1 in Supporting Information S1) for brevity. Figure 2a shows the multi-model-mean difference in 1-in-20-year winter cold extremes between futArcSIC and pd in the Northern Hemisphere. The largest warming, of over ∼2.5°C, is projected for northern and eastern Canada near Hudson Bay. The futArcSIC and pd winter minimum temperature distributions are statistically significantly different at the 5% level, indicating amplified warming in boreal winter cold extremes due to future Arctic sea-ice loss, as global average temperature is 1.4°C higher in futArcSIC than in pd (Smith et al., 2019). A statistically significant warming of ∼2°C is also projected for Alaska. These results are generally consistent across the models ( Figure S3 in Supporting Information S1), likely due the close proximity to imposed sea ice reductions in Hudson Bay, Labrador Sea, and Bering-Chukchi Seas (Smith et al., 2022).

Responses to Sea-Ice Loss
In the multi-model mean, ∼1°C warming is simulated in Greenland, across Scandinavia and in northern Russia. However, there is inconsistency in the sign between the models, with MIROC6 and TaiESM1 simulating some cooling in central Greenland, CanESM5, and CESM2 simulating cooling over Scandinavia, and four models simulating cooling in different parts of north Russia ( Figure S3 in Supporting Information S1). At the mid and low latitudes, cooling responses are seen for the United States, parts of Europe and central and eastern Asia. In some models, this cooling is up to about −1°C, suggesting intensified winter cold extremes. However, this response is not statistically significant and is less robust in terms of spatial extent and magnitude than the aforementioned higher-latitude warming response ( Figure S3 in Supporting Information S1). Figure 2b shows the multi-model-mean difference in 1-in-20-year winter warm extremes between futArcSIC and pd. Statistically significant changes are only simulated in the high latitudes, with northern Canada showing the strongest warming, of over ∼2.5°C, followed by northeastern Russia (∼2°C). These changes are generally consistent across the models ( Figure S4 in Supporting Information S1). The multi-model mean indicates widespread cooling of up to −0.4°C that is not statistically significant across most parts of North America, Eurasia, and central Africa. Individual models simulate a stronger cooling response in different parts of the continents, although the responses are not statistically significant ( Figure S4 in Supporting Information S1). A greater warming of cold extremes (Figure 2a and Figure S3 in Supporting Information S1) compared to warm extremes (Figure 2b and Figure S4 in Supporting Information S1) implies a reduced temperature variance.
Next, we examine the regionally averaged differences in 1-in-20-year winter cold and warm extremes between futArcSIC and pd in 14 selected regions over the mid-to-high northern latitudes, where the largest and most significant temperature responses are simulated. Figure 3 shows the results from the individual models (circles), as well as the multi-model-mean responses across the 10 models (yellow crosses) and the multi-model-mean responses across selected models (i.e., those that simulate regional present-day climates that are consistent with the ERA5 reanalysis; black squares).
Like in Figure 2a, the regional analysis for cold extremes (Figure 3a) reveals the largest average warming response in East Canada, with the models simulating regional-mean warming between 2 and 6°C. In the multi-model mean, all selected regions are projected to experience a warming of winter cold extremes due to Arctic sea-ice loss, with values ranging from 0.4°C (inter-model range: −0.9°C to 1.4°C) in West Siberia to 4.1°C (range: 2.5°C to 5.4°C) in East Canada. The mean results are similar across the subsets of models (note, no model is consistent with ERA5 in East and West Canada even after mean bias correction). Despite the general warming response seen in the multi-model mean, some models simulate intensified winter cold extremes in regions including West Canada, North EEA, Northwest and Southwest Russia, and West and Northeast Siberia. However, these cooling responses are not statistically significant ( Figure S3 in Supporting Information S1) and could be due to internal variability.
The regional results for winter warm extremes are shown in Figure 3b. Nunavut in northern Canada has the largest multi-model-mean warming response to future Arctic sea-ice loss, with individual models simulating 2°C-4°C warming. As shown in Figure 2b, North Pacific Russia has the second largest mean response at 1.5°C (range: 0.2°C to 2.9°C). For the rest of the regions, the multi-model-mean response is within about ±1°C, ranging from −0.2°C (range: −0.7°C to 0.4°C) in North EEA to 1.1°C (range: −0.2°C to 2.2°C) in Sakha. Eleven of the 14 selected regions (except Nunavut, East Canada and North Pacific Russia) have at least one model simulating a cooling response, showing a smaller signal-to-noise ratio than cold extremes. Overall, Figure 3 shows that the pd model simulations do not compare very well with reanalysis even after bias correction, partially because of the idealized nature of the experiments. However, this does not affect our main results. Figure 4a shows that warmer SSTs associated with 2°C global mean warming increase 1-in-20-year cold temperatures over land in the Northern Hemisphere in the multi-model mean. This warming is statistically significant at the 5% level. No cooling response is shown in the multi-model mean at any location. In general, individual models agree on a strong (∼3°C) warming signal in North America, particularly in the western parts ( Figure S5 in Supporting Information S1). The cold temperature response in Eurasia to future SST change is more variable, with IPSL-CM6A-LR showing strong warming in the northern parts, whereas FGOALS-f3-L shows cooling in those parts but relatively strong warming in east Asia ( Figure S5 in Supporting Information S1). These differences may be affected by sampling variability.

Responses to SST Change
Future SST change is also projected to warm winter warm extremes significantly in all Northern Hemisphere land grid cells in the multi-model mean, as shown in Figure 4b. However, both the multi-model mean and individual model results ( Figure S6 in Supporting Information S1) indicate that the warm extreme response is smaller compared to the cold extreme response in almost all places except northern Canada. Small inter-model differences are seen in the warm extreme response to SST change, with CanESM5 simulating cooling in Greenland and northeastern Russia that is not statistically significant, for example.
With previous evidence that responses to sea ice and greenhouse gas forcing are approximately linearly additive (McCusker et al., 2017), it may be reasonable to deduce the combined mean 1-in-20-year winter temperature responses to Arctic sea-ice loss and ocean warming from Figures 3 and 4. For cold extremes, even in places where Arctic sea-ice loss is simulated to intensify them (e.g., in southwestern United States, parts of Europe, central and eastern Asia, though not statistically significantly; Figure S3 in Supporting Information S1), warming due to SST change overwhelms this cooling effect, resulting in net warming (not shown).
Indeed, by comparing the multi-model mean of the 1-in-20-year cold temperature differences due to Arctic sea-ice change (x-axis) and SST change (y-axis) over the 14 selected regions in Figure 4c (navy markers), we find that the warming response to future SST change is larger than or equal to the response to future Arctic sea-ice loss in 11 regions (i.e., except Nunavut, East Canada, and North Pacific Russia). The three exceptions suggest that the response to sea-ice loss is by far the largest near the regions of sea-ice loss, whereas warming due to SST change is more spatially homogeneous. The ratio of SST change-induced response to sea ice loss-induced response ranges from 0.5 in East Canada to 7.5 in Southwest Russia. Since all selected regions are projected to experience multi-model-mean warming to both sea-ice and SST changes, an enhanced combined response is expected. For East Canada, this may mean a combined response of 5.8°C.
For warm extremes, warming from SST change also dominates over the small and non-statistically significant Arctic sea ice-loss induced cooling responses in North America, Eurasia, and Africa, resulting in net warming. Figure 4c shows this clearly, where all but one orange markers (i.e., except for Nunavut) are above the 1:1 identity line. The ratio of the magnitude of SST-induced response to the magnitude of sea ice-induced response ranges from 0.8 in Nunavut to 35 in North EEA (because of a near-zero response to sea ice). In Nunavut (northern Canada), where winter warm extreme is projected to become statistically significantly warmer due to Arctic sea-ice loss and SST change separately, the combined effect may mean intensification of warm extremes by 5.4°C, although we emphasize that our results are based on idealized atmosphere-only experiments.

Discussion and Conclusions
Arctic amplification has been a topic of interest in the literature, not only because it is one of the strongest anthropogenic climate change signals, but also because of its wide-reaching effects on the climate system (Labe et al., 2020;Screen et al., 2013). This study is the first to use targeted and coordinated PAMIP experiments to examine the response of rare (1-in-20-year) winter temperature extremes to Arctic sea-ice loss and SST change at 2°C global mean warming above pre-industrial levels. It is also the first to investigate winter warm extremes, and to bias-correct the PAMIP simulations and apply model selection based on reanalysis data.
We have shown a multi-model-mean warming response of winter cold extremes to future Arctic sea-ice loss across the mid and high latitudes. This is consistent with the projected decrease in the likelihood and severity of mid-and high-latitude cold extremes in previous studies (Ayarzagüena & Screen, 2016;Screen et al., 2015b). For 8 of the 14 selected regions (excluding West Canada, North EEA, Northwest and Southwest Russia, and West and Northeast Siberia), the sign of change is robust across 10 models. Where a local cooling response is simulated in some models, the location of this cooling is not robust across models, and may be a sign of internal variability. We Comparison between the multi-model mean temperature changes due to future Arctic sea-ice loss (x-axis) and the corresponding changes due to future SST change (y-axis). Navy points indicate changes in 1-in-20year DJF minimum of daily minimum temperature, whereas orange points indicate changes in 1-in-20-year DJF maximum of daily maximum temperature. Each point represents the regional mean in one particular region. The dashed line indicates a 1:1 relationship, whereas the dotted line indicates a 2:1 relationship.
cannot rule out a weak cooling response, as suggested by previous studies (Labe et al., 2020;Zappa et al., 2021), but it appears to be model dependent.
The winter cold extreme response to future SST change is more robust, with almost all of the Northern Hemisphere showing a warming response in all available models. Notably, this warming response exceeds the sea ice-induced cooling response in southwestern United States, Western Europe, and central and eastern Asia. Overall, our results imply that some of the adverse impacts of cold extremes on, for instance, human health (Mäkinen, 2007;Vasconcelos et al., 2013) and transport and power supply (Screen et al., 2015b) are expected to be lessened in the mid and high latitudes in the future. However, we stress that Arctic warming and sea-ice loss are already impacting the Arctic communities (Moerlein & Carothers, 2012), whose lifestyles and livelihoods were adapted to cold weather through generations of lived knowledge.
For winter warm extremes, we have shown that statistically significant responses to future Arctic sea-ice loss are limited to the high latitudes, primarily to northern Canada and northeastern Russia. Non-significant responses are found for the rest of the hemisphere, and overall the warm extreme response is weaker than the cold response. This suggests a reduced winter temperature variance due to Arctic sea-ice loss, which is consistent with the literature (Blackport et al., 2021;Collow et al., 2019;Screen, 2014). SST-induced warming is larger than the sea ice-induced changes in most places.
Warming of winter warm extremes in the high latitudes due to Arctic sea-ice loss and ocean warming can increase the chances of rain on snow events. Notable events have already occurred in Arctic Canada (Rennert et al., 2009) and Russia (Forbes et al., 2016), which led to declines in ungulate (e.g., reindeer and musk oxen) populations that persisted for years and herders losing food security and transportation (Serreze et al., 2021). Our results suggest that these communities are at an increased risk of these impacts in a 2°C warmer world, compared to the present day.
Model differences have been shown in this study. Future work is needed to understand these differences to reduce uncertainty. Furthermore, sea-ice loss does not happen in isolation, but considering it together with future ocean warming is not routinely done. Going forward, we recommend researchers place a stronger focus on the SST component or the net response. Moreover, the combined effect of Arctic sea-ice loss and SST change on winter temperature extremes has not been studied here. Potential non-linearities in their effects may mean that a combined future sea-ice and SST experiment in PAMIP is important. Future work should also quantify the resulting impacts on various aspects of society through coupled climate and impact modeling.
Aside from sea-ice concentration loss and SST change, PAMIP provides a range of experiments designed to investigate the impacts of sea-ice thickness changes and full ocean dynamics (Smith et al., 2019), which have not been studied here. Our estimates of the responses to sea-ice loss may be conservative because both ice thickness changes (Labe et al., 2018) and atmosphere-ocean coupling (Deser et al., 2015(Deser et al., , 2016Smith et al., 2017) have been suggested to strengthen the response. It is recommended that researchers fully exploit the PAMIP data to investigate the effects of these changes.

Conflict of Interest
The authors declare no conflicts of interest relevant to this study.

Data Availability Statement
The PAMIP data used in this study are available at Earth System Grid Federation (ESGF) via https://esgf-node. llnl.gov/search/cmip6/. A user guide for creating an ESGF account and downloading the data can be found at https://esgf.github.io/esgf-user-support/. PAMIP data information from each modeling center, including their contact information, can be found at https://www.cesm.ucar.edu/projects/CMIP6/PAMIP/. The ERA5 reanalysis data are available from the Copernicus Climate Data Store https://cds.climate.copernicus.eu/cdsapp#!/dataset/ reanalysis-era5-single-levels?tab=overview.