<figure><img src="image_1.png" /></figure>

HAL
open science

# Reducing uncertainty in local temperature projections

Saïd Qasmi, Aurélien Ribes

► To cite this version:

Saïd Qasmi, Aurélien Ribes. Reducing uncertainty in local temperature projections. Science Advances, 2022, 8 (41), 10.1126/sciadv.abo6872. hal-03977751

HAL Id: hal-03977751

https://hal.science/hal-03977751

Submitted on 8 Feb 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
CLIMATOLOGY

Reducing uncertainty in local temperature projections

Saïd Qasmi* and Aurélien Ribes

Planning for adaptation to climate change requires accurate climate projections. Recent studies have shown that the uncertainty in global mean surface temperature projections can be considerably reduced using historical observations. However, the transposition of these new results to the local scale is not yet available. Here, we adapt an innovative statistical method that combines the latest generation of climate model simulations, global observations, and local observations to reduce uncertainty in local temperature projections. By taking advantage of the tight links between local and global temperature, we can derive the local implications of global constraints. The model uncertainty is reduced by 30% up to 70% at any location worldwide, allowing to substantially improve the quantification of risks associated with future climate change. A rigorous evaluation of these results within a perfect model framework indicates a robust skill, leading to a high confidence in our constrained climate projections.

Copyright © 2022
The Authors, some
rights reserved;
exclusive licensee
American Association
for the Advancement
of Science. No claim to
original U.S. Government
Works. Distributed
under a Creative
Commons Attribution
License 4.0 (CC BY).

INTRODUCTION

As the global mean temperature keeps rising and climate change intensifies, there is a growing demand for local-scale monitoring of current and future climate change. Assessing and planning the adaptation to the expected unprecedented impacts of climate change on human activities, ecosystems, and the biosphere as a whole require an accurate local information with well-calibrated uncertainties. This need relates to estimates of warming to date and the future warming in response to set of scenarios of future greenhouse gas emissions.

It is unequivocal that human influence has warmed the atmosphere, ocean, and land since preindustrial times (1). Concurrently, the anthropogenic influence is not detected everywhere at the local scale (2, 3). Natural climate variability can blur the emergence of the anthropogenic signal for the next few years at high latitudes, while a substantial warming is already reported in several tropical regions (4, 5). Regarding climate projections, the Intergovernmental Panel on Climate Change (IPCC) concluded in its fifth assessment report (AR5) (6) that “Future [human-induced] warming trends cannot be predicted precisely, especially at local scales”.

In the IPCC AR6 (7), a new generation of climate models (8) has been used to provide a range of projections in response to different socioeconomic scenarios (9). On the basis of this new dataset, various studies have recently shown that uncertainty in global mean warming can be considerably reduced using the information provided by recent observed warming trends via so-called “constraint” methods (10–13). These studies consistently point toward a downward revision of the expected warming in all emission scenarios (10, 12), with a decrease in model uncertainty of nearly 40% for end of century projections (11), and even more at shorter lead times. This is an important result as, until then, observations have failed to provide clear evidence in reducing the range of climate projections (14).

The next challenge is to transpose these new findings on global warming to regional and local scales. At the regional scale, a few studies have adopted the partitioning from the Special Report on Managing the Risks of Extreme Events and Disasters to Advance Climate Change Adaptation (SREX) (15) and have attempted to narrow model uncertainty with sophisticated techniques with promising results (16, 17). However, the SREX regions are typically continental-wide and do not provide relevant information for local adaptation. At the local

