Substantial regional climate change expected following cessation of CO₂ emissions

Here we examine the full complexity ESMs that participated in ZECMIP and have submitted their results to Earth System Grid Federation. Nine ESMs have done so: ACCESS-ESM1-5 [16, 17], CanESM5 [18], CESM2 [19, 20], GFDL-ESM4 [21], GISS-E-2.1 [22], MIROC-ES2L [23], MPI-ESM1-2-LR [24, 25], NorESM2-LM [26], and UKESM1-0-LL [27]. Notably this list of nine ESMs is different from the nine that were analyzed by MacDougall et al [7] in the initial ZECMIP analysis paper. In the initial analysis output from GFDL-ESM2M was analyzed, while for the submission to Earth System Grid Federation these simulations were replaced with the new version of the model—GFDL-ESM4. An additional model, GISS-E2.1, completed the ZECMIP experiments too late to be analyzed by MacDougall et al [7]. Additionally, as of August 2022, the ZECMIP results from CNRM-ESM2 have yet to be submitted to Earth System Grid Federation.

Succinct model descriptions for aspects of the model relevant for understanding ZEC for the seven models examined here and in MacDougall et al [7] are detailed in the appendices MacDougall et al [7]. A full description of the GFDL-ESM4 model is given by Dunne et al [21] and of GISS-E2.1 in Kelley et al [22]. Because the biogeochemical components of ESM are computationally expensive, many ESMs use lower resolution variants of their ocean and atmosphere components. Thus MIROC-ES2L, MPI-ESM1-2-LR, NorESM2-LM, and UKESM1-0-LL all use lower resolution versions of their atmosphere and ocean components than the general circulation variants of these models (See appendices of [7]). Notably lower resolution model variants are expected to have somewhat different patterns of warming relative to high resolution versions of the same model (e.g. [28]).

The ZECMIP protocol [10] outlined six experiments to quantify ZEC and its dependence on total emissions of CO₂ and CO₂ emissions pathway. Only one of these experiments, esm-1pct-brch-1000PgC, was requested from all participating full complexity ESMs. The esm-1pct-brch-1000PgC experiment branches from either the standard 1pctCO2 experiment [11] after models diagnose that 1000 PgC has been emitted to the atmosphere from burning of fossil fuels, or from the esm-1pctCO2 which is emission driven from the start of the simulation [10]. Neither the 1pctCO2 nor the esc-1pctCO2 experiments simulate land-use-change. Thereafter the models switch from concentration driven to emissions driven configurations and CO₂ emissions are set to zero. Thus the atmospheric CO₂ pool is allowed to freely evolve in response to changes in oceanic and terrestrial carbon pools [10]. Although this experiment is highly idealized [29] the instantaneous cessation of emissions allows for clear distinction between ZEC and the transient climate response to cumulative CO₂ emissions [7, 30]. ESMs were requested to submit simulations of at least 100 years following cessation of emissions. While multiple ensemble members were not requested, two models—CESM2 and UKESM1-0-LL, did each submit three ensemble members to Earth System Grid Federation.

To quantify ZEC at large spatial scales 2 m air temperature anomalies following cessation of emissions were computed for the land and oceans, as well as zonal averages for tropical, extratropical, and polar regions. For all computations the protocol devised for computing the global value of ZEC was followed [7]. That is, ZEC is taken as relative to the 20 year mean temperature from the 1pctCO2 experiment centred about the point in time where the esm-1pct-brch-1000PgC branches from the 1pctCO2 experiment [7]. Due to the sudden halt in the rate of temperature change in the esm-1pct-brch-1000PgC experiment, using the 1pctCO2 experiment in place of the esm-1pct-brch-1000PgC experiment is necessary to avoid an underestimation of the global temperature at the time of emission cessation [7]. The land-ocean mask provided for each model is used to separate land and ocean 2 m air temperatures and compute areally weighted mean values. For zonal averages Earth is divided into a tropical region extending from 30 ^∘S to 30 ^∘N, two extratropical regions one in each hemisphere between 30 ^∘ to 60 ^∘, the Arctic region above 60 ^∘N, and an Antarctic region below 60 ^∘S. Note that we have used divisions approximating the extent of the Hadley, Farrell, and Polar atmospheric circulation cells [31], instead of the astronomical definitions of these regions.

For regional climate our principal aim is to compute ZEC at the grid-point scale and investigate whether the patterns observed are due to internal variability or are being generated by persistent dynamically driven changes to regional climate. Following MacDougall et al [7] our main focus is on 50 years after emissions cease, a metric called ZEC₅₀. To derive maps of ZEC₅₀ we take a 20 year average for each grid-point centred at 50 years after emissions cease, and find the difference between these grid-point averages and the 20 year mean grid-point temperature from the 1pctCO2 experiment centred around the year emissions cease. To explore whether the resultant patterns are distinguishable from internal variability we use the following statistical metric:

