Climate change impacts on linkages between atmospheric blocking and North American winter cold spells in CanESM2 and CanESM5

This study investigates changes in linkages between atmospheric blocking and winter (December–February) cold spells over the Pacific-North America region in two large-ensembles of Canadian Earth System Models (CanESM2 and CanESM5 under high-emission scenarios). The two ensembles show decreases in winter blocking frequency over the North Pacific from 1981–2010 baseline to 2071–2100, with larger decreases in CanESM5 (− 3.08%/decade) than CanESM2 (− 1.73%/decade). Using a time-invariant (stationary) threshold estimated from the baseline to define cold days, the two ensembles project a decline in cold spell events as future air temperature increases; the occasional occurrence of cold spell events is still projected to occur at the end of the century. Using a time-dependent (nonstationary) climatological threshold, CanESM2 and CanESM5 ensembles project modest decreases in cold spell days over North America (− 2.0 and − 2.3%/decade). With the nonstationary threshold, the two ensembles project decreases in winter cold spell frequency during blocking, with larger decreases in CanESM5 (13%) than CanESM2 (3%) for 2071–2100 period compared to the baseline. The two ensembles display similar blocking-cold spell linkages between the baseline and future periods; however, the linkage is weaker and exhibits larger uncertainty in the future. Moreover, temperature advection and net heat flux anomalies during blocking are generally weaker for the future period, resulting in weaker impacts on North American cold spells with larger uncertainty associated with increases in internal-variability.