scale (defined as the size of a global climate model grid box of about 200 km) and to the best of our knowledge, only a few studies have attempted to narrow climate model uncertainty, by using weighting methods to account for interdependencies between models (18, 19) or by focusing on specific and limited areas (20–22). In particular, although constrained projections of global mean temperature are now used in the IPCC AR6 (7), local climate projections are still solely based on a raw ensemble of available climate models (https://interactive-atlas.ipcc.ch/), derived from global warming levels.

Here, we assess how much uncertainty in local temperature projections can be reduced. We first take advantage of the tight links that exist between local climate and global mean temperature (23, 24). Specifically, we describe the local implications of the recent advances in the reduction of the uncertainties in global mean temperature projections. We then provide a set of local-scale temperature projections, which encapsulate another source of information: the observed local warming to date. If compared to the global mean temperature record, internal variability is larger in local observations. However, they still provide a useful source of information about both past and future trends, particularly over some specific regions. We discuss how much these two types of observations (global and local) narrow uncertainty on future warming ranges. Such a reduction is expected to provide more accurate information that becomes critical for policy-makers in the local climate risk management (25), as well as for the climate science community.

RESULTS

The Kriging for Climate Change (KCC) method used by Ribes et al. (11, 16) is one of the statistical techniques that have led to a significant reduction of uncertainty in probabilistic projections of global mean temperature by combining observations and models (7, 11). This Bayesian method involves three steps. First, the response to all external forcings for each climate model considered is estimated over the period 1850–2100, after filtering out internal variability as much as possible. Second, the sample of forced responses from available climate models is used as a prior of the real-world forced response for each grid point. This is done assuming that “models are statistically indistinguishable from the truth.” Third, observations are used to derive a posterior distribution, i.e., a constrained temperature response, of the past and future forced response given observations, in a Bayesian way.
Compared to many emergent constraint approaches in which
observations are often summarized into one single variable (26, 27),
the KCC method is able to take the full observed time series into
account. Here, we further extend this technique to account for mul-
tiple time series and potential dependencies between them (see
Materials and Methods). When applied to the global surface air
temperature (GSAT) time series simulated by the models from the
Coupled Models Intercomparison Project phase 6 (CMIP6) (8) models
and the Shared Socioeconomic Pathway (SSP) 5 to 8.5 scenario, the
amplitude of best estimate of the projected GSAT changes constrained
by the observations is revised downward by 0.5°C by 2100, with a
reduction in model uncertainty of nearly 40% (11). Here, we con-
sider global mean surface temperature (GMST; a blending of land
air and ocean sea surface temperatures; see Materials and Methods)
instead of GSAT, as GMST is more consistent with local observa-
tional record. Figure 1A offers an update of the GMST constraint:
the warming of 5.3°C projected by CMIP6 is in this case revised
downward by 0.4°C (best estimate) in 2100. Minor differences with
Ribes *et al.* (11) are explained by 2 years of additional observations,
by the addition of several CMIP6 models, which affect the prior dis-
tribution, and by the lower warming observed in GMST compared
to GSAT (28, 29).

Constrain local climate projections with global observations

Climate models exhibit a strong correlation between future GMST and local warming over most regions of the globe (Fig. 1B). To take such a relationship into account, we extend the KCC method to constrain local temperature projections. Beyond the simple correlation shown in Fig. 1B, this method uses all the information contained in the entire observed GSAT time series to derive local warming (considering that the annual time series provides useful additional information, e.g., to distinguish between greenhouse gases (GHG) and aerosol forcings). This is done by deriving the local warming conditional on the observed GMST record (hereafter the GMST-only case; see Materials and Methods). As an example, we consider the North American city of Dallas for which the simulated local temperature over the 2022–2100 period is significantly correlated with future GMST (see the corresponding point in Fig. 1B). Consistent with GMST results, the local temperature range constrained by GMST observations indicates a decrease in uncertainty of about 20% over the 2021–2040 period, up to 30% over the 2081–2100 period (Fig. 2A) in the GMST-only case. The best estimate of local warming is revised downward by 0.4°C by 2100 compared to the unconstrained projections (hereafter the unconstrained case).

When the method is applied to any location worldwide, the re-
sults in the projected mean temperature and in model uncertainty
show a clear relationship with the level of correlation between the
local temperature and GMST (Figs. 3, B and D, and 1B). Note that
the method is applied for each location separately, i.e., that it does
not take advantage of spatial correlations and does not provide con-
strained projection at a larger regional scale. The reduction of un-
certainty in local projections is the highest at the locations where
the correlation with GMST is the strongest. For these locations, e.g.,
over several continental regions, the North Pacific and the Indian
Ocean, a reduction of the ensemble spread of about 45% is obtained
over the 2021–2040 period, with a downward revision of the best
estimate warming between 0.2° and 1°C (Fig. 3, A and C). Conversely,
for locations where the correlation is low, such as in southern Indian
Ocean and the Barents Sea, the local temperature response is weakly

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** GMST time series and its correlation with local temperature. (A) GMST annual observations from the HadCRUT5 (48) dataset (black points) are used to constrain concatenated historical and SSP5-8.5 scenario simulations of GMST. The unconstrained (pink) and constrained (red) ranges stand for the 5 to 95% confidence interval of the forced response as estimated from 27 CMIP6 models. The thick pink (red) line stands for the ensemble mean (best estimate). All values are anomalies with respect to the 1850–1900 period. (B) Intermodel correlation between simulated GMST trends over the 2022–2100 period and local temperature trends over the 2022–2100 period. Stippling indicates regions with nonsignificant correlation (*P* > 0.05 based on a two-sided Student's *t* test).</figcaption></figure>

constrained, with a reduction of the model uncertainty of 10% and
a revision of the best estimate by 0.5°C or less. These revised ranges
lead to a warming pattern at +2°C of global warming, considerably
different from that assessed in the AR6 (Fig. 4, A and B, and fig.
S1C). For example, the constrained local temperatures over North
America are expected to be 0.4°C warmer compared to the uncon-
strained case. Note that, by construction, only the spatial pattern of
warming is affected by the observational constraint in this case,
since the +2°C warming level is let unchanged—so, the map of dif-
ferences necessarily mixes patches of positive and negative values
here, even if warming ranges at a given date are all revised down-
ward (see caption of Fig. 4).

Added value of local observations to the constraints

Beyond the useful information provided by the historical GMST
time series, it is insightful to assess the consistency between the ex-
pected local response (regardless whether observed GMST is accounted
for) and local historical observations. Current and past warming is
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Observational constraints on historical and SSP5-8.5 local temperature changes in Dallas. (A) Constrained local temperature for the grid point at (48.75°N; 11.25°E) (see blue point in Fig. 1B) in the GMST-only case. The 5 to 95% constrained spread of the simulated response to all external forcings is in red, the thick line stands for the best estimate. The spread in light blue and the associated thick line is relative to the unconstrained case. (B) Constrained local temperature in the Local + GMST case. The 5 to 95% constrained spread of the simulated forced response is in blue, the thick line stands for the best estimate. Black points are the observations. Note that because single-year observations are affected by internal variability, they often lie outside the 5 to 95% assessed ranges of the forced responses. The orange line stands for the smoothed observations time series using 6° of freedom. Unconstrained simulated ranges are identical between (A) and (B). All values are anomalies with respect to the 1850–1900 period.</figcaption></figure>

spatially heterogeneous, and some regions such as the Arctic are warming faster than others (Fig. 3A) (30). Therefore, it is relevant to account for both GMST and local observations to provide local projections consistent with all available observations. Using recent local observations could particularly affect short-term projections (typically over the 2021–2040 period) and could provide a different picture of the constrained temperature ranges.

To make such a calculation, we derive a posterior of the expected local warming given local historical observations in addition to the GMST observations (hereafter the Local + GMST case; see Materials and Methods). Following the example of Dallas considered in Fig. 2A, the constrained local temperature ranges become typically closer to long-term local observed variations (Fig. 2B compared to Fig. 2A). The constraint by local observations is found to influence the best estimate of future warming even on longer lead times. Regarding the uncertainty, the added value of local observations in the reduction of the model spread compared to the GMST-only case is limited in this example, with a decrease of about 10% of the confidence range width compared to the GMST-only case. Two reasons contribute to this limited impact and must be considered for any location. First, the local signal-to-noise ratio can be small. This may happen if local internal variability or measurement uncertainty is large (i.e., local observations provide little insight on the externally forced response). Second, the global and local responses can be highly correlated with each other so that they partly provide the same information, leading to a limited impact of local observations on uncertainty ranges. In both cases, the model uncertainty will be only marginally reduced by local historical data (see Materials and Methods).

The application of the Local + GMST constraint to all grid points worldwide results in a projected warming pattern, which remains quite close to the GMST-only case (Fig. 3, C and E) but with regional

differences. On the one hand, for several regions over the Arctic, the warming is revised downward compared to the unconstrained case, making the projections more consistent with recent observations and implying a reduced warming compared to the GMST-only case. On the other hand, an upward revision is obtained over Eastern Asia, Greenland, the East Siberian region, and Southern oceans. The added value of local observations in the reduction of model uncertainty is largest over these regions where the correlation in Fig. 1 is low (Fig. 3F). The estimated projections for the 2081–2100 period can also be derived and are consistent with these results (fig. S2). Note that even if each grid point is treated independently from the others (see Materials and Methods), the global mean of the constrained local ranges (for both the GMST-only and the Local + GMST cases) is very close to the constrained GMST ranges shown in Fig. 1A, suggesting a consistency between the different types of constraints.

The addition of local information can clearly modify the warming pattern at +2°C of global warming (Fig. 4, B and C; and fig. S1, B and C). For example, while a downward revision of the temperature change of -0.2°C is obtained over Europe in the GMST-only case, an upward revision of 0.3°C is obtained in the Local + GMST case. This change of sign is widespread over Eurasia. In the context of an urgent need of adaptation to the threat of climate change, our constrained warming pattern provides a revised and a more relevant information for local adaptation planning.

## Evaluation of the constrained projections

The robustness of these promising results is quantified within a so-called perfect model framework, using a leave-one-out cross validation (see Materials and Methods). Each member of each model is considered as pseudo-observations over the 1850–2021 period. These are subsequently used to constrain the temperature projections, using all other models as a prior. The constrained temperature
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Unconstrained and constrained local temperature projections over the 2021–2040 period. (A) Ensemble mean of the unconstrained local temperature changes. (B) Ensemble spread of the unconstrained local temperature changes, defined as the 5 to 95% confidence interval of the multimodel ensemble. (C) Difference of local temperature changes ensemble mean between the GMST-only case and the unconstrained case. (D) Relative difference of local temperature changes ensemble spread between the GMST-only case and the unconstrained case. (E) Same as (C) but for the Local + GMST case. (F) Relative difference of local temperature changes ensemble spread between the Local + GMST case and the GMST-only case, illustrating how much incorporating local observations narrows uncertainty. All values are anomalies with respect to the 1850–1900 period.</figcaption></figure>

range is then compared to the projected warming simulated by the model from which pseudo-observations were taken. As making this evaluation for all the grid points is computationally expensive, this procedure is applied to several locations, considered as representative of the diversity of the worldwide climate. As for the real observations, we assess both the GMST-only and Local + GMST constraints. The continuous ranked probability skill score (CRPSS) (31) is used to measure the accuracy of the method, taking the

unconstrained projections as a baseline as a first step (see Materials and Methods).

Figure 5A shows that the average of the CRPSS distribution based on all pseudo-observations is positive for every location in the GMST-only case with an improvement of about 30% over the 2021–2040 period. Depending on the location, the skill is remarkably improved by 10 to 50% (for the average) (Fig. 5 and figs. S5 to S9). In the Local + GMST case, the average skill is also positive (Fig. 5C), and for most
<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 4.** Warming pattern at +2°C of GMST warming. (A) Ensemble mean of the unconstrained local temperature changes. (B) Difference of local temperature changes between the GMST-only case (best estimate) and the unconstrained case (ensemble mean). (C) Difference of local temperature changes between the Local + GMST case (best estimate) and the unconstrained case (ensemble mean). By construction, global averaged differences are zero, as the GMST threshold remains at +2°C with or without constraint (contrary to Fig. 3). All values in (A) to (C) are anomalies with respect to the 1850–1900 period.</figcaption></figure>

locations, the CRPSS lies between 10 and 60% relative to the unconstrained case.

These results clearly demonstrate the performance of the method over short lead times (2021–2040). Moreover, the comparison between the Local + GMST and the GMST-only constraints indicates that the skill is slightly improved when adding local observations to constrain projections. The average of the CRPSS distribution is positive for 52 locations over 55 when comparing the Local + GMST and the GMST-only constraints (Fig. 5D and see boxplots in magenta in figs. S5 to S9) and is slightly negative for the remaining points. The significance of this result is assessed with a binomial test. Under the null hypothesis that adding local observations has no impact on the skill (i.e., that the GMST-only and Local + GMST cases are categories equally likely, such as a coin toss), the probability of getting this result by pure chance is less than 0.1%. This suggests that there is a clear added value in considering the constrained ranges derived from the Local + GMST case relative to the GMST-only case. A third case for which we only use local observations (Local-only case) to constrain projections indicates lower scores than in the Local + GMST case (Fig. 5B) and confirms that using the combination of global and local observations enhances the accuracy of the method.

Similar results are obtained when we consider the 2081–2100 period (fig. S11).

The rare cases where the CRPSS is negative are due to the large low-frequency variability in few CMIP6 models. This is a topic of concern, as several CMIP6 model are characterized by clear multidecadal and even centennial internal variability in GMST (32, 33). Figures S12 to S16 show that models associated with a strong decadal variability (e.g., CNRM-CM6-1, EC-Earth3, and IPSL-CM6A-LR) are those contributing negatively to the CRPSS in most cases. The assumptions used in the KCC method to statistically model internal variability might explain the failure to capture such a slow internal variability, resulting in overconfidence. Excluding those specific models from the evaluation process leads to a clear increase in the CRPSS by 15% both in the GMST-only and Local + GMST cases (fig. S17). Discussing the realism of these particular models requires further analyses on internal variability (34, 35) and climate sensitivity (36), which are beyond the scope of this study. Still, regardless of whether the real system contains such low frequency variability, the CRPSS distributions remain mostly positive across locations for the three types of constraint, even if models with large internal variability are included. Both distributions in fig. S17 (excluding models with
<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 5.** Spatial distribution of the CRPSS for the constrained temperature projections. The values indicate the average values from the CRPSS distributions (based on all pseudo-observations; see Materials and Methods) for the GMST-only (A), Local-only (B), and Local + GMST (C) constraints, respectively, compared to the unconstrained case for 55 locations (number in black over each point). (D) The added value of the local observations and stands for the CRPSS average in the Local + GMST case compared to the GMST-only case. For each location, the average is weighted by the number of members in each model. Calculations are made over the 2021–2040 period.</figcaption></figure>

large internal variability) are above the distributions (including all models) in fig. S10, which themselves show mostly positive values (see the median value). Last, the relevance and reliability of our method seems robust to including all CMIP6 models and would be even strengthened if low-variability models were proven less realistic.

This is supported by a second evaluation criterion of the method based on coverage probabilities, which lead to similar conclusions (see Supplementary Discussion). From all of these evaluation results, we retain the Local + GMST case to provide guidance in constraining local projections. The evaluation of the KCC method suggests that the constrained temperature ranges are reliable and demonstrate that relying on unconstrained projections to assess the local future climate is no longer the best approach.

## DISCUSSION

We have shown, using a statistical method combining the entire temperature observation records with model simulations, that uncertainty in local temperature projections can be substantially narrowed. Local projections constrained by both global and local observations exhibit a reduction of the uncertainty of 40% on average by 2100. This demonstrates the benefits of merging model simulations with observations to provide the best picture of future

climate change. Figure 6 offers a complementary perspective to the AR6 (7) conclusions that were solely based on raw (unconstrained) projections. For each location, a temporal evolution from 1850 to 2100 of the constrained temperature and its uncertainty can be derived, with revised projections for the near and the long-term time scales. In particular, the KCC method provides a way out of the concept of global warming levels, by estimating the uncertainty for a given date. This fills the gap in the IPCC atlas (only based on unconstrained projections) and provides a considerable revision of the local exposure to the consequences of the on-going climate change (37). An online tool that implements the method and illustrates the constrained temperature ranges for every point over a horizontal grid of 2.5° resolution is available via the following demonstrator: https://saidqasmi.shinyapps.io/KCC-shinyapp/.

Promising prospects exist to improve the constrained projections. CMIP5 (38) and CMIP6 (8) ensembles sample model uncertainty in a probabilistic way using all climate models as an “ensemble of opportunities” (39, 40). This approach has several limitations that can bias the estimation of climate uncertainty (41). First, they are not designed to sample uncertainty comprehensively (42), e.g., no models are being intentionally built to have low or high climate sensitivities. Forcing uncertainty is also poorly sampled in the CMIP6 ensemble, e.g., the magnitude of the aerosol forcing in 2014
<figure><img src="image_13.png" /><figcaption>**Fig. 6.** Mean temperature change at a +2°C GMST warming. Best estimate of the constrained local temperature changes in the Local + GMST case. Similarly to Fig. 2B, the constrained and unconstrained temperature ranges are shown for several world capitals cities over the 1850–2100 period. All values are anomalies with respect to the 1850–1900 period.</figcaption></figure>

(43). Second, each model output is considered as independent and contributes equally to the multimodel ensemble, although it is known that many CMIP6 models share common components and parametrizations (44). This “model democracy” paradigm has been largely used to summarize projection information in IPCC assessment reports (6), although it can be criticized (45). Therefore, using a subset of models qualified as independent a priori or weighting the models in this way (46, 47), before applying the observational constraint, may provide even more reliable results.

Our results demonstrate that available observations offer valuable information to reduce uncertainty in climate projections over the next decades. Even with the continuous improvement of climate

models, the intermodel spread is not necessarily reduced for several variables from one CMIP generation to the next (1). In this sense, the contribution of observations and constraint methods is expected to improve the reliability of the projections. Therefore, it is critical to account for this new source of information and to regularly bridge the gap between monitoring recent changes and predicting future changes.

The KCC method itself can also be improved. Although it can be used on larger areas to easily derive constrained projections, e.g., on the SREX regions (16), the current implementation does not take into account the spatial dependence in the climate variability between locations. Taking the spatial dimension fully into account
could bring additional useful information and would result in consistent uncertainties at all spatial scales. This requires dealing with much larger covariance matrices. Estimating those matrices is a key challenge, as the number of models is much smaller than the number of grid points and years. Reducing the dimension of the problem, with methods used in data assimilation or considering spatio-temporal principal components may help to regularize these high-dimensional matrices.

Generalizing the method to other variables of high societal impacts, e.g., extreme precipitation, droughts, and snow cover, some of which are also tightly related to GMST, would also be very relevant. For some of these variables, the observed data are sparse, which requires finding well-sampled covariables over the historical period. In addition, variables such as sea ice or precipitation do not necessarily follow Gaussian distributions, which makes it necessary to adapt the KCC method to other types of distributions. In this way, the climate science community could take a step forward toward a more accurate assessment of past and future human-induced climate change.

## MATERIALS AND METHODS

### Observational dataset and models

The temperature observations are from the HadCRUT5 (48) dataset over the 1850–2021 period. The temperature field comes from a blending of near-surface temperature and sea surface temperature using land sea mask and sea ice concentration. The measurement uncertainty of the HadCRUT5 dataset is estimated from an ensemble of 200 equiprobable realizations. Most of the other observational products are included in the temperature range estimated by this ensemble, which confirms our choice to consider the HadCRUT5 dataset as a reference.

CMIP6 models are selected according to the availability of the following data: at least 200 years of a preindustrial control (piControl) simulation; at least one member of a historical simulation and one member of a projection simulation for the SSP5-8.5 scenario. To constrain the simulated temperatures at a grid point scale in a consistent way, a blended temperature field $T_{\text{blend}}$ is computed in each CMIP6 model based on the formulation of Morice et al. (48)

$$ \begin{cases} w_{\text{air}} = (1 - f_{\text{ocean}}) + f_{\text{ocean}} f_{\text{ice}} \\ T_{\text{blend}} = w_{\text{air}} T_{\text{air}} + (1 - w_{\text{air}}) T_{\text{ocean}} \end{cases} \quad (1) $$

where $T_{\text{air}}$, $T_{\text{ocean}}$, $f_{\text{ice}}$, and $f_{\text{ocean}}$ are for each grid point near-surface air temperature, sea surface temperature, sea ice concentration, and sea area fraction. The 27 models for which these variables are available and which satisfy the above criteria are listed in Table 1.

We define the GSAT as the global mean of $T_{\text{air}}$ and the GMST as the global mean of $T_{\text{blend}}$. Several studies have shown that GMST and GSAT clearly differ as GMST warms less than GSAT (28, 29).

Models are interpolated on a common horizontal grid of 2.5° resolution before calculating blended temperatures and applying the constraining method. This choice is motivated by a compromise between the different resolutions of the CMIP6 models (between 1.5° and 2.5°). Note that the KCC method can be applied to finer resolutions if observations are available at this scale. For temperature, for which the spatial autocorrelation is high, the reduction in the uncertainty is expected to be similar as for the 2.5° resolution.

### Statistical method

The statistical method is based on the same one used by Ribes et al. (11), whose formulation and principle is similar to kriging, which is a method originally developed to interpolate geophysical data based on prior covariances. In Ribes et al. (11), this method is applied to the analysis of time series from climate simulations of CMIP5 and CMIP6 models and is used for several purposes: (i) reducing model uncertainty on past and future global warming estimated by CMIP and ScenarioMIP (9) simulations, (ii) reducing uncertainty on warming attributed to several external forcings via the Detection and Attribution Model Intercomparison Project (DAMIP) (49) models, and (iii) complete or statistically reconstruct missing simulations from other physically relevant simulations (e.g., using the so-called 1% CO₂ simulations in which the CO₂ concentration increases by 1% each year, to reconstruct DAMIP historical simulations in which greenhouse gases follow their historical concentrations, while other forcings are kept constant). Here, we apply this method of KCC to

<table><caption>Table 1. List of the available CMIP6 models and the associated number of members in the historical and SSP5-8.5 simulations used to constrain temperature projections.</caption><thead><tr><td>Model simulation</td><td>Historical + SSP5-8.5</td></tr></thead><tbody><tr><td>ACCESS-CM2</td><td>1</td></tr><tr><td>ACCESS-ESM1-5</td><td>3</td></tr><tr><td>BCC-CSM2-MR</td><td>1</td></tr><tr><td>CanESM5-CanOE</td><td>3</td></tr><tr><td>CanESM5</td><td>50</td></tr><tr><td>CESM2</td><td>2</td></tr><tr><td>CESM2-WACCM</td><td>1</td></tr><tr><td>CNRM-CM6-1</td><td>6</td></tr><tr><td>CNRM-CM6-1-HR</td><td>1</td></tr><tr><td>CNRM-ESM2-1</td><td>5</td></tr><tr><td>EC-Earth3</td><td>7</td></tr><tr><td>EC-Earth3-Veg</td><td>3</td></tr><tr><td>FGOALS-g3</td><td>1</td></tr><tr><td>FIO-ESM-2-0</td><td>3</td></tr><tr><td>GFDL-ESM4</td><td>1</td></tr><tr><td>HadGEM3-GC31-LL</td><td>1</td></tr><tr><td>INM-CM4-8</td><td>1</td></tr><tr><td>INM-CM5-0</td><td>1</td></tr><tr><td>IPSL-CM6A-LR</td><td>6</td></tr><tr><td>MIROC6</td><td>3</td></tr><tr><td>MIROC-ES2L</td><td>1</td></tr><tr><td>MPI-ESM1-2-HR</td><td>2</td></tr><tr><td>MPI-ESM1-2-LR</td><td>10</td></tr><tr><td>MRI-ESM2-0</td><td>1</td></tr><tr><td>NESM3</td><td>1</td></tr><tr><td>NorESM2-MM</td><td>1</td></tr><tr><td>UKESM1-0-LL</td><td>5</td></tr><tr><td>27 models</td><td>121 members</td></tr></tbody></table>
reduce the model uncertainty in the past and future temperature changes simulated by CMIP6 models at each grid point. Note that a confusion can be made with techniques based on so-called emergent constraints methods (26, 27). Emergent constraints would usually consider the sole observed global warming trend (a single scalar); e.g., over the 1980–2021 period, to constrain the simulated temperature changes in the future. The KCC method has several advantages compared to this approach. Instead of simply constraining a trend over a subperiod, it uses the entire observed time series of temperature, which avoids ignoring useful information. In addition, the method takes into account the model temporal pattern uncertainty and provides confidence ranges specifically for the forced response, while many other studies also include internal variability.

For a given grid point, we define $\mathbf{y}^*_{\text{loc}}$ as the yearly time series of the real (and unknown) temperature response to external forcings over the 1850–2021 period and $\mathbf{y}_{\text{loc}}$ as the observed yearly temperature time series over the same period. Similarly, we define for the GMST the vectors $\mathbf{y}^*_{\text{glo}}$ and $\mathbf{y}_{\text{glo}}$ as the unknown response to external forcings and the observed time series, respectively. They constitute the following $\mathbf{y}^*$ and $\mathbf{y}$ vectors, both of size $2n_y$ where $n_y = 170$

$$ \mathbf{y}^* = \begin{pmatrix} \mathbf{y}^*_{\text{loc}} \\ \mathbf{y}^*_{\text{glo}} \end{pmatrix}, \mathbf{y}^* = \begin{pmatrix} \mathbf{y}^*_{\text{loc}} \\ \mathbf{y}^*_{\text{glo}} \end{pmatrix} \quad (2) $$

Assuming that the observed temperature total variability can be decomposed as the sum of a term of forced variability and a term including both internal variability and measurement errors, $\mathbf{y}$ takes the following form

$$ \mathbf{y} = \mathbf{y}^* + \boldsymbol{\epsilon} \quad (3) $$

where $\boldsymbol{\epsilon} = (\boldsymbol{\epsilon}_{\text{loc}}, \boldsymbol{\epsilon}_{\text{glo}})$ is a vector of size $2n_y$ and corresponds to the local and global terms of measurement errors and observed internal variability. Further assuming that models are indistinguishable from the truth, i.e., that observations and models are exchangeable (50–52), observations $\mathbf{y}$ can be rewritten

$$ \left\{ \begin{array}{l} \mathbf{x} = \begin{pmatrix} \mathbf{x}_{\text{loc}} \\ \mathbf{x}_{\text{glo}} \end{pmatrix} \\ \mathbf{y} = \mathbf{H}\mathbf{x} + \boldsymbol{\epsilon} \end{array} \right. \quad (4) $$

where $\mathbf{x}_{\text{loc}}$ and $\mathbf{x}_{\text{glo}}$ are the yearly time series over the 1850–2100 period of the local and global temperature responses to external forcings estimated in CMIP6 models, respectively, i.e., vectors of size $n_x = 251$. $\mathbf{H}$ is an observation operator of size $2n_y \times 2n_x$, which extracts the part of $\mathbf{x}$ that is observed in $\mathbf{y}$, i.e., the forced response from 1850 to 2021, and whose form depends on the type of the applied constraint (using only GMST observations or both GMST and local observations; see eq. S21). Note that the assumption of exchangeability between observations and models has been suggested as well supported by observations, especially for temperature (50, 53).

For a given CMIP6 model $m$ listed in Table 1, we choose to estimate the simulated response to all external forcings $\mathbf{x}_{m,glo}$, by decomposing the simulated GMST over 1850–2100 into an anthropogenic response $\mathbf{x}_{m,ant,glo}$, and a natural response $\mathbf{x}_{m,nat,glo}$. Therefore, after averaging all available members of the model $m$, the simulated GMST time series over 1850–2100 $\mathbf{x}_{m,glo}$ writes

$$ \mathbf{x}_{m,glo} = \mathbf{x}_{m,nat,glo} + \mathbf{x}_{m,ant,glo} + \boldsymbol{\epsilon}_m \quad (5) $$

where $\boldsymbol{\epsilon}_m$ is a random term for internal variability.

To estimate $\mathbf{x}_{m,nat,glo}$ and $\mathbf{x}_{m,ant,glo}$ in the model $m$, we use a generalized additive model (GAM) to compute the response to all external forcings, $\mathbf{x}_{m,glo}$ (recall that $\mathbf{x}_{m,glo}$ follows Eq. 4)

$$ \left\{ \begin{array}{l} \mathbf{x}_{m,glo} = \mathbf{x}_{m,all,glo} + \boldsymbol{\epsilon}_m \\ \mathbf{x}_{m,all,glo} = \frac{\boldsymbol{\beta}_m \mathbf{e} + f(\mathbf{t})}{\mathbf{x}_{m,nat,glo} + \mathbf{x}_{m,ant,glo}} \end{array} \right. \quad (6) $$

where $\boldsymbol{\beta}_m$ is an unknown scaling factor. $\mathbf{e}$ is a vector of size $n_x$ and is the temperature response to the natural forcings computed from a two-layer (atmosphere-ocean) energy balance model (EBM) following equations 1 and 2 of Geoffroy et al. (54), using the historical and SSP5-8.5 natural forcings between 1850 and 2100 estimated by the Priestley Center (55) as a radiative term. Here, $\mathbf{e}$ is calculated using an average of EBM parameters fitted to the CMIP6 ensemble and aims at estimating rapid year-to-year variations of natural forcings. $f(\mathbf{t})$ is a time series [with $\mathbf{t} = (1850, \dots, 2100)$] and refers to an assumed smoothed response of GMST to the anthropogenic forcings (i.e., a smoothed response of $\mathbf{x}_{m,glo}$). The function $f$ corresponds to smoothing splines to filter out part of internal variability, with 6° of freedom (a value that was selected as a bias-variance trade-off). The total forced responses as estimated by this procedure are illustrated in fig. S18.

We apply the exact same procedure to estimate the local forced responses as simulated by each CMIP model. For each grid point from the model $m$, we consider $\mathbf{x}_{m,loc}$, the average of all available members, to estimate the local forced response, $\mathbf{x}_{m,all,loc}$. We assume that the local natural response scales linearly with the globally averaged natural forcings time series, as the EBM response $\mathbf{e}$ used to calculate $\mathbf{x}_{m,nat,glo}$ is also used when we fit the GAM to compute the local natural response $\mathbf{x}_{m,nat,loc}$. Thus, $\mathbf{x}_{m,nat,glo}$ and $\mathbf{x}_{m,nat,loc}$ only differ by their scaling factor $\boldsymbol{\beta}_m$. We believe that our results are not sensitive to this choice given the reduced strength of, and uncertainty in, the natural response compared to the anthropogenic response.

The multimodel ensemble of the local and global simulated responses to all external forcings is used to derive a distribution of $\mathbf{x}$, noted $\Pi(\mathbf{x}) \sim \mathcal{N}(\boldsymbol{\mu}, \Sigma_{\text{mod}})$, built from all $\mathbf{x}_{m,glo}$ and $\mathbf{x}_{m,loc}$. $\boldsymbol{\mu} = (\boldsymbol{\mu}_{\text{loc}}, \boldsymbol{\mu}_{\text{glo}})$ is a vector of size $2n_x$ and is the multimodel ensemble mean of the concatenated local and global forced responses. $\Sigma_{\text{mod}}$ is a variance-covariance matrix of size $2n_x \times 2n_x$ that describes the model spread, with the following form

$$ \Sigma_{\text{mod}} = \begin{bmatrix} \Sigma_{\text{mod,loc}} & \Sigma_{\text{mod,dep}} \\ \Sigma_{\text{mod,dep}} & \Sigma_{\text{mod,glo}} \end{bmatrix} \quad (7) $$

where $\Sigma_{\text{mod,loc}}$ and $\Sigma_{\text{mod,glo}}$ are the sample covariance matrices of size $n_x \times n_x$ modeling local and global model spread within $\mathbf{x}_{\text{loc}}$ and $\mathbf{x}_{\text{glo}}$, respectively. $\Sigma_{\text{mod,dep}}$ is the covariance matrix modeling the dependence between $\mathbf{x}_{\text{loc}}$ and $\mathbf{x}_{\text{glo}}$.

In our Bayesian framework, $\Pi(\mathbf{x})$ is a first (probabilistic) estimate of $\mathbf{x}$, which makes no use of observations, and is only based on climate models. We want to update this estimate by incorporating the observational evidence provided by $\mathbf{y}$. Following the Bayesian theory, the calculation of the posterior distribution $p(\mathbf{x} | \mathbf{y})$ is required. A prerequisite is to define the observational uncertainty, i.e., the covariance matrix associated with $\mathbf{y}$.
## Modeling of observational uncertainty

Given Eq. 4, we assume that $\epsilon \sim N(0, \Sigma_{obs})$, where $\Sigma_{obs} = \Sigma_{meas} + \Sigma_{iv}$ is the observation error covariance matrix. $\Sigma_{meas}$ and $\Sigma_{iv}$ are both of size $2n_y \times 2n_y$ and describe the measurement error and internal variability, respectively. $\Sigma_{meas}$ is estimated as the sample covariance matrix over the 200-member ensemble of the HadCRUT5 dataset (48).

$\Sigma_{iv}$ is estimated using observed annual time series of global and local temperature over the 1850–2021 period. First, we compute the global observational residuals by subtracting the CMIP6 response to all external forcings $\mu_{glo}(1,...,n_y)$ from the observations $y_{glo}$. Similarly, we derive local residuals by subtracting $\mu_{loc}(1,...,n_y)$ from $y_{loc}$. These residuals constitute a first estimate of global and local internal variability, noted $\hat{\epsilon}_{iv,loc,1}$, $\hat{\epsilon}_{iv,glo,1}$, respectively.

We define $\Sigma_{iv}$ as a matrix of size $2n_y \times 2n_y$ of the following form

$$ \Sigma_{iv} = \begin{bmatrix} \Sigma_{iv,loc} & \Sigma_{iv,dep} \\ \Sigma'_{iv,dep} & \Sigma_{iv,glo} \end{bmatrix} \quad (8) $$

where $\Sigma_{iv,loc}$ and $\Sigma_{iv,glo}$ are the covariance matrices of size $n_y \times n_y$ modeling local and global internal variability within $y_{loc}$ and $y_{glo}$, respectively. $\Sigma_{iv,dep}$ is the covariance matrix modeling the dependence between local and global internal variability, i.e., $\epsilon_{iv,loc}$ and $\epsilon_{iv,glo}$.

To compute $\Sigma_{iv}$, we take into account decadal internal variability that exists in the global (56), regional (57), and even local (58) observations, using a mixture of two autoregressive processes or order 1 (AR1), hereafter mixture of autoregressive processes (MAR), as done by Ribes et al. (11). The MAR formulation includes a fast (f) and a slow (s) components such that global internal variability $\epsilon_{iv,glo}$ within the GMST residuals writes at a time t

$$ \left\{ \begin{aligned} \epsilon_{iv,glo}(t) &= \epsilon_{iv,f,glo}(t) + \epsilon_{iv,s,glo}(t), \\ \epsilon_{iv,f,glo}(t) &= \alpha_{f,glo} \epsilon_{iv,f,glo}(t-1) + Z_{f,glo}(t), \\ \epsilon_{iv,s,glo}(t) &= \alpha_{s,glo} \epsilon_{iv,s,glo}(t-1) + Z_{s,glo}(t) \end{aligned} \right. \quad (9) $$

where the parameters $\alpha_{s,glo}$ and $\alpha_{f,glo}$ are the lag 1 coefficients of the AR1 processes and $\alpha_{s,glo} \ge \alpha_{f,glo}$ by convention. $Z_{s,glo}(t) \sim N(0, \sigma_{s,glo}^2)$ and $Z_{f,glo}(t) \sim N(0, \sigma_{f,glo}^2)$ are white noises associated with the two AR1. The slow component is able to generate a dependence on time scales of typically one decade, while the fast component accounts for interannual variability. Following the principle of parsimony, only four coefficients ($\sigma_{f,glo}^2, \alpha_{f,glo}, \sigma_{s,glo}^2$, and $\alpha_{s,glo}$) are thus needed to characterize internal variability at the global scale and to make $\Sigma_{iv,glo}$ invertible. We fill the covariance matrix $\Sigma_{iv,glo}$ following the calculations of each of its coefficients, as detailed in eq. S8. In practice, we apply a maximum likelihood procedure to the local and global residuals according to the statistical model from Eq. 9. Uncertainty related to these coefficients is not taken into account. Then, we make the same assumptions, and estimate four other parameters, ($\sigma_{f,loc}^2, \alpha_{f,loc}, \sigma_{s,loc}^2$, and $\alpha_{s,loc}$), to characterize fast and slow components in local internal variability $\epsilon_{iv,loc}$ and to compute $\Sigma_{iv,loc}$.

The initial estimate of $\Sigma_{iv}$, noted $\hat{\Sigma}_{iv,1}$ is solely based on the residuals $\hat{\epsilon}_{iv,loc,1}$ and $\hat{\epsilon}_{iv,glo,1}$ derived from the unconstrained forced response. This first estimate is likely flawed as the real (and unknown) forced response $y^*$ is not necessarily consistent with the unconstrained forced response estimated by $\mu$. In addition, as $\mu$ can be by construction different from the best estimate of the constrained forced response $\hat{\mu}_1$ (the mean of the posterior distribution $p(x|y)$,

see section below), the residuals $\hat{\epsilon}_{iv,loc,1}$, $\hat{\epsilon}_{iv,glo,1}$ before constraint are not always coherent with the residuals $\hat{\epsilon}_{iv,loc,2}$, $\hat{\epsilon}_{iv,glo,2}$ computed as the $y - \hat{\mu}_1$ difference.

Hence, to ensure an accurate estimation of internal variability in the constraint procedure, an iterative algorithm is applied to find the MAR parameters that fit the residuals from the constrained forced response

$$ \begin{aligned} y - \mu &\xrightarrow{\text{residuals}, \hat{\epsilon}_{iv,1}} \text{constraint} \xrightarrow{\hat{\mu}_1} \\ y - \hat{\mu}_1 &\xrightarrow{\text{residuals}, \hat{\epsilon}_{iv,2}} \text{constraint} \xrightarrow{\hat{\mu}_2} \\ y - \hat{\mu}_{n-1} &\xrightarrow{\text{residuals}, \hat{\epsilon}_{iv,n}} \text{constraint} \xrightarrow{\hat{\mu}_n} \end{aligned} \quad (10) $$

where, for each iteration $n$, $\hat{\mu}_n$ and $\hat{\epsilon}_{iv,n}$ are estimates of the forced response and internal variability, respectively. The termination criterion is based on the Frobenius norm $\| \cdot \|_F$. Hence, we consider that the algorithm converges at the iteration $n$, i.e., that $\hat{\epsilon}_{iv,n} \rightarrow \epsilon_{iv}$ when the relative difference between $\|\hat{\Sigma}_{iv,n}\|_F$ and $\|\Sigma_{iv,n-1}\|_F$ is inferior to 1%, meaning that the MAR parameters values have also converged. In practice, $n$ varies between 2 and 4 depending on the location.

The autocorrelations from this MAR model suggest that our statistical representation of internal variability effectively captures decadal variability (typically between lag 5 and lag 10) in the GMST and local temperature time series, e.g., for the Atlantic, African, and South American regions (Fig. 7). We are aware that initial condition large ensembles and long piControl simulations provide a nice sampling of internal variability and could also be used to estimate this variability. However, we choose to not directly rely on it because of the huge discrepancies between models in terms of their simulated internal variability (56). Figures S20 to S30 illustrate this aspect with the piControl simulations from the CMIP6 models, including those used to build large ensembles. In all cases, the models do not converge to a consistent estimate of internal variability. For instance, over the Atlantic ocean, many models exhibit clear pseudo-periodic low frequency variability, while other models do not simulate decadal variability.

## Modeling of the dependence between local and global internal variability

As impacts from Pacific and Atlantic decadal variability (and potential other modes of variability) on GMST have been reported over the historical period (59, 60), we need to allow a potential dependence between global and local internal variability in $\Sigma_{iv,dep}$. Therefore, finding a simple and parsimonious dependence model that is compatible with the MAR structure is required. Allowing the covariances $\text{Cov}[\epsilon_{s,glo}(t), \epsilon_{s,loc}(t)]$ and $\text{Cov}[\epsilon_{f,glo}(t), \epsilon_{f,loc}(t)]$ to be nonzero is not trivial, and these terms need to be quantified to fill the covariance matrix $\Sigma_{iv,dep}$. Note that the fast and slow components remain always independent and that $\Sigma_{iv}$ is computed for each location separately, as the spatial dependence among various locations is not considered in the method. To compute $\Sigma_{iv,dep}$, we introduce a ninth parameter $\lambda$ accounting for some correlation between the local versus global components in the MAR modeling. The formulation of the covariances is slightly different in this case, and the calculations are detailed in the Supplementary Materials.

## Calculation of $p(x|y)$

As $\Pi(x)$ and $\epsilon$ are assumed to follow normal distributions, the Gaussian conditioning theorem is applicable to derive the posterior or the
<figure><img src="image_14.png" /><figcaption>**Fig. 7.** Estimation of observed internal variability. (A) Autocorrelation function (ACF) of the GMST observations (bars) after applying the iterative algorithm (see Materials and Methods) to the CMIP6 multimodel mean forced response and the GMST observed time series. The ACF of the MAR processes fitted to the same residual time series according to Eq. 9 is in red. (B to L) Same as (A) but for local residuals from 11 selected points (see Fig. 5).</figcaption></figure>

"constrained" distribution $p(x|y)$. Its formulation (detailed in eq. S23) indicates that the method is conservative: The uncertainty in $p(x|y)$ is never larger than that in $\Pi(x)$. Therefore, if observed internal variability is very large, then the model uncertainty in $p(x|y)$ will remain very close to that in $\Pi(x)$.

## Perfect model evaluation

We evaluate the performance of the KCC method within a perfect model framework, following a leave-one-out cross-validation:

1. For a given model, we consider a single member as pseudo-observations $y$ over the 1850–2021 period (the historical simulation is extended by the SSP5-8.5 simulation over the 2015–2021 period).

2. We use the other 26 models to derive the prior $\Pi(x) \sim N(\mu, \Sigma_{\text{mod}})$.

3. As there is no measurement uncertainty in models, $\Sigma_{\text{meas.}}$ is null; therefore, $\Sigma_{\text{obs}} = \Sigma_{\text{iv}}$. As done with the real observations, internal variability within the pseudo-observations is estimated from the difference between the pseudo-observations time series and the forced temperature response estimated by the ensemble mean of the 26 other models. $\Sigma_{\text{iv}}$ is then derived from the MAR fitted on the obtained residuals.

4. We apply the KCC method using the inputs $y, \Sigma_{\text{obs}}, \mu, \Sigma_{\text{mod}}$ to calculate projected changes constrained by pseudo-observations.

5. These four steps are repeated for each available member of the considered model and for all available models.
Continuous ranked probability score

We use the CRPS (31) to quantify the performance of the KCC method. It is defined as the quadratic measure of discrepancy between (i) $1(x \ge y_{pobs})$, the empirical cumulative distribution function (CDF) of a scalar pseudo-observation $y_{pobs}$ simulated by one model and averaged over the future period, and (ii) the projected CDF $G_{cons}$ of $p(x|y)$ (derived from all of the other models) over the same period

$$CRPS_{\text{cons}}(G_{\text{cons}}, y_{\text{pobs}}) = \int_{\mathbb{R}} [G_{\text{cons}}(x) - 1(x \ge y_{\text{pobs}})]^2 dx \quad (11)$$

where 1 is the indicator function (note that $x$ is here a bound variable in the integral, different from the vector $\mathbf{x}$ in Eq. 4). Similarly, we define a reference CRPS, $CRPS_{ref}$ based on $G_{ref}$, the CDF of $\Pi(\mathbf{x})$, the unconstrained distribution, and $y_{pobs}$. We can compute the CRPSS, which quantifies the performance of the KCC method if compared to the reference

$$\text{CRPSS} = 1 - \frac{\text{CRPS}_{\text{cons}}}{\text{CRPS}_{\text{ref}}} \quad (12)$$

The CRPSS is computed over all available pseudo-observations (121 values; see Table 1). $CRPS_{cons}$ is calculated in both GMST-only and Local + GMST cases. Therefore, the quantity $1 - \frac{CRPS_{\text{cons}}(\text{Local} + \text{GMST})}{CRPS_{\text{cons}}(\text{Global} - \text{only})}$ allows quantifying the added value from local observations compared to the sole use of GMST observations. A positive (negative) value, indicates an improvement (deterioration). The higher the CRPSS (bounded at 1), the better the performance.

SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at https://science.org/doi/10.1126/sciadv.abo6872

REFERENCES AND NOTES

1. Intergovernmental Panel on Climate Change, *Climate Change 2021: The Physical Science Basis*. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change (Cambridge Univ. Press, 2021).

2. E. Hawkins, R. Sutton, Time of emergence of climate signals. *Geophys. Res. Lett.* **39**, L01702 (2012).

3. E. Hawkins, D. Frame, L. Harrington, M. Joshi, A. King, M. Rojas, R. Sutton, Observed emergence of the climate change signal: From the familiar to the unknown. *Geophys. Res. Lett.* **47**, e2019GL086259 (2020).

4. C. Deser, R. Knutti, S. Solomon, A. S. Phillips, Communication of the role of natural variability in future North American climate. *Nat. Clim. Change* **2**, 775–779 (2012).

5. I. Mahlstein, R. Knutti, S. Solomon, R. W. Portmann, Early onset of significant local warming in low latitude countries. *Environ. Res. Lett.* **6**, 034009 (2011).

6. Intergovernmental Panel on Climate Change, *Climate Change 2013: The Physical Science Basis*. IPCC Working Group I Contribution to AR5 (Cambridge Univ. Press, 2013).

7. J. Lee, et al., *Climate Change 2021: The Physical Science Basis*. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change (Cambridge Univ. Press, 2021).

8. V. Eyring, S. Bony, G. A. Meehl, C. A. Senior, B. Stevens, R. J. Stouffer, K. E. Taylor, Overview of the Coupled Model Intercomparison Project Phase 6 (CMIP6) experimental design and organization. *Geosci. Model Dev.* **9**, 1937–1958 (2016).

9. B. C. O'Neill, C. Tebaldi, D. P. van Vuuren, V. Eyring, P. Friedlingstein, G. Hurtt, R. Knutti, E. Kriegler, J.-F. Lamarque, J. Lowe, G. A. Meehl, R. Moss, K. Riahi, B. M. Sanderson, The Scenario Model Intercomparison Project (ScenarioMIP) for CMIP6. *Geosci. Model Dev.* **9**, 3461–3482 (2016).

10. K. B. Tokarska, M. B. Stolpe, S. Sippel, E. M. Fischer, C. J. Smith, F. Lehner, R. Knutti, Past warming trend constrains future warming in CMIP6 models. *Sci. Adv.* **6**, eaaz9549 (2020).

11. A. Ribes, S. Qasmi, N. P. Gillett, Making climate projections conditional on historical observations. *Sci. Adv.* **7**, eabc0671 (2021).

12. L. Brunner, A. G. Pendergrass, F. Lehner, A. L. Merrifield, R. Lorenz, R. Knutti, Reduced global warming from CMIP6 projections when weighting models by performance and independence. *Earth Syst. Dynam.* **11**, 995–1012 (2020).

13. F. J. M. M. Nijssse, P. M. Cox, M. S. Williamson, Emergent constraints on transient climate response (TCR) and equilibrium climate sensitivity (ECS) from historical warming in CMIP5 and CMIP6 models. *Earth Syst. Dynam.* **11**, 737–750 (2020).

14. M. Schlund, A. Lauer, P. Gentine, S. C. Sherwood, V. Eyring, Emergent constraints on equilibrium climate sensitivity in CMIP5: Do they hold for CMIP6? *Earth Syst. Dynam.* **11**, 1233–1258 (2020).

15. C. B. Field, V. Barros, D. J. Dokken, Managing the Risks of Extreme Events and Disasters to Advance Climate Change Adaption: Special report of the Intergovernmental Panel on Climate Change (Cambridge Univ. Press, 2012).

16. L. Brunner, C. McSweeney, A. P. Ballinger, D. J. Befort, M. Benassi, B. Booth, E. Coppola, H. de Vries, G. Harris, G. C. Hegerl, R. Knutti, G. Lenderink, J. Lowe, R. Nogherotto, C. O'Reilly, S. Qasmi, A. Ribes, P. Stocchi, S. Undorf, Comparing methods to constrain future European climate projections using a consistent framework. *J. Climate* **33**, 8671–8692 (2020).

17. L. Brunner, R. Lorenz, M. Zumwald, R. Knutti, Quantifying uncertainty in European climate projections using combined performance-independence weighting. *Environ. Res. Lett.* **14**, 124010 (2019).

18. G. Abramowitz, C. H. Bishop, Climate model dependence and the ensemble dependence transformation of CMIP projections. *J. Climate* **28**, 2332–2348 (2015).

19. C. H. Bishop, G. Abramowitz, Climate model dependence and the replicate Earth paradigm. *Climate Dynam.* **41**, 885–900 (2013).

20. A. Borodina, E. M. Fischer, R. Knutti, Emergent constraints in climate projections: A case study of changes in high-latitude temperature variability. *J. Climate* **30**, 3655–3670 (2017).

21. F. Lehner, A. W. Wood, J. A. Vano, D. M. Lawrence, M. P. Clark, J. S. Mankin, The potential to reduce uncertainty in regional runoff projections from climate models. *Nat. Clim. Chang.* **9**, 926–933 (2019).

22. R. Lorenz, N. Herger, V. Sedláček, V. Eyring, E. M. Fischer, R. Knutti, Prospects and caveats of weighting climate models for summer maximum temperature projections over North America. *J. Geophys. Res. Atmos.* **123**, 4509–4526 (2018).

23. R. Sutton, E. Suckling, E. Hawkins, What does global mean temperature tell us about local climate? *Philos. Trans. R. Soc. A.* **373**, 20140426 (2015).

24. C. Tebaldi, J. M. Arblaster, Pattern scaling: Its strengths and limitations, and an update on the latest model simulations. *Clim. Change* **122**, 459–471 (2014).

25. D. Viner, M. Ekstrom, M. Hulbert, N. K. Warner, A. Wreford, Z. Zommers, Understanding the dynamic nature of risk in climate change assessments—A new starting point for discussion. *Atmos. Sci. Lett.* **21**, e958 (2020).

26. A. Hall, X. Qu, Using the current seasonal cycle to constrain snow albedo feedback in future climate change. *Geophys. Res. Lett.* **33**, L03502 (2006).

27. A. Hall, P. Cox, C. Huntingford, S. Klein, Progressing emergent constraints on future climate change. *Nat. Clim. Chang.* **9**, 269–278 (2019).

28. K. Cowtan, Z. Hausfather, E. Hawkins, P. Jacobs, M. E. Mann, S. K. Miller, B. A. Steinman, M. B. Stolpe, R. G. Way, Robust comparison of climate models with observations using blended land air and ocean sea surface temperatures. *Geophys. Res. Lett.* **42**, 6526–6534 (2015).

29. M. Richardson, K. Cowtan, R. J. Millar, Global temperature definition affects achievement of long-term climate goals. *Environ. Res. Lett.* **13**, 054004 (2018).

30. J. A. Screen, I. Simmonds, The central role of diminishing sea ice in recent Arctic temperature amplification. *Nature* **464**, 1334–1337 (2010).

31. T. Gneiting, A. E. Raftery, A. H. Westveld, T. Goldman, Calibrated probabilistic forecasting using ensemble model output statistics and minimum CRPS estimation. *Mon. Weather Rev.* **133**, 1098–1118 (2005).

32. A. Voldoire, D. Saint-Martin, S. Sénési, B. Decharme, A. Alias, M. Chevallier, J. Colin, J. F. Guérémy, M. Michou, M. P. Moine, P. Nabat, R. Roehrig, D. Salas y Mélia, R. Séférían, S. Valcke, I. Beau, S. Belamari, S. Berthet, C. Cassou, J. Cattiaux, J. Deshayes, H. Douville, C. Ethé, L. Franchistéguy, O. Geoffroy, C. Lévy, G. Madec, Y. Meurdesoiff, R. Msadek, A. Ribes, E. Sanchez-Gomez, L. Terray, R. Waldman, Evaluation of CMIP6 DECK Experiments With CNRM-CM6-1. *J. Adv. Modeling Earth Syst.* **11**, 2177–2213 (2019).

33. O. Boucher, J. Servonnat, A. L. Albright, O. Aumont, Y. Balkanski, V. Bastrikov, S. Bekki, R. Bonnet, S. Bony, L. Bopp, P. Braconnot, P. Brockmann, P. Cadule, A. Caubel, F. Cheruy, F. Codron, A. Cozic, D. Cugnet, F. D'Andrea, P. Davini, C. de Lavergne, S. Denvil, J. Deshayes, M. Devilliers, A. Ducharne, J.-L. Dufresne, E. Dupont, C. Ethé, L. Fairhead, L. Falletti, S. Flavoni, M.-A. Foujols, S. Gardoll, G. Gastineau, J. Ghattas, J.-Y. Grandpeix, B. Guenet, E. G. Lionel, E. Guilyardi, M. Guimberteau, D. Hauglustaine, F. Hourdin, A. Idelkadi, S. Joussaume, M. Kageyama, M. Khodri, G. Krinner, N. Lebas, G. Levavasseur, C. Lévy, L. Li, F. Lott, T. Lurton, S. Luyssaert, G. Madec, J.-B. Madeleine, F. Maignan, M. Marchand, O. Marti, L. Mellul, Y. Meurdesoiff, J. Mignot, I. Musat, C. Ottlé, P. Peylin, Y. Planton, J. Polcher, C. Rio, N. Rochetin, C. Rousset, P. Sepulchre, A. Sima, D. Swingedouw, R. Thiéblemont, A. K. Traore, M. Vancoppenolle, J. Vial, J. Vialard, N. Viovy, N. Vuichard, Presentation and evaluation of the IPSL-CM6A-LR climate model. *J. Adv. Modeling Earth Syst.* **12**, e2019MS002010 (2020).
34. M. E. Mann, B. A. Steinman, D. J. Brouillette, S. K. Miller, Multidecadal climate oscillations during the past millennium driven by volcanic forcing. *Science* **371**, 1014–1019 (2021).

35. K. Haustein, F. E. L. Otto, V. Venema, P. Jacobs, K. Cowtan, Z. Hausfather, R. G. Way, B. White, A. Subramanian, A. P. Schurer, A limited role for unforced internal variability in twentieth-century warming. *J. Climate* **32**, 4893–4917 (2019).

36. K. Wyser, T. van Noije, S. Yang, J. von Hardenberg, D. O'Donnell, R. Döscher, On the increased climate sensitivity in the EC-Earth model from CMIP5 to CMIP6. *Geosci. Model Dev.* **13**, 3465–3474 (2020).

37. C.-F. Schleussner, T. K. Lissner, E. M. Fischer, J. Wohland, M. Perrette, A. Golly, J. Rogelj, K. Childers, J. Schewe, K. Frieler, M. Mengel, W. Hare, M. Schaeffer, Differential climate impacts for policy-relevant limits to global warming: The case of 1.5 °C and 2 °C. *Earth Syst. Dynam.* **7**, 327–351 (2016).

38. K. E. Taylor, R. J. Stouffer, G. A. Meehl, An overview of CMIP5 and the experiment design. *Bull. Am. Meteorol. Soc.* **93**, 485–498 (2012).

39. C. Tebaldi, R. Knutti, The use of the multi-model ensemble in probabilistic climate projections. *Philos. Trans. R. Soc. A.* **365**, 2053 (2007).

40. B. M. Sanderson, R. Knutti, P. Caldwell, A representative democracy to reduce interdependency in a multimodel ensemble. *J. Climate* **28**, 5171–5194 (2015).

41. V. Eyring, P. M. Cox, G. M. Flato, P. J. Gleckler, G. Abramowitz, P. Caldwell, W. D. Collins, B. K. Gier, A. D. Hall, F. M. Hoffman, G. C. Hurtt, A. Jahn, C. D. Jones, S. A. Klein, J. P. Krasting, L. Kwiatkowski, R. Lorenz, E. Maloney, G. A. Meehl, A. G. Pendergrass, R. Pincus, A. C. Ruane, J. L. Russell, B. M. Sanderson, B. D. Santer, S. C. Sherwood, I. R. Simpson, R. J. Stouffer, M. S. Williamson, Taking climate model evaluation to the next level. *Nat. Clim. Change* **9**, 102–110 (2019).

42. R. Knutti, R. Furrer, C. Tebaldi, J. Cermak, G. A. Meehl, Challenges in combining projections from multiple climate models. *J. Climate* **23**, 2739–2758 (2010).

43. C. J. Smith, R. J. Kramer, G. Myhre, K. Alterskjær, W. Collins, A. Sima, O. Boucher, J.-L. Dufresne, P. Nabat, M. Michou, S. Yukimoto, J. Cole, D. Paynter, H. Shiogama, F. M. O'Connor, E. Robertson, A. Wiltshire, T. Andrews, C. Hannay, R. Miller, L. Nazarenko, A. Kirkevåg, D. Olivie, S. Fiedler, A. Lewinschal, C. Mackallah, M. Dix, R. Pincus, P. M. Forster, Effective radiative forcing and adjustments in CMIP6 models. *Atmos. Chem. Phys.* **20**, 9591 (2020).

44. R. Knutti, D. Masson, A. Gettelman, Climate model genealogy: Generation CMIP5 and how we got there. *Geophys. Res. Lett.* **40**, 1194–1199 (2013).

45. B. M. Sanderson, R. Knutti, On the interpretation of constrained climate model ensembles. *Geophys. Res. Lett.* **39**, L16708 (2012).

46. R. Knutti, J. Sedláček, B. M. Sanderson, R. Lorenz, E. M. Fischer, V. Eyring, A climate model projection weighting scheme accounting for performance and interdependence. *Geophys. Res. Lett.* **44**, 1909–1918 (2017).

47. B. M. Sanderson, M. Wehner, R. Knutti, Skill and independence weighting for multi-model assessments. *Geosci. Model Dev.* **10**, 2379–2395 (2017).

48. C. P. Morice, J. J. Kennedy, N. A. Rayner, J. P. Winn, E. Hogan, R. E. Killick, R. J. H. Dunn, T. J. Osborn, P. D. Jones, I. R. Simpson, An updated assessment of near-surface temperature change from 1850: The HadCRUT5 Data Set. *J. Geophys. Res. Atmos.* **126**, e2019JD032361 (2021).

49. N. P. Gillett, H. Shiogama, B. Funke, G. Hegerl, R. Knutti, K. Matthes, B. D. Santer, D. Stone, C. Tebaldi, The Detection and Attribution Model Intercomparison Project (DAMIP v1.0) contribution to CMIP6. *Geosci. Model Dev.* **9**, 3685–3697 (2016).

50. J. D. Annan, J. C. Hargreaves, Reliability of the CMIP3 ensemble. *Geophys. Res. Lett.* **37**, L02703 (2010).

51. A. Ribes, F. W. Zwiers, J.-M. Azaïs, P. Naveau, A new statistical approach to climate change detection and attribution. *Climate Dynam.* **48**, 367–386 (2017).

52. J. Rougier, M. Goldstein, L. House, Second-order exchangeability analysis for multimodel ensembles. *J. Am. Stat. Assoc.* **108**, 852–863 (2013).

53. G. J. van Oldenborgh, F. J. D. Reyes, S. S. Drijfhout, E. Hawkins, Reliability of regional climate model trends. *Environ. Res. Lett.* **8**, 014055 (2013).

54. O. Geoffroy, D. Saint-Martin, D. J. L. Olivie, A. Voldoire, G. Bellon, S. Tytéca, Transient climate response in a two-layer energy-balance model. Part I: Analytical solution and parameter calibration using CMIP5 AOGCM experiments. *J. Climate* **26**, 1841–1857 (2013).

55. N. J. Leach, S. Jenkins, Z. Nicholls, C. J. Smith, J. Lynch, M. Cain, T. Walsh, B. Wu, J. Tsutsui, M. R. Allen, FalRv2.0.0: A generalized impulse response model for climate uncertainty and future scenario exploration. *Geosci. Model Dev.* **14**, 3007–3036 (2021).

56. L. A. Parsons, M. K. Brennan, R. C. J. Wills, C. Proistosescu, Magnitudes and spatial patterns of interdecadal temperature variability in CMIP6. *Geophys. Res. Lett.* **47**, e2019GL086588 (2020).

57. S. Qasmi, C. Cassou, J. Boé, Teleconnection between Atlantic multidecadal variability and European temperature: Diversity and evaluation of the coupled model intercomparison project phase 5 models. *Geophys. Res. Lett.* **44**, 11,40–11,49 (2017).

58. T. Laepple, P. Huybers, Ocean surface temperature variability: Large model–data differences at decadal and longer periods. *Proc. Natl. Acad. Sci.* **111**, 16682–16687 (2014).

59. K. E. Trenberth, J. T. Fasullo, An apparent hiatus in global warming? *Earth's Future* **1**, 19–32 (2013).

60. H. Douville, A. Voldoire, O. Geoffroy, The recent global warming hiatus: What is the role of Pacific variability? *Geophys. Res. Lett.* **42**, 880–888 (2015).

61. C. J. Smith, G. R. Harris, M. D. Palmer, N. Bellouin, W. Collins, G. Myhre, M. Schulz, J.-C. Golaz, M. Ringer, T. Storelvmo, P. M. Forster, Energy budget constraints on the time history of aerosol forcing and climate sensitivity. *J. Geophys. Res. Atmos.* **126**, e2020JD033622 (2021).

**Acknowledgments:** We thank H. Douville for discussions about this work. We thank the climate modeling groups involved in CMIP6 exercises for producing and making available their simulations. We thank K. Cowtan for providing the blended temperature data and the ETH-Zürich for providing CMIP6 data through their cmip6-ng interface (http://dx.doi.org/10.5281/zenodo.3734128). The analyses and figures were produced with the R software (https://R-project.org/) and the NCAR Command Language software (http://dx.doi.org/10.5065/D6WD3XHS). **Funding:** This work was supported by the European Union's Horizon 2020 Research and Innovation Programme in the framework of the EUCP project (grant agreement no. 3776613), CONSTRAIN project (grant agreement no. 820829), and Météo-France. **Author contributions:** S.Q. and A.R. contributed to the design of the study, the interpretation of the results, and the writing of the manuscript. S.Q. processed the data and produced the figures. **Competing interests:** The authors declare that they have no competing interests. **Data and materials availability:** The datasets generated during and/or analyzed during the current study are available at https://doi.org/10.5281/zenodo.7037597. The HadCRUT5 dataset is available at https://metoffice.gov.uk/hadobs/hadcrut5/index.html. The CMIP6 datasets are available at https://esgf-node.llnl.gov/projects/esgf-llnl/. Effective radiative forcing time series corresponding to the CMIP6 exercise (61) are available at https://doi.org/10.5281/zenodo.3515338. All required programs to run the statistical method are in the associated KCC R package, which is available under a GNU General Public License, version 3 (GPLv3), at https://doi.org/10.5281/zenodo.5233947. All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials.

Submitted 18 February 2022

Accepted 26 August 2022

Published 12 October 2022

10.1126/sciadv.abo6872
# Reducing uncertainty in local temperature projections

Sad Qasmi and Aurlien Ribes

Sci. Adv., 8 (41), eabo6872.

DOI: 10.1126/sciadv.abo6872

## View the article online

https://www.science.org/doi/10.1126/sciadv.abo6872

## Permissions

https://www.science.org/help/reprints-and-permissions