$$\begin{equation} A = \frac{Z_{L}-Z_{G}}{\sigma_{L}}, \end{equation} \tag{ 1 }$$

where, A is a climate anomaly expressed as units of standard deviations (formally dimensionless), Z_L is the local ZEC, and Z_G is the globally averaged ZEC value of the respective model. σ_L is the local standard deviation of differences between pairs of 20 year average temperature anomalies. That is, σ_L quantifies how much variance one would expect between any two (non-overlapping) 20 year average temperatures at a given location. σ_L was calculated from the pre-industrial control simulation (pi-control) for each model. The control simulations were first detrended at each grid-point and then discreetly averaged into non-overlapping 20 year mean temperatures. Thereafter differences were found between every combination of 20 year mean temperatures at that grid-cell location and the standard deviation taken of these differences to compute σ_L . Providing a minimum 500 year pre-industrial control simulation (pi-control) wherein models are evolving under constant pre-industrial forcing was a requirement for entry into CMIP6 [11], and thus such a simulation exists for all models that participated in ZECMIP. Note that we have used that standard pi-control simulation where CO₂ concentration is held constant. Evaluated over the whole planetary surface the metric can give a sense for whether a dynamic shift has taken place or that internal variability has increased relative to pre-industrial times. Generally the higher fraction of the surface exhibiting large absolute values of the metric indicates that it is more likely that such a shift has occurred. For example assuming a Gaussian (normal) distribution we would expect that 5% of the surface would have absolute values over 2, if we instead find that 25% of the surface has absolute values above 2 then it is likely that either dynamic changes have occurred or the internal variability has increased relative to pre-industrial times. The metric is a slight variant of the commonly used Standardized Mean Difference statistical metric. For conciseness we call this metric 'Regional Standard Difference'. The metric can be computed for any given time after emissions cease.

A more conventional way of removing internal variability and establishing the robustness of dynamic changes in ESM output is to average together many ensemble members from the same model conducting the same numerical experiment (e.g. [32]). However only two of the models submitted more than one ensemble member, and these models submitted only three members each. Sample sizes of three are not kind to the reliability of formalized statistical metrics, however the range of values between these ensemble members is nonetheless informative about the persistence of dynamic changes.

Figure 1 shows ZEC for the land and ocean area for each model. In general the response over land is more extreme and more variable than over the ocean. Over land three of the models (CESM2, MPI-ESM1-2, and MIROC-ES2L) have cooled by more than 0.5 ^∘C by the end of their simulations, after 150 years, while over the ocean no model breaches that cooling threshold. Other models such as ACCESS-ESM1-5 and GISS-E2.1, show a more moderate response with both land and ocean ZEC staying close to zero. One model (UKESM1-0-LL) shows substantial warming of similar magnitude over both the land and the ocean, although with much higher variability over land. The ZEC₅₀ values for land and ocean are shown in table 1. Seven of the models have negative ZEC₅₀ over the land and the ocean, with ZEC₅₀ over the land being on average 0.15 ^∘C cooler than over the ocean. Note that due to a slightly different collection of ESMs than that assessed in MacDougall et al [7] the multi-model, global mean value of ZEC is lower than that previously assessed. As this is largely due to one of the higher ZEC models being slow to upload their data to Earth System Grid Federation this should not be interpreted as a meaningful change.

Zoom In Zoom Out Reset image size

Simulations with general circulation models long ago established that land areas are expected to warm more than the ocean even when the climate system has come into full equilibrium with a change in radiative forcing [33]. This persistence anomaly was diagnosed to be caused by inter-linked changes in humidity, lapse-rate and moisture transport affecting local radiative transfer in the atmosphere over the continents [33]. How these processes will interact with changes in ocean heat uptake in a scenario of zero emissions is unclear and likely a fruitful avenue for further research.

Figure 2 show the zonal averages for ZEC for each model. The figure shows clear zonal differences in ZEC, with ZEC being closer to zero, less variable, and exhibiting smaller trends in the tropical region compared to other regions (recall that the tropics cover half the planetary surface). The southern hemisphere mid-latitudes show remarkably little inter-annual variability and ZEC values close to zero, however, continuing warming trends in CESM2 and UKESM1-0-LL are clear. Note, that the southern hemisphere mid-latitudes are dominated by the Southern Ocean, and contain very little land. The northern hemisphere mid-latitudes are more variable and are generally consistent with ZEC over land (figure 1(a)), which is unsurprising given the extent of continents in this region. Both Arctic and Antarctic exhibit high variability, both in time and between models, with more variability in the Arctic relative to the Antarctic. It is noteworthy that the three models experiencing strong negative ZECs over land after 150 years, are the ones showing the strongest cooling over the Arctic region.