Introduction
An overall decrease in cold extremes has been seen at the global and continental scales in observations collected since 1950; this decrease is mostly attributable to global warming driven by anthropogenic influences (IPCC 2012;Kim et al. 2016;Van Oldenborgh et al. 2019). However, despite the decreasing trends, cold extremes still lead to adverse impacts on human and natural systems (e.g., WMO 2012; Lormée et al. 2013;Roland and Matter 2013;Planchon et al. 2015;Smith and Sheridan 2019). Meanwhile, it is well established over the past decades that atmospheric blocking has been a main driver of persistent winter cold extremes, particularly over Europe (e.g., Hoskins and Sardeshmukh 1987;Demirtaş 2017) and North America (e.g., Quiroz 1984;NOAA-NCEI 2018). Blocking is a quasi-stationary highpressure system occurring at mid-to high-latitudes that interrupts the predominant westerly flow and persists from days up to few weeks (Trigo et al. 2004;Schwierz et al. 2004).
The blocking life cycle (onset, maintenance, and decay), which is generally driven by interactions between background conditions, planetary waves, and synoptic-scale eddies, is an important and ongoing research topic (Woollings et al. 2018;Luo and Zhang 2020;Lupo 2021). Blocking onset is known to be connected with the dynamics of long waves of Rossby type and their nonlinear instability, and breaking is caused by interaction with large-scale motions, synoptic-scale eddies, cyclones, and anticyclones (Luo and Zhang 2020;Lupo 2021). A key feature in blocking onset is a rapid poleward displacement of the subtropical air mass, following up a large-scale extended ridge on the mid-latitude 1 3 jet stream (Woollings et al. 2018). Luo et al. (2019) showed that that cyclonic wave breaking or eddy straining is a concomitant phenomenon of the blocking establishment and maintenance.
Blocking events result in anomalous atmospheric circulation and cold air advection that lead to cold extremes at mid-latitudes (Sillmann et al. 2011;Pfahl and Wernli 2012;Whan et al. 2016). Over the Pacific-North American sector,  reported that a negative-phase of the North Atlantic Oscillation (NAO) that coincides with North Pacific blocking enhances a warm west/cold east dipole over midlatitude North America and also leads to a warm north/cold south dipole anomaly in eastern North America. Che et al. (2021) showed that strong (weak) wintertime sea surface temperature fronts around the Kuroshio extension tend to inhibit (promote) the genesis of North Pacific blocking due to intensified (weakening) background mid-latitude zonal wind; thus, this leads to less (more) cold spells in the midlatitudes of North America.
The Intergovernmental Panel on Climate Change (IPCC 2013) reported that decreases in cold temperature extremes will occur over most land areas through the twenty-first century as global mean temperatures increase, based on projections of a multi-model ensemble of the Coupled Model Intercomparison Project Phase 5 (CMIP5) global climate models (GCMs). Nevertheless, cold extremes are projected to continue to occur in the future, even though their frequency will decline (e.g., Kodra et al. 2011;IPCC 2013;Gao et al. 2015). Therefore, understanding the physical mechanisms underlying changes in cold extremes in future climate is essential for impacts assessments and policy decisions regarding mitigation and adaptation strategies. Accordingly, it is of importance to investigate changes in the behavior of atmospheric blocking and the linkages between blocking and regional cold extremes in a changing climate (e.g., Brunner et al. 2018). Over the Euro-Atlantic sector, Brunner et al. (2018) reported that blocking will continue to play an important role in the development of winter cold spells in the future; however, evolution of blocking-cold spell linkages under a changing climate has not yet been investigated widely over the Pacific-North American sector.
CMIP5 GCMs tend to underestimate winter blocking frequency with small negative biases over the Pacific and Ural regions and significant underestimation over the Euro-Atlantic sector (Davini and D'Andrea 2016;Woollings et al. 2018). These results are an improvement from CMIP Phase 3 (CMIP3) and are likely due to higher horizontal and vertical resolution in the atmosphere-ocean and better representation of the Earth system in CMIP5 GCMs (Anstey et al. 2013;Davini and D'Andrea 2016). The recent CMIP Phase 6 (CMIP6) GCMs also underestimate winter blocking frequency, despite their improved representation of blocking frequency and persistence in the Atlantic and Pacific over CMIP5 models (Davini and D'Andrea 2020;Schiemann et al. 2020). Meanwhile, future projection from the three generations of CMIP GCMs are in agreement that winter blocking frequency will decrease over the Atlantic and Pacific as greenhouse gas concentrations increase (IPCC 2013;Kennedy et al. 2016;Matsueda and Endo 2017;Brunner et al. 2018;Davini and D'Andrea 2020). However, CMIP5 and CMIP6 GCMs generally project weaker decreases in Northern Hemisphere blocking frequency and a more complex response relative to those from CMIP3 (Dunn-Sigouin and Son 2013; Masato et al. 2013;Davini and D'Andrea 2020). Nevertheless, the future intensity and persistence of blocking and the accompanying impacts on regional climate remain unclear (Sillmann and Croci-Maspoli, 2009;IPCC 2013), particularly over the Pacific-North American sector. Jeong et al. (2021) reported that the Canadian Earth system model version 5 (CanESM5) yields an improved representation of winter (December-February) blocking climatology at the Pacific-North American sector relative to its predecessor version 2 (CanESM2), based on overall improvements in simulating large-scale climate patterns, while both models yield similar results in terms of the blocking-cold spell linkages. Nevertheless, the two generations of CanESMs are not yet evaluated for future winter blocking climatology and its impacts on North American winter cold spells. Therefore, here we investigate changes in blocking activity over Pacific-North America and the impacts on North American cold spells in winter throughout the twentyfirst century using two large initial condition ensembles of GCMs: CanESM2 participated in the CMIP5 and CanESM5 contributed to CMIP6. We consider two different ways of defining cold spells (colder than 10th percentiles of daily minimum temperature): stationary (no-adaptation), using a time-invariant cold day threshold calculated from the 1981-2010 baseline, and nonstationary (full-adaptation), using time-variant cold day thresholds calculated over future 30-year periods. Uncertainties due to internal variability in the baseline and future blocking-cold spell linkages are quantified using the large ensembles. We also investigate changes in other synoptic-scale circulation characteristics relevant to blocking-cold spell linkages, including air temperature advection, surface energy fluxes, and 1000-500 hPa geopotential thickness. Note that this is a complement study to the work conducted by Jeong et al. (2022) for the climate change impacts on the linkages between atmospheric blocking and summer heatwaves over the Pacific-North American sector.

Data
CanESM2 and CanESM5 are developed by the Canadian Centre for Climate Modelling and Analysis (CCCma) of Environment and Climate Change Canada (ECCC). They are global scale Earth System Models used to produce seasonal and decadal climate predictions, simulate the historical climate, and make centennial-scale projections of the future climate. CanESM2 has fully-coupled atmosphere, ocean, sea-ice, land, and terrestrial and ocean carbon cycle components (Arora et al. 2011). CanESM5, which is CCCma's latest version of its coupled GCM, incorporates incremental updates to the atmosphere, land-surface, and terrestrial ecosystem models, completely new components for the ocean, sea-ice, and marine ecosystems, as well as a new coupler (Swart et al. 2019).
The CanESM2 ensemble analyzed here consists of 50 initial-condition simulations under historical  and RCP 8.5 (2006-2100) forcing scenarios (Fyfe et al. 2017), while the CanESM5 ensemble also has 50 initialcondition simulations under historical (1981-2014) and SSP5 8.5 (2015-2100) forcing scenarios. Members of each CanESM ensemble share the same external forcings, differing only in initial conditions of the atmosphere and ocean state at the beginning of the historical simulation period. The initialization strategy involves a strategic combination of different oceanic initial states and different atmospheric perturbations; therefore, each member simulates a different realization of internal-variability. RCP 8.5 and SSP5 8.5 are high emission (no mitigation policy) scenarios used for the CMIP5 and CMIP6, respectively, with a radiative forcing of 8.5 W/m 2 . CanESM2 under RCP8.5 and CanESM5 under SSP5 8.5 simulate roughly 5.5 and 8 °C increases in global mean temperature by the end of the twenty-first century compared to the pre-industrial level, respectively.
Daily minimum air temperature at 2-m, geopotential heights at 500, 850, and 1000 hPa, U and V components of wind at 850 hPa, surface short and long wave radiation, and surface sensible and latent heat fluxes for the period 1981-2100 are used. All of these variables are available from the 50-member CanESM2 ensemble and a 20-member subset of the CanESM5 ensemble, which is used in this assessment. Both CanESM2 and CanESM5 employ the spectral transform method with T63 atmospheric resolution in the horizontal (approximately 2.8° × 2.8°), and all variables are interpolated to a grid spacing of 2.5° × 2.5°. The two ensembles are compared in their ability to simulate blocking events and their connection to cold spells over Pacific-North America with observationally-constrained estimates from the National Centers for Environmental Prediction-Department of Energy reanalysis 2 (NCEP-DOE-R2) (Kanamitsu et al. 2002) and ERA-Interim (Dee et al. 2011) over the 1981-2010 baseline period. NCEP-DOE-R2 provides pressure level variables on a 2.5° grid but the other variables on a 1.9° grid; the latter are interpolated onto the 2.5° grid prior to analysis. All climate variables listed above are available for the ERA-Interim at a horizontal grid spacing of 0.75°. As Jeong et al. (2021) has already reported on the historical performance of CanESM2 and CanESM5 in representing the blocking-cold spell linkage in the baseline period, this study mainly focuses on future projections.

Atmospheric blocking
Existing blocking indices can be classified into three broad types: those based on absolute fields (Davini et al. 2012), anomaly fields (Schwierz et al. 2004), or combination of the two (Diao et al. 2006;Barriopedro et al. 2010). Woollings et al. (2018) compared winter blocking frequency diagnosed by the three types for 1958-2012 using three reanalysis datasets. All types display reasonably similar spatial patterns of blocking frequency, particularly over the Pacific-North American sector which shows more frequent blocking over the North Pacific than the North American continent (Fig. 2 of Woollings et al. (2018)).
A blocking event is then identified from the 2-dimensional IB field by taking account the spatial extent and temporal persistence. A spatial IB (SIB) event is identified if contiguous IB extends over at least 15° longitude, while allowing for meridional movement within ± 5° latitude (Davini et al. 2012;Woollings et al. 2018;Brunner et al. 2018). A temporally persistent SIB (TSIB), known as a blocking event, is then defined when a SIB event stays in place, within 10° longitude, for at least five consecutive days (Woollings et al. 2018;Brunner et al. 2018). Blocking events are further diagnosed into five different regions within the 50°N-75°N latitudes from the North Pacific to northern North America to investigate the remote and co-located linkages between blocking location and North American cold spell. The five regions are referred to as B1 (210°W-180°W), B2 (180°W-150°W), B3 (150°W-120°W), B4 (120°W-90°W), and B5 (90°W-60°W). A blocking day in each region is defined when more than half of the blocking region is covered by a blocking event (Brunner et al. 2018;Jeong et al. 2021).
Anomalous patterns of other atmospheric circulation variables, such as temperature advection at 850 hPa, net surface heat flux, and 1000-500 hPa geopotential thickness, are investigated during blocking episodes in winter for the baseline and 2071-2100 future periods. Temperature advection anomalies due to anomalous wind flow are estimated by − T x u � − T y v � (Lehtonen and Karpechko 2016;Sousa et al. 2018), where T represents climatological mean temperature, and u � and v� are anomalous zonal and meridional wind components during blocking relative to climatology for the baseline and future periods. Net surface radiative and turbulent heat fluxes are calculated on a daily timescale via R S + R L − LE − H (Muñoz et al. 2010;Pettenuzzo et al. 2010), where the terms represent downward shortwave and longwave solar radiation and upward surface latent and sensible heat fluxes, respectively. Therefore, daily net surface heat flux roughly represents the summation of the ground heat flux and total heat energy stored in the surface layer for a day. Anomalies of the net surface heat flux with respect to their climatology are calculated during blocking days at the five regions for the baseline and future periods. Four subregions of North America are defined as northwest (NW; 170°W-102°W and 50°N-70°N), northeast (NE; 102°W-52°W and 50°N-70°N), southwest (SW; 140°W-102°W and 30°N-50°N), and southeast (SE; 102°W-52°W and 30°N-50°N) to summarize regional differences.

Cold spell
A cold spell is defined as a period of at least three consecutive cold days with daily minimum temperatures below a threshold. The threshold is calculated at each grid point for each calendar day as the 10th percentile of daily minimum temperatures in a 31-day window centered on the calendar day for a 30-year period (Fischer and Schär, 2010;Pereira et al. 2017). Future cold spells using a stationary cold day definition are identified using time-invariant thresholds estimated from the baseline, while those under a nonstationary definition are diagnosed using time-variant thresholds calculated for each 30-year future period at 10-year time-step from the baseline to 2071-2100. The blocking-cold spell linkage is mainly investigated using the nonstationary cold day definition, as it is difficult to evaluate with the stationary definition due to the significant decrease in the number of cold spells in a warming climate. For each blocking region from B1 to B5, the winter blocking-cold spell relationship is quantified by the cold spell frequency anomaly P ano = (N CS|B ∕N B )∕(N CS ∕N) (Brunner et al. 2018), where N , N CS , N B , and N CS|B represent total number of days, number of cold spell days, number of blocking days, and number of cold spell days during blocking in winter for an analysis period. Consequently, this anomaly indicates the ratio of the conditional probability of cold spell days during blocking to the unconditional probability. It is larger than one when the conditional probability is larger than the unconditional (climatological) occurrence probability.
Statistical significance of the cold spell frequency anomaly is tested for each ensemble member using a block bootstrap resampling approach that is suggested to deal complex correlated structure of dataset (Efron and Tibshirani, 1994). Autocorrelation in blocking and cold spells is conserved using contiguous 15-day bootstrap resamples that could include a long blocking episode as well as a cold spell. The resampling is repeated 100 times from the series of daily blocking in each blocking region and cold days at each North American grid point. Statistical significance of the cold spell frequency anomaly is defined when the 5th-95th percentiles range quantified by the block bootstrap approach excludes one. The anomaly is therefore either statistically significantly higher than climatology when the 5th percentile is larger than one or lower than climatology when the 95 th percentile is smaller than one. Anomalies of surface heat fluxes and 1000-500 hPa geopotential thickness are also tested their statistical significance from zero using the bootstrap resampling approach.

Results
It should be noted again that abilities of CanESM2 and CanESM5 ensembles in representing historical blockingcold spell linkages has been reported by Jeong et al. (2021) through comparison with ERA-Interim and NCEP-DOE-R2 for the 1981-2010 baseline period. Therefore, this section mainly focuses on describing future projections of the two ensembles for blocking-cold spell linkages and relevant atmospheric circulation characteristics.

Changes in atmospheric blocking and cold spell
CanESM2 and CanESM5 ensembles generally agree with ERA-Interim and NCEP-DOE-R2 for the spatial variability of winter blocking frequency for the baseline, with larger percentages of blocking days from B1 to B3 (i.e., North Pacific, Alaska, and western Canada) compared to B4 and B5 (i.e., central and eastern Canada) (Fig. 1a, b). This spatial variability of the Pacific-North American blocking was also well documented in previous studies (e.g., Whan et al. 2016;Woollings et al. 2018). However, the two ensembles underestimate the blocking frequency of the two reanalyses over the North Pacific at B1 and B2, consistent with previous findings using CMIP5 models including CanESM2 (e.g., Anstey et al. 2013; Lee and Ahn 2017) as well as using CMIP6 models including CanESM5 (Davini and D'Andrea 2020). In particular, ERA-Interim and NCEP-DOE-R2 yield 30 and 26 blocking days per winter, while CanESM2 and CanESM5 ensembles simulate 15 and 20 blocking days on average in the North Pacific. CanESM5 shows slightly better performance than CanESM2 for blocking frequency, which is consistent with other reported improvements by CanESM5 to reproduce large-scale circulation patterns in the climate system (Swart et al. 2019;Jeong et al. 2021).
The two ensembles are generally consistent in projected changes of winter blocking days over the North Pacific, as they show a gradual decrease in the North Pacific from the baseline to 2071-2100, with a larger decrease in CanESM5 than CanESM2 (− 0.62 and − 0.25 day/decade) (Fig. 1b). The two ensembles project smaller changes at regions B3-B5. In terms of signal-to-noise ratio, the ratio of ensemble mean change to the ensemble spread (defined by the standard deviation), the signal in CanESM2 and CanESM5 ensembles exceeds noise at B1 (since the periods 2051-2080 and 2001-2030, respectively) and B2 (since the periods 2071-2100 and 2021-2050, respectively), and the ensemble mean changes increase as the projection horizon increases. However, the ensemble spread in blocking frequency of both CanESM2 and CanESM5 also decreases by about 18% over the North Pacific for the 2071-2100 period compared to the baseline. The decrease of future winter blocking frequency in the North Pacific is in line with previous findings from projections of CMIP5 and CMIP6 GCMs, under conditions of enhanced warming, both in the tropics and over the Arctic (e.g., Dunn-Sigouin and Son, 2013; Kennedy et al. 2016;Woollings et al. 2018;Davini and D'Andrea 2020).
The two ensembles also project similar spatial patterns in the changes for IB, SIB, and TSIB for the 2071-2100 period relative to the baseline (From 1st to 3rd panels of Fig. 1c), with statistically significant decreases with regard to internal variability (tested by two-sample t test at 90% temporally persistent SIB (TSIB), and meridional potential vorticity gradient (PVy) at 500 hPa of CanESM2 and CanESM5 ensemble averages (EAs) in winter for the 2071-2100 period relative to the baseline. The values of IB, SIB, and TSIB are percentages of the differences between the future and baseline periods relative to the baseline in frequency. Grid points are marked by black dots when the projected changes are statistically significant based on the two-sample t test at 10% significant level confidence level), in the North Pacific, western Alaska, and southwestern and eastern coastal regions of Canada. However, the two ensembles show statistically significant increases for IB, SIB, and TSIB relative to the internal variability in central and southeastern parts of Canada for the 2071-2100 period, with larger increases for CanESM2. This indicates that the change in blocking (i.e., TSIB) is mainly driven by the changes in IB and associated changes in spatial extent and temporal persistence. The spatial patterns of projected changes are consistent with findings of Kennedy et al. (2016), Woollings et al. (2018), and Davini and D'Andrea (2020) from CMIP5 and CMIP6 GCMs. Accordingly, it is important to understand how the climate change influences the background condition and subsequent blocking frequency. Recent theoretical studies (e.g., Luo et al. 2019;Luo and Zhang, 2020) indicated that the persistence or duration of blocking is directly related to the magnitude of background meridional potential vorticity gradient (PVy). In particular, when the PVy shows an increase in the future relative to the baseline, the blocking persistence or duration shows a decrease that result in a decrease in future blocking frequency (see Luo et al. (2019) for more details). The changes in the PVy at 500 hPa of the two ensembles (bottom panel of Fig. 1c) generally show an increase over the Pacific region but a decrease over the North American continental region in the future (2071-2100) relative to the baseline, which might support the changes in the blocking frequency over those regions.
The two ensembles show an agreement with the reanalyses for regional averages of cold spell days, with larger values at NW and SW regions compared to NE and SE regions (Fig. 2a). The two ensembles also show an agreement for the spatial distribution of cold spell frequency with larger values at NW and central parts (i.e., Great Plains) of North America for the baseline (Fig. 2b). In particular, CanESM2 and CanESM5 yield 0.10 day (2.1%) and − 0.39 day (− 8.1%) of biases with ERA-Interim, respectively, and 0.50 day (11.1%) and 0.002 day (0.05%) of biases with NCEP-DOE-R2 for cold spell days per winter, on average over North America for the baseline. Climatological cold spell days are determined by surface temperature anomalies that occur in response to the persistence of synoptic-scale circulation anomalies and land-atmosphere interactions (Bieli et al. 2015;Grotjahn et al. 2016). It is notable that the climatological probability of the cold spell days represents the denominator (unconditional probability) of the cold spell frequency anomaly P ano .
Under the stationary cold day definition, the two ensembles exhibit a gradual but significant decrease in cold spell days in four regions of North America as the projection Unfilled and filled box-plots represent stationary and nonstationary cold day threshold conditions. b Spatial plots of CanESM2 and CanESM5 ensemble averages (EAs) for winter cold spell days over North America for the baseline and their projected changes under nonstationary cold day threshold for the 2071-2100 period with respect to the baseline. Grid points are marked by black dots when the projected changes are statistically significant based on the twosample t test at 10% significant level horizon increases; faster rates of convergence to zero are seen in the north (i.e., NW and NE) as compared to the south (i.e., SW and SE) (Fig. 2a). These decreases are due to increases in future air temperature above the baseline cold day threshold in response to increased greenhouse gas concentrations. In addition, larger increases in future air temperature over northern North America under the highemission scenarios is a reflection of polar amplification (IPCC 2013; Bush and Lemmen 2019), which results in larger decreases in cold spell days in the north. However, the two ensembles still project the occasional occurrence of cold spells until 2071-2100, which is in line with the previous findings (e.g., Kodra et al. 2011;IPCC 2013;Gao et al. 2015). CanESM5 shows a faster decrease (e.g., 8% larger decreases during two decades from the baseline) in cold spell days compared to CanESM2. This is because CanESM5 has a higher climate sensitivity and warms more rapidly than CanESM2 in response to the external forcings (Swart et al. 2019). For instance, CanESM5 ensemble show 1.4 °C larger increase of global mean temperature than CanESM5 ensemble for 2071-2100 period compared to the baseline.
Under the nonstationary cold day threshold, the two ensembles project gradual decreases in cold spell days at the four regions as the projection horizon increases (Fig. 2a).
The two ensembles show a consistent regional response with, on average, larger decreases in the NE and NW (0.15 and 0.13 day/10-year) than those in the SE and SW (0.08 and 0.05 day/10-year). They also display an agreement in their spatial patterns of projected changes between the future (2071-2100) and baseline periods (Fig. 2b). These decreases are likely due to a reduction in the persistence of future cold spells associated with changes in the persistence of synopticscale circulation anomalies. It is remarkable that the total number of cold days per winter in each future period is constant based on the nonstationary cold day threshold.

Changes in blocking-cold spell linkage
In the baseline period, CanESM5 shows better agreement with ERA-Interim and NCEP-DOE-R2 than CanESM2 in terms of the percentages of cold spell events that coincide with blocking occurring over all regions considered (Fig. 3a). This is mainly due to better representation of blocking frequency by CanESM5, particularly in the North Pacific (B1 and B2) (Fig. 1a, b). However, the two ensembles do still exhibit negative biases relative to the reanalyses for the baseline, with the largest bias at the SE region, mainly due to underestimated blocking frequency over the North Pacific (Fig. 1a, b). The two ensembles generally agree in Fig. 3 a Regional averages of the percentage of number of cold spell events coincident with any blocking days occurring from B1 to B5 relative to total number of cold spell events for ERA-Interim (black dots), NCEP-DOE-R2 (green dots), CanESM2 ensemble (red boxplot), and CanESM5 ensemble (blue box-plot) in the four regions (NW, NE, SW, and SE) for the 1981-2010 baseline period, and projections of the two ensembles for future 30-year periods. Unfilled and filled box-plots represent stationary and nonstationary cold day threshold conditions. b Spatial plots of CanESM2 and CanESM5 ensemble averages (EAs) for the percentage for the 1981-2010 baseline period and their projected changes (nonstationary threshold) for the 2071-2100 period relative to the baseline. Grid points are marked by black dots when the projected changes are statistically significant based on the two-sample t test at 10% significant level the spatial distribution of cold spell events during blocking, with larger values over southern NW and NE, northern SW, and most of the SE (Fig. 3b).
Under the stationary cold spell definition, the two ensembles project the percentage of cold spell events that occur during blocking to remain similar to the baseline until certain future periods (i.e., 2021-2050at NW, 2001-2030at NE, 2061-2090at SW, and 2021-2050 (Fig. 3a). This indicates that, despite significant decreases in the number of cold spells due to future warming, blocking occurring over Pacific-North America will result in similar percentages of cold spells for a few decades. In particular, the historical cold spell-blocking percentage in CanESM2 persists longer into the future than does that for CanESM5.
Under the nonstationary cold day definition, the two ensembles project decreases in the cold spell-blocking percentages, with larger decreases in CanESM5 than in CanESM2 (− 8.0% and − 1.8% points, respectively) for 2071-2100 period relative to the baseline (Fig. 3).
Regionally, CanESM5 ensemble projects statistically significant decreases (tested by two-sample t-test at 90% confidence level in term of internal-variability) in the percentage at NW, NE, and SE regions, while CanESM2 ensemble projects some insignificant decreases at NE and SE regions but also increases at western North America for the future period (Fig. 3b). These decreases in the percentages of cold spell events coincident with blocking in the future are driven both by decreases in blocking frequency over the North Pacific and decreases in cold spells over North America. In particular, CanESM5 shows larger decreases at NW and NE regions, based on larger decreases in the blocking frequency and cold spell events at these regions compared to the other regions (Figs. 1, 2).
Blocking-cold spell linkages (i.e., cold-spell frequency anomalies on blocking occurring in the five regions) of the two ensembles are presented in maps in Fig. 4 and are compared with those of ERA-Interim and NCEP-DOE-R2 using regional averages at four regions of North America for the Fig. 4 Winter cold-spell frequency anomaly ( P ano ) over North America during blocking in the five regions (B1-B5; black boxes) of CanESM2 and CanESM5 EAs for the baseline (left panels) and their projected changes (nonstationary threshold) for the 2071-2100 period relative to the baseline (right panels). Black dots of left panels represent statistically significant grid points at the 10% significant level tested by the bootstrap approach and the larger dot size represents the larger number of ensemble members yield statistical significance (percentages of the number are defined in the legend), while those of right panels show statistical significance based on the two-sample t test at 10% significant level baseline in Fig. 5. It is notable that the linkages are based on frequencies of blocking and cold spell days as defined by P ano . The two ensembles display a statistically significant connection between blocking and above-normal cold spell frequency at remote locations on eastern and/or southern flanks from each blocking region for the baseline, which is in line with previous findings (e.g., Masato et al. 2014;Brunner et al. 2018). These connections between blocking and regional above-normal cold spell frequency are statistically significant for both ensembles, particularly over NW and NE during blocking on B1 and B2 and over SW and SE during blocking on B2-B4, as the regional frequency values of the ensembles are generally larger than one (i.e., climatological condition) (Fig. 5). Meanwhile, the two ensembles exhibit a statistically significant connection between below-normal cold spell frequency at the core and surrounding areas of blocking where high-pressure anomalies and blocking-like circulation patterns are dominant. The two ensembles overestimate the above-normal cold spell frequency at NW and NE regions during blocking at regions B1-B3, due to overrepresented negative temperature advection at NW and negative surface heat flux at NW and NE during the blocking (Jeong et al. 2021). The two ensembles, however, display a consistency with the two reanalysis datasets for regional averages in the other cases, while CanESM5 shows a better agreement with the reanalyses than does CanESM2.
Under the nonstationary cold day threshold, projected changes of the two ensembles for the blocking-cold spell linkages are compared using maps (Fig. 4) and regional averages at four regions of North America between the baseline and future 2071-2100 period (Fig. 5). Note that this assessment is not considered under the stationary cold day threshold due to the significant decrease in cold spell frequency with future warming (Fig. 2a). Basically, the two ensembles project similar spatial distributions of the blocking-cold spell linkages for the future period relative to those for the baseline. This is consistent with the finding over the Euro-Atlantic sector that blocking will continue to play an important role in the development of winter cold spells in the future (Brunner et al. 2018). The two ensembles project fewer cold spell events than climatology over northwestern North America (i.e., Alaska and/or northwestern Canada) during blocking in the North Pacific (i.e., region B2) in the future; however, the linkage is weaker than in the baseline. This is due to weaker future blocking activity in the North Pacific as shown in Fig. 1b. For the cold spell frequency anomalies, the two ensembles generally show weaker anomalous signals (Fig. 4) and larger spreads in regional averages (Fig. 5) for the future period compared to the baseline. This indicates that the two ensembles simulate larger uncertainty associated with internal-variability in the future blockingcold spell linkages over the Pacific-North America.

Changes in temperature advection, radiation and turbulent heat fluxes during blocking
The two ensembles show strong negative temperature advection anomalies at 850 hPa on the east and south flanks of the blocking for the baseline period (Fig. 6a). These negative anomalies are driven by northerly and/or north-easterly wind anomalies due to blocking events and associated anomalous anti-cyclonic circulations. This spatial behavior of the negative advection anomalies during blocking in Pacific-North America is consistent with the patterns reported in the Euro-Atlantic sector (e.g., Trigo et al. 2004;Sousa et al. 2018). On the other hand, the two ensembles yield positive advection anomalies on the west flank of blocking regions B3-B5 for the baseline, induced by the anomalous anti-cyclonic circulation and associated southerly wind anomalies. In addition, the two ensembles show positive net heat flux anomalies in the northern parts of blocking regions B2-B5 and their east and west sides (Fig. 6b), due to the anti-cyclonic high-pressure of blocking that bring positive surface radiation anomalies with clear-sky conditions at the core and surrounding areas of the blocking. Meanwhile, the two ensembles generally display negative net heat flux anomalies over the other land areas of North America. Consequently, the combined influence of temperature advection and net heat flux anomalies over North America during blocking in regions B1-B5 can explain the blocking-cold spell linkages for the baseline period shown in Figs. 4 and 5. The two ensembles generally agree with ERA-Interim and NCEP-DOE-R2 in terms of regional averages of these anomalies during blocking for the baseline (Fig. 7). However, relative to the reanalysis datasets, the two ensembles over-represent regional averages of negative temperature advection at NW and negative net surface heat flux at NW and NE during blocking at regions B1-B3 for the baseline period. This explains the Fig. 6 a Temperature (850 hPa) advection anomalies with wind vector anomalies (arrows) (during blocking in regions B2, B3, and B5) of CanESM2 and CanESM5 EAs for the baseline and 2071-2100 future periods. b Net surface heat flux anomalies (during blocking in regions B2, B3, and B5) of CanESM2 and CanESM5 EAs for the baseline (left panels) and their projected changes (nonstationary threshold) for the 2071-2100 period relative to the baseline (right panels). Black dots of left panels represent statistically significant grid points at the 10% significant level tested by the bootstrap approach and the larger dot size represents the larger number of ensemble members yield statistical significance (percentages of the number are defined in the legend), while those of right panels show statistical significance based on the two-sample t test at 10% significant level. Cases for blocking regions B1 and B4 are not shown overestimation of the linkage in the two ensembles between blocking in regions B1-B3 and above-normal cold spell frequencies over NW and NE (Fig. 5).
The two ensembles project similar spatial patterns of temperature advection and net heat flux anomalies during blocking for the 2071-2100 period compared to those for the baseline (Fig. 6). Note that the future anomalies during blocking are taken relative to the corresponding future 30-year climatologies. Similar temperature advection and heat flux anomalies in the baseline and future periods supports their similar spatial patterns in the blocking-cold spell linkages between the two periods shown in Fig. 4. However, spatially distributed temperature advection and net heat flux anomalies of the two ensembles during blocking (Fig. 6) and their regional averages (Fig. 7) are generally weaker for the future period compared to those for the baseline, particularly over NW and NE for blocking B1-B3 and over SW and SE for blocking B4 and B5. These weaker anomalies in the future result in the weaker future blocking-cold spell linkages with larger uncertainty associated with internal-variability (Figs. 4, 5). The weaker net surface heat flux anomalies over North America in the future are mainly due to the weaker anomalous signals of downward radiative components during blocking events (Fig. 8a). The two ensembles also project weaker anomalies for 1000-500 hPa geopotential thickness during blocking in the future (Fig. 8b). Thinner (thicker) geopotential thickness than climatology generally implies colder (warmer) air temperature than climatology. Therefore, weaker geopotential thickness anomalies support weaker temperature anomalies and associated weaker blocking-cold spell linkages during blocking over North America.

Summary and discussion
This study is the first to investigate future changes in the linkages between atmospheric blocking and winter (December-February) cold extremes over Pacific-North America. This complements past studies for the Euro-Atlantic sector (e.g., Sillmann and Croci-Maspoli 2009;Brunner et al. 2018). This study also addresses the mechanisms responsible for these changes by investigating future blocking activities and synoptic-scale variables that are directly related to temperature anomalies, including temperature advection and surface heat fluxes. This investigation has been conducted using two large initial condition ensembles of successive generations of Canadian Earth System Models (CanESM2 and CanESM5) under high-emission scenarios. Uncertainties in changes of the cold spell-blocking linkages that stem from internal variability are evaluated using the largeensemble of simulations.
The two ensembles project gradual decreases in winter blocking frequency from the 1981-2010 baseline to 2071-2100 over the North Pacific with larger decreases in CanESM5 (− 3.08%/decade) than CanESM2 (− 1.73%/ decade) (Fig. 1b). The decreases of future winter blocking frequency over the North Pacific are consistent with previous Fig. 8 a Anomalies of the radiation components of the net surface heat flux (i.e., shortwave plus longwave downward solar radiations) and b those of geopotential thickness (between 500 and 1000 hPa) (during blocking in regions B2, B3, and B5) of CanESM2 and CanESM5 EAs for the baseline and 2071-2100 future periods. Cases for blocking regions B1 and B4 are not shown. Black dots of left pan-els represent statistically significant grid points at the 10% significant level tested by the bootstrap approach and the larger dot size represents the larger number of ensemble members yield statistical significance (percentages of the number are defined in the legend), while those of right panels show statistical significance based on the twosample t test at 10% significant level findings from CMIP3 to CMIP6 GCMs (e.g., Dunn-Sigouin and Son 2013; Lee and Ahn 2017;Woollings et al. 2018;Davini and D'Andrea 2020). Based on enhanced warming both in the tropics and over the Arctic, blocking activity is projected to be weaker over the North Pacific with an increase in background meridional potential vorticity gradient (PVy) at 500 hPa, while stronger over the northern North America with a decrease in PVy in the future (Fig. 1c).
Under a stationary cold day definition, the two ensembles project gradual but statistically significant decreases in cold spell frequency with increasing projection horizon due to increases in future daily minimum temperatures with warming but fixed cold day threshold defined from the baseline period. The two ensembles show faster rates of decreases in cold spell frequency over northern than southern parts of North America, due to larger increases in future air temperature over the region under the high-emission scenarios (IPCC 2013; Bush and Lemmen 2019). Nevertheless, the two ensembles continue to project the occasional occurrence of cold spells over North America out to 2071-2100 as reported by previous studies (e.g., Kodra et al. 2011;IPCC 2013;Gao et al. 2015). CanESM5 yields a faster decrease in the number of cold spell days, as it warms more rapidly in response to external forcing compared to CanESM2 (Swart et al. 2019). Despite the significant decreases in cold spell frequency, blocking occurring in the Pacific-North America region plays a similarly important role in the development of cold spells over North America for 4-8 decades from the baseline period (Fig. 3a).
Under the nonstationary cold day definition, the CanESM2 and CanESM5 ensembles project modest decreases in winter cold spell days over North America (− 2.0 and − 2.3%/decade, respectively) with larger decreases over the north (− 2.87%/decade) than the south (− 1.56%/decade). This could be caused by decreases in the persistence of future cold days associated with the changes in the persistence of synoptic-scale atmospheric circulation anomalies, given the same number of cold days for each future period based on the nonstationary cold day threshold. The two ensembles project decreases in winter cold spells coincident with blocking, with larger decreases in CanESM5 (13%) than those in CanESM2 (3%) for the 2071-2100 period relative to the baseline (Fig. 3b). These decreases are accompanied by both decreases in blocking frequency in north Pacific and cold spells over North America. The two ensembles project similar patterns in the blocking-cold spell linkages for 2071-2100 period compared to the baseline (Fig. 4). This is in line with the finding over Euro-Atlantic sector by Brunner et al. (2018) that blocking will continue to play an important role in the development of winter cold spells under a nonstationary cold day definition. However, the two ensembles generally show weaker significance in cold spell frequency anomalies during blocking and larger spreads in their regional averages for the future period (Figs. 4,5), indicating larger uncertainty in the future blocking-cold spell linkages. Meanwhile, the two ensembles project similar spatial patterns between the baseline and future periods for temperature advection and net heat flux anomalies during blocking (Fig. 6). This supports the similar spatial distributions in the blocking-cold spell linkages of the two ensembles in the baseline and the future periods under the nonstationary cold day definition. Nevertheless, these anomalies of the synoptic-scale circulation variables during blocking are also weaker for the future; hence, this results in weaker blocking-cold spell linkages with larger uncertainty due to internal variability in the future.
CanESM5 shows better agreement than CanESM2 for the blocking frequency and the blocking-cold spell linkages with reanalysis datasets at the Pacific-North America (see also Jeong et al. 2021); this is consistent with overall improvements of CanESM5 in simulating large-scale circulation and climate patterns (Swart et al. 2019). Larger projected changes for the blocking-cold spell linkages are suggested by CanESM5, although the two generation GCMs show a consistency in projected changes for the linkages and associated synoptic-scale atmospheric circulations during blocking. Ongoing research using global climate models under different future emission scenarios is needed to evaluate the impacts of blocking on other climate extremes over North America, such as heat waves and dry spells, which are usually co-located with blocking.