Zoom In Zoom Out Reset image size

Table 2 shows the ZEC₅₀ values for the zonal regions. The changes in the Arctic are most striking with some models showing a ZEC close to zero (MIROC-ES2L and CanESM5), others showing almost 1 ^∘C of warming (UKESM1-0-LL), and another over 1 ^∘C of cooling (CESM2). The bounds of Arctic ZEC extend between $-2\,^{\circ}$C and 2 ^∘C (figure 2). The tropics in contrast exhibit a range of ZEC₅₀ from −0.31 to $0.21\,^{\circ}$C. Five of the nine models agree in the sign of change in ZEC₅₀ for all zonal regions (CanESM5, GFDL-ESM4, MPI-ESM1-2-LR, NorESM2-LM, and UKESM1-0-LL), while the other four show broad zonal regions of warming and cooling.

Maps of ZEC₅₀ for each ESM are shown in figure 3. The patterns suggest that many regions can expect substantial climate change following cessation of emissions. At the extremes, UKESM1-0-LL has a grid-cell with a ZEC of 5.10 ^∘C and GISS-E2.1 has a grid-cell with a ZEC₅₀ of $-4.35\,^{\circ}$C. Each ESM exhibits a unique pattern of regional warming and cooling with few consistent features between models. These range from relatively uniform cooling in NorESM2 to patchy warming and cooling in a number of models, to the strong inter-hemispheric gradient in CESM2. General patterns that can be spotted are a stronger cooling over land compared to the ocean, as well as a high ZEC response in the Arctic region, confirming the signals found in the zonal means of section 3.1.

Zoom In Zoom Out Reset image size

Figure 4 shows changes in temperature 50 years after emissions cease relative to the change in temperature at the time emissions cease. The figures show how climate change after emissions cease compares to the climate change already experienced before emissions cease. These relative changes show somewhat more consistent patterns with all of the models showing anomalies of various sizes and magnitudes in the North Atlantic, with these anomalies being most extensive in CESM2 and GISS-E2.1. The other region where large relative changes occur in all models after cessation of emissions is the Southern Ocean, with an expanse in the South Pacific being a particularly anomalous region in most models. Overall very little of the planetary surface is expected to return to pre-industrial temperatures (median 1.7%, range of 0.1%–3.78% of the planetary surface cooling to pre-industrial temperature), and likewise only a small fraction of the planetary surface is expected to warm more in the 50 years after emission cease than warmed by time emission ceased (median 0.5%, range of 0%–1.9% of the planetary surface warming more in the 50 years after emissions cease than before).

Zoom In Zoom Out Reset image size

The diverse patterns seen in the spatially resolved ZEC₅₀ in the model output leads to the question of whether the patterns that are seen are persistent or are being driven by internal variability. The best way to answer such a question would be to examine a large ensemble of simulations from each model conducting the experiment. However, only a small ensemble (three members each) is available for two of the models. Figure 5 shows the ensemble mean and range for CanESM5 and UKESM1-0-LL. Both ensemble means show persistent patterns between the three ensemble members. CanESM5 exhibits warming in the North Atlantic and cooling in the Southern Ocean and UKESM1-0-LL exhibits strong Arctic warming. The range between the ensemble members shows that both models exhibit more internal variability over the Arctic, Antarctic, and Southern Ocean, and smaller and consistent variability between 45 ^∘ North and 30 ^∘ South. These results suggest that the patterns observed in figure 3 are created by a mix of internal variability and persistent dynamically driven changes in regional climate.

Zoom In Zoom Out Reset image size

Figure 6 shows maps of the Regional Standard Difference (see section 2.3) for each of the models. Two of the models, ACCESS ESM1-5 and GFDL-ESM4, show values within the −2 to 2 standard deviation range that would be expected for Gaussian (normal) distributions. The other ESMs show extensive regions outside this range, particularly in the Southern Ocean and North Atlantic. In addition, UKESM1-0-LL has values far outside the −2 to 2 standard deviation range in the Arctic, and CESM2 has extensive regions outside this range of opposite signs in both hemispheres. Figure 7 shows the fraction of the planetary surface above a given absolute value of regional standard difference for each model, along with a reference value for a Gaussian distribution of the same metric. Figure 7(a) shows this for an example from the pre-industrial climate and figure 7(b) for 50 years after emissions cease. For the pre-industrial climate values are close to the Gaussian expectation, with some models above and some below this line. For the ZEC₅₀ climate all models are either equal to or above the Gaussian expectation. This indicates that for most models the ZEC climate is exhibiting either an increase in variability relative to the pre-industrial climate, large dynamically driven changes in regional climate, or a combination of both.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size