Testing the assumptions in emergent constraints: why does the “emergent constraint on equilibrium climate sensitivity from global temperature variability” work for CMIP5 and not CMIP6?

It has been shown that a theoretically derived relation between annual global mean temperature variability and climate sensitivity held in the CMIP5 climate model ensemble (Cox et al., 2018a, hereafter CHW18). This so-called emergent relationship was then used with observations to constrain the value of equilibrium climate sensitivity (ECS) to about 3 °C. Since this study was published, CMIP6, a newer ensemble of climate models has become available. Schlund et al. (2020) showed that many of the emergent constraints found in CMIP5 were much weaker in the newer ensemble, including that of CHW18. As the constraint in CHW18 was based on a relationship derived from reasonable physical principles, it is of interest to find out why it is weaker in CMIP6. Here, we look in detail at the assumptions made in deriving the emergent relationship in CHW18 and test them for CMIP5 and CMIP6 models. We show one assumption, that of low correlation and variation between ECS and the internal variability parameter, a parameter that captures chaotic internal variability and sub-annual (fast) feedbacks, that while true for CMIP5 is not true for CMIP6. When accounted for, an emergent relationship appears once again in both CMIP ensembles, implying the theoretical basis is still applicable while the original assumption in CHW18 is not. Unfortunately, however, we are unable to provide an emergent constraint in CMIP6 as observational estimates of the internal variability parameter are too uncertain.

1 Introduction

Since the first general circulation climate models were introduced in the 1960s (Manabe and Bryan, 1969; Manabe and Wetherald, 1975), an ever-increasing amount of effort has been spent developing and improving these models to produce simulations that are increasingly more realistic and feature more of the processes and interactions present in the real world. The progress and understanding of the processes governing the Earth's climate as a result has been impressive. However, even after decades of research, the range of predictions of some key characteristics of the Earth's future climate coming from these models are actually increasing rather than narrowing with time, one particular characteristic being the amount of warming due to doubling of CO₂ at equilibrium, known as equilibrium climate sensitivity (ECS, Sherwood et al., 2020). Even though the latest state-of-the-art climate models in the Coupled Model Intercomparison Project 6 (CMIP6, Eyring et al., 2016; CMIP6_database, 2021) have a larger range of ECS values ([1.84; 5.68 K]) than previous CMIP model ensembles (Forster et al., 2020), the latest IPCC estimates have actually narrowed. For decades the IPCC “likely” range for ECS was between 1.5 and 4.5 K. In the latest report (IPCC, 2021) this was reduced to between 2.5 and 4 K with a best estimate of 3 K.

There have been numerous attempts (Knutti et al., 2017) to constrain ECS using the historical warming record, paleoclimate data and climate model experiments. Researchers have also used the emergent constraint technique (Hall et al., 2019; Brient, 2020; Williamson et al., 2021) to constrain ECS (see, for example, Covey et al., 2000; Knutti et al., 2006; Masson and Knutti, 2011a; Hargreaves et al., 2012; Sherwood et al., 2014; Caldwell et al., 2018, and the many references listed in Williamson et al., 2021). The basic idea of emergent constraints is to identify an observable of the climate x that varies significantly across a climate model ensemble and that exhibits a statistically significant relationship f(x) with another variable y describing an aspect of the climate model's future state. The relationship $y=f\left(x\right)+\mathit{\epsilon }$, is referred to as an “emergent relationship” where ε is a relatively small departure from f. Since x is observable, it can be measured in the real world. As such, f may then place a useful constraint on y, provided that the measurement uncertainty in x is small compared to the range of simulated values. This constraint is “emergent” because the emergent relationship f cannot be diagnosed from a single climate model. It becomes apparent only when the full ensemble is analysed.

There are pitfalls with the emergent constraint approach that must be guarded against particularly when the emergent relationships are not founded on well-understood physical processes. For example, data-mining outputs from climate models could lead to spurious correlations (Caldwell et al., 2014) and less than robust constraints on future changes (Bracegirdle and Stephenson, 2013). Care is also needed drawing statistical inferences from ensembles of small numbers of models. The problem is compounded if models within the ensemble share common components giving a smaller effective ensemble size (Pennell and Reichler, 2010; Masson and Knutti, 2011b; Herger et al., 2018). Observations used to guide model development also may lead to dependencies (Masson and Knutti, 2012) and common structural inaccuracies (Sanderson et al., 2021).

One way of guarding against spurious correlations between x and y is to use analytical solutions of simplified models of the full-complexity climate models to predict the emergent relationship f; f can then be tested against the results from the complex models. This approach was used in Cox et al. (2018a) (CHW18) where the analytical solution of the one-box or Hasselmann model (Hasselmann, 1976) provided an emergent relationship between the statistics of historical global annual mean temperature variability (x) and ECS (y; see Sect. 2 for further details). This emergent relationship was tested and found to hold in the CMIP5 (Taylor et al., 2011; CMIP5_database, 2021) models, although this was not without some debate regarding the applicability of the theory (Po-Chedley et al., 2018; Brown et al., 2018; Rypdal et al., 2018; Cox et al., 2018b; see Sect. 3 for a discussion of these points). However, since these works were published, the newer CMIP6 ensemble has become available. Schlund et al. (2020) showed that many of the emergent constraints found in CMIP5 were much weaker in the newer ensemble, including that of CHW18.

As the constraint in CHW18 was based on a relationship derived from reasonable physical principles, it is of interest to find out why it got weaker in CMIP6. Some possible reasons are as follows.

• The simple theory is not applicable to climate models and the real world. However, simple models (particularly two-box models) are regularly used to reproduce the annual global mean temperature response of climate models, and they do it well (see Caldeira and Myhrvold, 2013; Geoffroy et al., 2013b, a; Gregory, 2000; Held et al., 2010; MacMynowski et al., 2011).

• Estimates of the temperature variability observable (x) are uncertain enough to mask the relationship with ECS (y). This is unlikely as historical observations are long (>100 years) and relatively un-autocorrelated in time (a few years) giving good estimators of the true values.

• The assumptions made in deriving the emergent relationship that held for CMIP5 no longer hold for CMIP6. This is something we test in this paper.

The central interest of this paper is to test the assumptions that go into the derivation of the emergent relationship in CHW18. These assumptions are outlined in Sect. 3 and then tested in the CMIP5 and CMIP6 model ensembles with the aim of understanding why the emergent relationship in CHW18 is weaker for the CMIP6 model ensemble. Of course all assumptions will be ultimately wrong if perfect agreement is expected (the often used quote “all models are wrong” applies). However, “some models are useful”, and we look for agreement “for all practical purposes (FAPP)”, a term coined by John Bell (Bell, 1990). We will largely not be interested in the final step of obtaining the emergent constraint that results from combining the emergent relationship with observations for reasons we will outline later in the paper.

The structure of the rest of the paper is as follows. In Sect. 2 we review the methodology of CHW18 and how it is used in this study. In Sect. 3 we explicitly list, discuss and test the assumptions in CHW18 and show which assumption fails for the CMIP6 model ensemble. In Sect. 4 we show how to recover a robust emergent relationship in both CMIP5 and CMIP6 ensembles by including the forcing parameter in the predictor x. In Sect. 6 we perform a rigorous test of the emergent relationship theory by numerical simulation and show it does a reasonable job (FAPP) reproducing the results seen in each of the ensembles of the full-complexity CMIP climate models. We discuss and conclude in Sect. 7.

2 CHW18 methodology

The response of the global mean surface air temperature anomaly T(t) with time t to forcing Q(t) is assumed to be well modelled by the one-box or Hasselmann model (Hasselmann, 1976, hereafter H76) in CHW18. Forcing in this model comes from random, short-timescale weather noise and other external sources such as solar radiation and changes in greenhouse gas concentrations. Air temperature sensitivity to forcing is parameterised by λ, a term that lumps all the effects of the Earth system's feedbacks together. The single box has heat capacity C. In this model, T(t) evolves according to

$$\begin{array}{}\text{(1)}& C{\displaystyle \frac{\mathrm{d}T}{\mathrm{d}t}}=Q\left(t\right)-\mathit{\lambda }T\left(t\right).\end{array}$$

Solving this model results in a linear relation between ECS and a metric of temperature variability Ψ, which is a form of fluctuation–dissipation theorem (Kubo, 1966; Leith, 1975). Explicitly

$$\begin{array}{}\text{(2)}& \mathrm{ECS}=\sqrt{\mathrm{2}}{\displaystyle \frac{Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}}{{\mathit{\sigma }}_Q}}\mathrm{\Psi }.\end{array}$$

Here $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ is the radiative forcing resulting from doubling the atmospheric CO₂ concentration and σ_Q is the standard deviation of a zero mean white-noise process designed to model the fast (sub-annual), chaotic weather forcing on the slower Earth system components. Ψ can be measured from temperature observations and is defined as

$$\begin{array}{}\text{(3)}& \mathrm{\Psi }={\displaystyle \frac{{\mathit{\sigma }}_T}{\sqrt{-\log {\mathit{\alpha }}_{\mathrm{1}T}}}},\end{array}$$

where σ_T is the standard deviation and α_1T is the autocorrelation at 1-year lag of annual global mean temperature. Details of this derivation can be found in CHW18 and Williamson et al. (2018).

CHW18 calculated the pair of values (Ψ_i, ECS_i) for each of the n=16 CMIP5 climate models labelled by $i\in \mathit{\{}\mathrm{1},\mathrm{2},\mathrm{\dots },n\mathit{\}}$ performing a simulation of the historical period 1880–2016. Plotting the n pairs confirmed the theoretically expected Ψ vs. ECS linear “emergent relationship” with good correlation (r=0.77, r in this paper denotes Pearson's correlation coefficient). Combining this resulting emergent relationship with Ψ from observational records of the same period gave an emergent constraint on ECS of 2.8±0.6 °C (plus minus values are 66 % confidence intervals).

Although there were more CMIP5 models available than the n=16 used in CHW18, the choice of one model per modelling centre was made to avoid biasing the emergent constraint towards similar models. Where multiple models were available from the same centre, the model with the lowest root-mean-square error to the observational temperature record was chosen. Po-Chedley et al. (2018) and Schlund et al. (2020) repeated the analysis of CHW18 including these additional models and thus had a larger CMIP5 ensemble (larger n). They found the emergent relationship got slightly weaker, although it was still “highly significant” in the language of Schlund et al. (2020).

In this paper we use the CHW18 methodology (further detailed in the original paper) and apply it to CMIP5 and CMIP6 models with the following differences. Here we look at the historical period 1880–2005 for both CMIP5 and CMIP6 ensembles following Schlund et al. (2020) rather than 1880–2016 as in CHW18. This is because the standard CMIP5 historical experiment ends in 2005. Increasing the time period to present day by concatenating with one of the CMIP rcp or ssp future projection experiments slightly increases the strength of the correlation in the emergent relationship. We also use a different ensemble of 15 CMIP5 models corresponding with those analysed in Geoffroy et al. (2013b). Geoffroy et al. (2013b) also lists FGOALS-s2; however, we leave this model out as it does not have a historical simulation with which to calculate Ψ. We use the Geoffroy et al. (2013b) ensemble as their published parameter values are used in Sect. 6 to run simulations of box models. These simulations are used to compare the theory with the full-complexity CMIP5 models. To make a fair comparison limits us to the same set. For the CMIP6 ensemble we use all models that have the necessary simulations for our analysis (piControl, historical and abrupt-4xCO2), a set of n=33 models. For both CMIP ensembles we use one run for each model, preferably the one labelled r1i1p1 (CMIP5) or r1i1p1f1 (CMIP6) where it exists. We look at the results for different runs of the same model in Sect. 5, although we find no qualitative changes to the findings with the r1i1p1 (CMIP5) or r1i1p1f1 choices. A list of models used and their parameter values is given in Appendix B.

Figure 1 Ψ vs. ECS emergent relationships in the CMIP5 (a) and CMIP6 (b) model ensembles running the historical experiment. The period 1880–2005 of each model's time series is used to calculate Ψ. ECS is determined from the abrupt4xCO2 experiment using the standard Gregory plot method. Individual models are plotted as circles (CMIP5 models are blue, while CMIP6 models are red). The best-fit line in the ordinary least-squares sense is shown in black along with the standard deviation of the prediction error (dotted black line). Pearson correlation r and p value are given for each emergent relationship in each subplot title.

The (Ψ, ECS) emergent relationships for CMIP5 and CMIP6 ensembles are shown in Fig. 1. The CMIP5 ensemble shows good correlation between Ψ and ECS, $r(\mathrm{\Psi },\mathrm{ECS})=\mathrm{0.66}$; however, for CMIP6 this is weaker, $r(\mathrm{\Psi },\mathrm{ECS})=\mathrm{0.31}$, confirming the results of CHW18 (although with slightly different historical period and set of CMIP5 models) and Schlund et al. (2020) (CMIP6).

Schlund et al. (2020) use the following definitions for significance based on p value. An emergent relationship is called “highly significant” if p<0.02, “barely significant” if $\mathrm{0.02}\le p<\mathrm{0.05}$, “almost significant” if $\mathrm{0.05}\le p<\mathrm{0.1}$ and “far from significant” if p≥0.1. We adopt their definitions in this paper. We find the (Ψ, ECS) emergent relationship highly significant for CMIP5 and almost significant for CMIP6.

3 Assumptions in CHW18

The following assumptions are made in the CHW18 methodology to obtain the emergent relationship between Ψ and ECS.

A1. The T(t) response to Q(t) is modelled well by H76 for timescales greater than 1 year and less than the detrending window length (55 years in CHW18).

A2. H76 is solved with a random, white-noise forcing Q(t) of zero mean and standard deviation σ_Q. This is designed to parameterise internally generated variability (from weather, for example, Hasselmann, 1976). It is assumed that the response from all other sources of forcing in the historical period such as (but not limited to) greenhouse gases (GHGs), solar irradiance and volcanoes can be removed via detrending to a good approximation so that Eq. (2) applies to this period in both observations and CMIP models.

A3. The forcing parameters, $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ and σ_Q in Eq. (2), are uncorrelated to ECS and their variation is small relative to the variation in Ψ. This requirement makes Ψ a good predictor of ECS.

There are further assumptions concerning the quantification of sources of uncertainty (structural, observational, etc.) in deriving the emergent constraint in CHW18. These are considered in more detail in Williamson and Sansom (2019) and Williamson et al. (2021). However, as we look only at the emergent relationship here, these will not discussed further.

3.1 Testing the assumptions

To summarise this subsection, assumption A3 is violated for the CMIP6 models. However, all other assumptions still apply FAPP for CMIP5 and CMIP6. In particular, it is the assumption of no correlation between ECS and the forcing parameter, σ_Q, that is no longer true for CMIP6. In CMIP6 significant correlation exists. Each assumption is discussed below in order.

Assumption A1 was studied in detail in Williamson et al. (2018) and Cox et al. (2018b). To summarise, H76 only really has any physical justification when the timescales of interest are dominated by the well-mixed atmosphere and ocean surface layer (a few years to decades). It is well known that H76 does a poor job of reproducing T(t) on longer timescales (e.g. Caldeira and Myhrvold, 2013; Schwartz, 2007, 2008; Foster et al., 2008; Kirk-Davidoff, 2009; Knutti et al., 2008; Scafetta, 2008). This led some to question the use of H76 in CHW18, e.g. Rypdal et al. (2018). However, one can show analytically (Williamson et al., 2018) that a near-linear emergent relationship is also expected between ECS and Ψ for the more realistic and widely used two-box (Gregory, 2000; Held et al., 2010) and diffusion models (MacMynowski et al., 2011). Both two-box and diffusion models are known to do a good job of reproducing the global annual mean temperature response of CMIP climate models (Caldeira and Myhrvold, 2013; Geoffroy et al., 2013b). As the T(t) solutions of CMIP6 models qualitatively have the same structural form as CMIP5 models to stepped and linearly increasing forcing (abrupt-4xCO2 and 1pctCO2 experiments, respectively), we expect that two-box and diffusion models also emulate the CMIP6 models well. We fit two-box models to the CMIP6 ensemble (as Geoffroy et al. (2013b) did for CMIP5) later in the paper and can confirm this is indeed the case. The reason the Ψ vs. ECS linear relationship still holds to a good degree in the more complete two-box and diffusion models is because Ψ is a statistic that is dominated by fast timescale processes of a few years, a feature H76 does capture well.

A2 assumes the response to all external forcing (GHGs, volcanoes, etc.) in the historical period can be removed to a good approximation by linearly detrending T(t) in a 55-year moving window, leaving just the internally generated random variability parameterised as the response to random “forcing” in H76.

This was the procedure introduced in CHW18, and we continue with the same procedure here for consistency and comparison. The reasons for using a 55-year window have been discussed in the original paper (Cox et al., 2018a) and subsequent publications (Cox et al., 2018b; Williamson et al., 2018). The reason for linear detrending is to remove the response due to the slow timescale in the climate. It turns out that when fitting two-box models to the CMIP models, i.e. a fast timescale (∼4 years) and a slow timescale (∼200 years) response result; see Geoffroy et al. (2013b), for example, or the tables in the Appendix B. Linear detrending with a 55-year timescale fits nicely between the short and fast timescale and removes the slow response component. It also minimises the uncertainty in the resulting emergent constraint (Cox et al., 2018a). Removing the slow timescale response leaves a signal that is more like the H76 (one-box) model and therefore more like the underlying simple theory of the emergent relationship.

Assumption A2 seeks to make the derivation of Eq. (2) (which is a derivation that applies to the piControl experiment) applicable to the historical simulations. Several works (Po-Chedley et al., 2018; Brown et al., 2018) showed this assumption to be false. In particular they showed that the detrending procedure in CHW18 does not remove the response to all external forcing. They also showed that better methods of removing forced variability slightly weakened the emergent relationship. Cox et al. (2018b) acknowledged this to be true; however, they also showed that external forcing, provided it is common for all models in the ensemble, would actually be helpful and improve the emergent relationship. This was demonstrated using an ensemble of H76 and two-box models tuned to mimic the CMIP5 models running a variety of experiments with and without common and random forcing. A sketch of the reason is as follows. Ψ is linearly proportional to sensitivity, and given an ensemble of models with a range of sensitivities, more sensitive models will respond with a larger Ψ (or response) if all models in the ensemble are given the same (common) forcing, providing a natural way of ordering the model's sensitivities. The common forcing in the historical simulations comes from volcanoes, anthropogenic trends, solar cycles, etc.

Equation (2) predicts a linear relationship between ECS and Ψ provided $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ and σ_Q can be treated as “constants” across the model ensemble (assumption A3). A looser definition of constant for Eq. (2) is stated in A3. In Fig. 2a we plot $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ against ECS and compute their correlation in both CMIP5 and CMIP6 ensembles. For both ensembles $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ is uncorrelated to ECS ($r=-\mathrm{0.17}$ for CMIP5, and $r=-\mathrm{0.07}$ for CMIP6; both p values, p≥0.1, are far from significant). $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ is determined in the standard way for each model running an abrupt-4xCO2 experiment via a Gregory plot (Gregory et al., 2004).

In Fig. 2b we plot the other forcing constant σ_Q against ECS. σ_Q is estimated from the detrended temperature residual of each climate model's historical run. The standard deviation of white noise forcing σ_Q is fitted for each model from the global annual mean temperature time series. This time series is linearly detrended with a rolling 55 year window. This is to isolate the T(t) response to internal variability, analogous to how Ψ is determined in the CHW18 methodology, to leave the noisy T(t) response to white noise with standard deviation σ_T. The theoretical formula is given by (see Williamson et al., 2018 for example)

$$\begin{array}{}\text{(4)}& {\mathit{\sigma }}_T^{\mathrm{2}}={\displaystyle \frac{{\mathit{\sigma }}_Q^{\mathrm{2}}}{\mathrm{2}\mathit{\lambda }C}}.\end{array}$$

We rearrange this relation to get σ_Q in terms of the observable σ_T and the parameters λ and C (given in Tables B1, B2, B3 and B4; see Sect. 6 for details on how the H76 model parameters are fitted). Values of σ_Q in both the historical and piControl runs are also reported in these tables.

Consistent with CHW18, σ_Q is uncorrelated to ECS in CMIP5 ($r=-\mathrm{0.09}$); however, in CMIP6 there is highly significant anti-correlation ($r=-\mathrm{0.58}$, p<0.001). We could equally estimate σ_Q from piControl simulations. We choose the historical experiment for consistency with estimation of Ψ. Whichever simulation is used, the correlation with ECS remains largely invariant (piControl $r({\mathit{\sigma }}_Q,\mathrm{ECS})=-\mathrm{0.09}$ and $r({\mathit{\sigma }}_Q,\mathrm{ECS})=-\mathrm{0.58}$ in CMIP5 and CMIP6, respectively).

Figure 2 Individual models are plotted as circles (CMIP5 models are blue, while CMIP6 models are red). The best-fit line in the ordinary least-squares sense is shown in black along with the standard deviation of the prediction error (dotted black line). The Pearson correlation r and p value are given for each emergent relationship in each subplot title. Shown are $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ against CMIP5 (a) and CMIP6 (b) model ensembles. $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ is inferred from the abrupt4xCO2 experiment using the standard Gregory plot method. σ_Q against ECS in the CMIP5 (c) and CMIP6 (d) model ensembles. σ_Q is calculated from the period 1880–2005 of each model's historical experiment time series.

When plotting the combination of constants, $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}/{\mathit{\sigma }}_Q$, multiplying Ψ in Eq. (2), (figure not shown), CMIP6 still has highly significant correlation between $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}/{\mathit{\sigma }}_Q$ and ECS (r=0.74, p<0.001). CMIP5 shows some anti-correlation, although it is far from significant ($r=-\mathrm{0.21}$, p=0.46).

4 Recovering an emergent relationship

We have confirmed that Ψ is a good predictor of ECS for CMIP5 models, albeit not for CMIP6 models. In Fig. 3, when σ_Q is included in the x-axis predictor variable, a good emergent relationship is recovered for both CMIP ensembles, and both have highly significant p values of p<0.001. One can also include $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ (although it is uncorrelated to ECS in both CMIP ensembles) in the predictor, i.e. ECS $\propto \frac{Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}}{{\mathit{\sigma }}_Q}\mathrm{\Psi }$, to get a similarly skilful emergent relationship (figure not shown). We restrict this to $\mathrm{\Psi }/{\mathit{\sigma }}_Q$ as minimal degrees of freedom are preferred.

Figure 3 $\mathrm{\Psi }/{\mathit{\sigma }}_Q$ against ECS in the CMIP5 (a) and CMIP6 (b) model ensembles running the historical experiment. The period 1880–2005 of each model's time series is used to calculate Ψ and σ_Q. Individual models are plotted as circles (CMIP5 models are blue, while CMIP6 models are red). The best-fit line in the ordinary least-squares sense is shown in black along with the standard deviation of the prediction error (dotted black line). Pearson correlation r and p value are given for each emergent relationship in each subplot title.

Where does the skill in predicting ECS using $\mathrm{\Psi }/{\mathit{\sigma }}_Q$ come from? In CMIP5, it came from Ψ (an observable). There is no skill in σ_Q (it is uncorrelated with ECS). In CMIP6 the converse is roughly correct. There is limited correlation with ECS from Ψ but good correlation from σ_Q, which, to our knowledge, is unfortunately not directly observable.

Theoretically, these findings should hold equally well in the piControl run, although the emergent relationships should have slightly weaker correlation for reasons outlined in Sect. 3.1 and Cox et al. (2018b). Again, we find this is roughly true (see Fig. 4). For the piControl experiments we analyse the longest common period simulated in the CMIP5 and CMIP6 ensembles, which is 200 years.

Figure 4 Individual models are plotted as circles (CMIP5 models are blue and CMIP6 models are red). The best-fit line in the ordinary least-squares sense is shown in black along with the standard deviation of the prediction error (dotted black line). Pearson correlation r and p value are given for each emergent relationship in each subplot title. The first 200 years of each model's time series is used to calculate Ψ and σ_Q running the piControl experiment (no external forcing). (a) Ψ against ECS in the CMIP5 (a, b) and CMIP6 (c, d) model ensembles. (b) $\mathrm{\Psi }/{\mathit{\sigma }}_Q$ against ECS in the CMIP5 (a, b) and CMIP6 (c, d) model ensembles.

H76 in Eq. (2) predicts ECS $=k\frac{Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}}{{\mathit{\sigma }}_Q}\mathrm{\Psi }$ with a constant of proportionality $k=\sqrt{\mathrm{2}}$. We investigate the empirical value of k for the full-complexity models next. As this relation should hold for all models running any experiment (provided you can remove the forced signal), we have plotted ECS against $\frac{Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}}{{\mathit{\sigma }}_Q}\mathrm{\Psi }$ for an ensemble composed of all CMIP5 and CMIP6 models running both piControl and historical experiments (a total of 96 data points) to determine the empirical k (Fig. 5a). While the proportionality holds with a high correlation value and significance (r=0.74, p<0.001), the empirical constant is $k\sim \mathrm{2}\sqrt{\mathrm{2}}$ rather than the $\sqrt{\mathrm{2}}$ predicted by H76, i.e. the theoretical prediction of ECS from H76 is lower than the full-complexity models suggest.

Williamson et al. (2018) showed the two-box and diffusion models also shared the linear ECS–Ψ proportionality of H76 FAPP, albeit with slightly different variables and constants. We have therefore also compared the more realistic two-box model theoretical predictions of k (Eq. 23 in Williamson et al., 2018) to the full-complexity models (Fig. 5b). Using the two-box values in Tables B5 and B6 brings the empirically determined $k\sim \mathrm{1.3}\sqrt{\mathrm{2}}$ closer to the two-box theoretical value ($k=\sqrt{\mathrm{2}}$) with similar high correlation and significance (r=0.76, p<0.001). The two-box model adds a second, longer timescale to H76, mimicking the full-complexity models more closely; however, the theoretical k is still slightly low. The lower prediction of k than the empirical results suggest this could be due to the full-complexity models having other timescales that the conceptual models do not. Although the conceptual models predict a linear ECS-Ψ proportionality also seen in the full-complexity CMIP models, they do not predict the constant of proportionality well. This is why the empirically determined k should be used to obtain an emergent constraint as in CHW18.

Figure 5 ECS against the theoretical predictor it is proportional to in (a) the one-box or H76 model and (b) the two-box model. Individual models are plotted as circles (CMIP5 models are blue, while CMIP6 models are red) running both historical and piControl experiments between 1880–2005 and the first 200 years of the simulation respectively. The best-fit line in the ordinary least-squares sense is shown in black along with the standard deviation of the prediction error (dotted black line). Pearson correlation r and p value are given for each emergent relationship in each subplot title. The empirically determined constant of proportionality between the x-axis variable and ECS, k, is given in each subplot title. The H76 and two-box theoretical values of $k=\sqrt{\mathrm{2}}$ are plotted as the dashed black line.

5 Robustness to choice of model run

For both CMIP ensembles we have used one run for each model, preferably the one labelled or r1i1p1 (CMIP5) or r1i1p1f1 (CMIP6) where it exists; however, we could have equally chosen any r*i*p* (CMIP5) or r*i*p*f* (CMIP6) for each model provided multiple runs of the same model exist. In this section we show the results for r1i1p1 or r1i1p1f1 are representative of a typical random run choice.

For models with multiple runs, we have drawn at random one run (r*i*p* or r*i*p*f* for CMIP5 and CMIP6 respectively) for each model and repeated the analysis in the previous sections multiple times. For the historical runs, in both CMIP5 and CMIP6 ensembles, many models do repeated runs, sometimes multiple times. For example, the CMIP6 model CanESM5 has the most runs, performing the historical experiment 50 times. There are 1.6×10¹² and 5.4×10²⁰ unique permutations for the same set of CMIP5 and CMIP6 models, respectively, performing the historical experiment. These numbers are clearly too large to search exhaustively. We have therefore drawn 1000 unique permutations for the historical experiment and repeated the analysis in this paper, i.e. calculated the Pearson correlation, r, for every one of these 1000 permutations. The results are shown in the upper half of Tables 1 (CMIP5) and 2 (CMIP6). We find that the results reported for r1i1p1 and r1i1p1f1 where they exist are fairly typical of a randomly chosen set of runs; i.e. they fall within 1 standard deviation of the mean value in the CMIP5 ensemble. The CMIP6 historical experiment r(Ψ,ECS) with r1i1p1f1 is slightly higher than would be expected (mean value for a randomly chosen permutation is $r=\mathrm{0.18}\pm \mathrm{0.11}$, far from significant, compared to r=0.31, almost significant). We have also listed the outermost values (min and max) found in the 1000 random run choice distribution for completeness. As one would expect with a large enough sample, there is a chance of finding non-representative correlations; i.e. for the CMIP5 historical ensemble there is a small possibility that you might find an $r(\mathrm{\Psi },\mathrm{ECS})=\mathrm{0.14}$ (far from significant) or even $r(\mathrm{\Psi },\mathrm{ECS})=\mathrm{0.82}$ (highly significant) in a random pick of runs. These are, however, not typical results.

For the piControl runs, there are 15 and 24 unique permutations for the same set of CMIP5 and CMIP6 models, respectively. The relatively low number of unique permutations is due to the low number of repeated runs for the piControl experiment. We have run the same analysis in this paper; i.e. we calculated the Pearson correlation, r, for every one of these permutations. The results are shown in the lower half of Tables 1 (CMIP5) and 2 (CMIP6). Because of the low number of unique permutations, the range of results is much narrower than the historical experiment. Again, we find that the results reported for r1i1p1 and r1i1p1f1 where they exist are fairly typical of a randomly chosen set of runs; i.e. they fall within 1 standard deviation of the mean value in the CMIP5 and CMIP6 ensembles with the exception of CMIP5 $r(\mathrm{\Psi },\mathrm{ECS})=\mathrm{0.64}$ (highly significant) for r1i1p1 compared to $r=\mathrm{0.61}\pm \mathrm{0.02}$ (highly significant) for a randomly chosen permutation. This is actually the highest value of r (max) found in that experiment.

Table 1 Details of the CMIP5 experiment showing the robustness of correlation to choice of run, including a comparison of Pearson correlation r for model ensembles using the r1i1p1 run (third column) with randomly chosen runs r*i*p*. There are 15 unique permutations of runs for the piControl experiment and 1.6×10¹² unique permutations for the historical experiment. We calculate r for all unique permutations for the piControl experiment and 1000 randomly chosen unique permutations for the historical experiment. We report the mean value of r with its standard deviation (±) in the fourth column and the most extreme values of r in the fifth and sixth columns.

Table 2 Details of the CMIP6 experiment showing the robustness of correlation to choice of run, including a comparison of Pearson correlation r for model ensemble using the r1i1p1f1 run if available (third column) with randomly chosen runs r*i*p*f*. There are 24 unique permutations of runs for the piControl experiment and 5.4×10²⁰ unique permutations for the historical experiment. We calculate r for all unique permutations for the piControl experiment and 1000 randomly chosen unique permutations for the historical experiment. We report the mean value of r with its standard deviation (±) in the fourth column and the most extreme values of r in the fifth and sixth columns.

6 Can theory simulate the CMIP model results?

In Sect. 4 we found that by including the forcing parameter σ_Q in the predictor for ECS an emergent relationship could be recovered for both CMIP5 and CMIP6 ensembles. These relationships are present in both the historical and piControl experiments, giving confidence in the underlying theoretical basis FAPP.

In this section we make a more demanding test of the theoretical basis by asking if theory alone can simulate the full-complexity CMIP model ensemble r(Ψ,ECS) and $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ results. To do this, we create a H76 model emulator of each of the $i\in \mathit{\{}\mathrm{1},\mathrm{2},\mathrm{\dots },n\mathit{\}}$ full-complexity CMIP5 and CMIP6 models used in the preceding figures. With the emulator H76 models we can build emulator H76 CMIP5 and CMIP6 ensembles and run analogous historical and piControl experiments with them. This will allow us to compare the results of the pure theory used in CHW18 with that of the full-complexity CMIP model ensembles. We also fit the more complete two-box model in addition for comparison (see Sect. A).

6.1 Methodology

The H76 model fitted to each of the full-complexity CMIP models is given by Eq. (1). Parameters are fitted from the full-complexity abrupt-4xCO2 CMIP model experiments: λ and $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ are determined from Gregory plots (Tables B1, B2), C is found using a modification of the methodology of Geoffroy et al. (2013b) (Tables B3 and B4). Geoffroy et al. (2013b) published parameter values for two-box models fitted to CMIP5 models. The two-box model is H76's well-mixed upper-ocean and atmosphere box extended by coupling it to a large heat capacity deep-ocean box (see Sect. A). This gives the two-box model a fast and a slow e-folding timescale of adjustment with typical values of ∼4 and ∼200 years when fitted to CMIP models (see Geoffroy et al., 2013b and Tables B5 and B6), which is known to do a good job reproducing the global annual mean temperature response of climate models.

As H76 only has one box and therefore one timescale, it cannot capture both fast and slow responses of CMIP models. Because Ψ is a statistic that is dominated by fast timescale processes of a few years, a feature H76 does capture well, we choose to fit the fast response with H76 using the fast timescale fitting methodology of Geoffroy et al. (2013b) (see Eq. 18 in that paper). When modified for H76, this equation becomes

$$\begin{array}{}\text{(5)}& \mathit{\tau }=-{\displaystyle \frac{t}{\log (\mathrm{1}-\frac{T\left(t\right)\mathit{\lambda }}{Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}})}}.\end{array}$$

We fit H76's timescale parameter $\mathit{\tau }=\frac{C}{\mathit{\lambda }}$ (and therefore the heat capacity C) by averaging over the first 5 years of the abrupt-4xCO2 experiment. We choose the average over 5 years rather than the first 10 in the two-box fits of Geoffroy et al. (2013b). This is because the H76 fit gets worse as the number of years in the average increases (in the root-mean-square error of the fit). For the two-box fits we use the methodology of Geoffroy et al. (2013b) unmodified (see Sect. A for complete details).

Fitted values of λ, τ, C and σ_Q are reported for CMIP5 and CMIP6 ensembles in Tables B3 and B4 respectively.

6.2 Emulator piControl experiments

We perform analogous piControl experiments with the H76 and two-box CMIP5 and CMIP6 ensembles by integrating each of the individual CMIP H76 (and two-box) emulators (Eqs. 1 and A1 respectively) numerically with forcing Q_i(t), a zero mean random variable with model specific standard deviation σ_Qi and Gaussian pdf. We write this as

$$\begin{array}{}\text{(6)}& Q_i\left(t\right)={\mathit{\sigma }}_{Qi}{\mathit{\eta }}_i\left(t\right),\end{array}$$

where η_i(t) is the Gaussian random variable with unit standard deviation. The equations are integrated with a time step of 0.1 years using the Euler–Maruyama method. The T_i(t) time series that result are then analysed in the same way as the full-complexity CMIP model time series to produce a pair of values (Ψ_i,ECS_i). For the full set of n two-box models in each CMIP emulator ensemble r(Ψ,ECS) and $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ are calculated. Because Q_i(t) is a random variable of finite length, repeating the same experiment results in slightly different values of $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ for each run due to the properties of statistical estimators (estimation converges as $\mathrm{1}/\sqrt{N}$ where N is the number of points in the time series). The same applies to different initial value runs in the full-complexity models due to the chaotic weather variability the random forcing captures in the H76 and two-box models. We therefore repeat each piControl emulator experiment 250 times and compare the distribution of emulator r(Ψ,ECS) values with the single full-complexity CMIP piControl experiment.

Results are shown in Fig. 6. Agreement between the H76 emulator CMIP r(Ψ,ECS) and the full-complexity CMIP ensembles is reasonable. Full-complexity CMIP5 ensemble r(Ψ,ECS) results (Fig. 6a, dotted blue line) fall in the upper end of the distribution of r(Ψ,ECS) H76 emulator values. Although full-complexity r(Ψ,ECS) CMIP6 results (in red) were shown to be lower in correlation (red dotted line), they can still be simulated reasonably well by the H76 emulator ensemble, falling like CMIP5, in upper end of simulated r(Ψ,ECS) values. In Fig. 6b where Ψ is now normalised by the mean amplitude of the random forcing σ_Q both CMIP5 and CMIP6 results are much more similar, with histograms of the H76 emulator ensembles and the full-complexity results having much more overlap, although simulated values are still on average slightly lower than the full-complexity ensembles.

Analogous figures simulated with two-box emulator CMIP ensembles are shown in Fig. A1. The two-box ensembles do an even better job of simulating the full-complexity CMIP results.

Figure 6 Probability of obtaining r(Ψ,ECS) (a) and $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ (b) in the H76 CMIP emulator ensembles performing a piControl simulation of 200 years. The CMIP5 emulator ensemble histogram is given in blue, and the equivalent data for CMIP6 are shown in red. The full-complexity CMIP ensemble results performing the same experiment are shown as vertical dotted lines.

6.3 Emulator historical experiments

The analogous historical experiments are performed in the same way to the piControl experiments but with a common external forcing component Q_i(t) in addition to the random forcing. This comes from GHGs, volcanoes, solar cycles and others. For this common external forcing component we use Meinshausen et al. (2011) reconstructed historical forcing (Q_IPCC(t)). Explicitly, this is written as

$$\begin{array}{}\text{(7)}& Q_i\left(t\right)=Q_{\mathrm{IPCC}}\left(t\right)+{\mathit{\sigma }}_{Qi}{\mathit{\eta }}_i\left(t\right),\end{array}$$

in the historical simulations. We integrate the H76 and two-box ensembles between the years 1765 and 2005 but calculate Ψ and r(Ψ,ECS) between 1880–2005 to correspond to the full-complexity model analysis. Results are shown in Fig. 7. As with the piControl experiments in Sect. 6.2, agreement between the pure-theory H76 emulator ensembles and the full-complexity ensembles is reasonable, giving confidence that the underlying theory used in CHW18 is good FAPP. The analogous figure simulated with two-box CMIP emulator ensembles does an even better job (Fig. A2).

Figure 7 Probability of obtaining r(Ψ,ECS) (a) and $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ (b) in the H76 CMIP emulator ensembles performing a historical simulation of the period 1880–2005. The CMIP5 emulator ensemble histogram is given in blue, and the equivalent data for CMIP6 are shown in red. The full-complexity CMIP ensemble results performing the same experiment are shown as vertical dotted lines.

6.4 Emulator experiments with constant σ_Q

We have shown that by taking into account model-specific σ_Q in the CHW18 theory we can both understand r(Ψ,ECS) correlation results and can recover good emergent relationships for both CMIP5 and CMIP6 in piControl and historical runs. In CHW18 and Cox et al. (2018b) it was assumed that σ_Q was constant for each model in the CMIP5 ensemble (assumption A3). We now test this assumption with the H76 and two-box CMIP ensembles. Instead of fitting σ_Q to each CMIP model we fix it to be a constant, σ_Q=0.25 W m⁻² following Cox et al. (2018b). This value was chosen (even though it is lower than the values given in the tables) as it was the mean value of the standard deviation of net top-of-the-atmosphere radiation, which was thought to be a good proxy for σ_Q at that time. Results with constant, model-independent σ_Q for r(Ψ,ECS) are shown in Figs. 8 (H76) and A3 (two-box) for both piControl and historical experiments. $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ results are not shown as they are identical to r(Ψ,ECS). This is because the predictors, the set of {Ψ_i}, are all divided by the same constant.

Figure 8 Probability of obtaining r(Ψ,ECS) in the H76 CMIP emulator ensembles performing a piControl (a) and historical (b) simulation if each two-box emulator is given the same value of σ_Q=0.25 W  m⁻². The CMIP5 emulator ensemble histogram is given in blue, and the equivalent data for CMIP6 are shown in red. The full-complexity CMIP ensemble results performing the same experiment are shown as vertical dotted lines.

The constant σ_Q assumption can be seen to be good for the CMIP5 ensemble (blue) with full-complexity models (dotted blue line) agreeing well with likely values of the H76 and two-box CMIP5 emulator ensembles (blue histogram). However, the full-complexity CMIP6 ensemble (dotted red line) correlations are generally much lower than the CMIP6 emulators (red histogram). This is again supporting evidence that the underlying theory in CHW18 is sound FAPP. The similarity in the CMIP5 and CMIP6 histograms also suggests there is no real difference in the parameters of the emulator ensembles. The difference can be attributed to the amount of correlation between σ_Q and ECS in the CMIP5 and CMIP6 ensembles.

7 Discussion and conclusion

The aim of this paper was to understand why the strong emergent relationship from CHW18 found in the CMIP5 model ensemble weakened in the newer CMIP6 ensemble. This emergent relationship was based on reasonable, albeit simple, physical principles, and thus it is interesting (and important) to understand the differences between the theory and full-complexity models. A number of assumptions (Sect. 3) were made in deriving the theoretical emergent relationship between the predictor Ψ, a metric based on annual global mean temperature variability and ECS, the predictand in CHW18. We have shown the “no correlation between forcing and ECS” assumption no longer holds for the CMIP6 ensemble. In particular, the parameter σ_Q describing random forcing from internally generated variability, is correlated to ECS in CMIP6, and when this parameter is incorporated into the predictand, a good emergent relationship is recovered for both CMIP ensembles.

Assumption A3 stated that the forcing parameters, $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ and σ_Q, could be treated as constants across a model ensemble. While this is a fair assumption for $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ for both CMIP ensembles and σ_Q in the CMIP5 ensemble, we have shown that σ_Q is correlated to ECS in the CMIP6 ensemble. We have also shown that when the predictor of ECS is changed to $\mathrm{\Psi }/{\mathit{\sigma }}_Q$, good emergent relationships are recovered in both CMIP ensembles for both piControl and historical experiments. We also showed that pure theory could reproduce the full-complexity CMIP model results using H76 and two-box CMIP emulator ensembles. Although the proportionality between ECS and the predictor $\frac{Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}}{{\mathit{\sigma }}_Q}\mathrm{\Psi }$ has a high correlation and significance, simple pure theory underestimates the constant of proportionality. Aside from this, these results give us confidence the theoretical basis of CHW18 still applies to CMIP6 models as it did for CMIP5 FAPP. Testing the theoretical basis was the underlying aim of our study.

However, several questions remain. Can we estimate σ_Q from observations and therefore get an emergent constraint on ECS from the CMIP6 ensemble? Why is σ_Q correlated to ECS in CMIP6 and not CMIP5? σ_Q is a parameter designed to reproduce the observed global annual mean temperature variability, σ_T, in the non-chaotic H76 and two-box models. In the full-complexity models and the real world, this parameter attempts to capture chaotic internal variability and sub-annual (fast) feedbacks. It is fitted in this study using σ_T (an observable) as well as the unobservable two-box parameters. The reliance on these unobservable two-box parameters makes it appear that getting an estimate of σ_Q in the real world and thus an emergent constraint may be tough. However, there may be observable proxies for σ_Q that we have not yet found.

An obvious place to start looking for a proxy for σ_Q is in basic theory. The simplest example that one can imagine is

$$\begin{array}{}\text{(8)}& Q\left(t\right)=N\left(t\right)+\mathit{\lambda }T\left(t\right),\end{array}$$

where N(t) is the net top-of-the-atmosphere radiative flux. However this still requires knowledge of λ. Even given knowledge of λ it is well known (Forster, 2016) that N is poorly correlated to T where most of the change in N and T is driven by internal variability (although this relation works very well for large forced trends; e.g. the Gregory method works well applied to large stepped increases in CO₂). There are several models (Winton et al., 2010; Geoffroy et al., 2013a) and methods (Dessler et al., 2018; Bloch-Johnson et al., 2020) that get much better correlations between N and T when most of the changes are driven internally by taking into account the spatial distributions (the so-called pattern effect; see Armour et al., 2012). We leave this to a future study.

The question of why σ_Q is correlated to ECS in CMIP6 and not CMIP5 is also left unanswered. However, one can speculate why this may be the case. As previously mentioned σ_Q is a fitting parameter that is designed to capture the effect of chaotic internal variability and sub-annual (fast) feedbacks on global mean temperature variability. Zelinka et al. (2020) showed that the increased range of ECS in the CMIP6 models could be explained by the increased range in cloud feedbacks (see also Bock and Lauer, 2024). As σ_Q is fitted to annual temperature time series, some of this fast (sub-annual) cloud feedback effect could be included in σ_Q, correlating it to ECS. We leave concrete answers to a future study.

We have understood what assumption in the theoretical emergent relationship for CHW18 was responsible for the weakened correlation in CMIP6, namely that σ_Q is correlated with ECS. When accounted for, good emergent relationships are recovered. Although the information that the simple theory holds FAPP is useful and that σ_Q is correlated in CMIP6 with ECS is interesting, it is disappointingly not useful in constraining ECS due to the unobservable nature (we think) of σ_Q. In this sense, the method in CHW18 does not produce a useful emergent constraint on CMIP6 because the extra degree of freedom in σ_Q needs to be incorporated.

Schlund et al. (2020) tested 11 emergent constraints found in CMIP5, and nearly all of these got weaker in CMIP6. We do not know whether they failed for similar reasons. Indeed, many of them do not have a simple theoretical model as a basis for their emergent relationship, which means that assumption testing, the approach we follow in this paper, would be difficult to do. This is why we argue that emergent constraints should be based on a testable, falsifiable theoretical model. This aids understanding and lifts emergent constraint research from looking for strong correlations between variables to a more scientific approach of testing hypotheses of how the Earth system works. However, looking at all these other emergent constraints and identifying why they got weaker in CMIP6 would be very beneficial to understanding and useful to the community.

Emergent constraints based on theory with minimal degrees of freedom are most likely to be the most robust and useful. Constraints such as Hall and Qu (2006) on snow albedo feedback where the predictor (seasonal cycle snow albedo feedback) and predictand (climate change snow albedo feedback) are the essentially the same variable have been shown to be robust through three CMIP generations (Thackeray et al., 2021). Other constraints of this type that are likely to be more robust are the transient climate response constraints of Nijsse et al. (2020) and Tokarska et al. (2020) where near-term historical warming is the predictor of future longer-term warming.

Even if emergent relationships based on sound theoretical principles do fail there is still information to be gleaned on understanding why. Today, there is even more of an opportunity for the top-down insights of specific conceptual models to meet and complement the comprehensive, bottom-up approach from state-of-the-art climate models; there are many more high-quality observations; the global warming signal has also become clearer over time; and there is also a large archive of past and present climate model simulations.

Appendix A: Two-box CMIP emulators

The two-box model is H76's low thermal inertia atmosphere and well-mixed ocean surface layer with heat capacity C extended with a large-heat-capacity C₀ deep-ocean box coupled to the surface box by flux γ. This gives the model two timescales of adjustment, a fast (τ_f) and a slow e-folding time (τ_s). When fitted to CMIP models, typical values for the timescales are τ_f∼4 years and τ_s∼200 years (see Tables B5 and B6).

Each CMIP model labelled with i is “mimicked” by the two-box equations

$$\begin{array}{}\text{(A1)}& \begin{array}{rl}{C}_i{\displaystyle \frac{\mathrm{d}{T}_i}{\mathrm{d}t}}& =Q_i\left(t\right)-{\mathit{\lambda }}_i{T}_i\left(t\right)-{\mathit{\gamma }}_i\left(T_i\left(t\right)-T_{\mathrm{0}i}\left(t\right)\right),\\ C_{\mathrm{0}i}{\displaystyle \frac{\mathrm{d}{T}_{\mathrm{0}i}}{\mathrm{d}t}}& ={\mathit{\gamma }}_i\left(T_i\left(t\right)-T_{\mathrm{0}i}\left(t\right)\right).\end{array}\end{array}$$

T_0i is the annual global-mean deep-ocean temperature anomaly of model i. Parameters are fitted from the full-complexity abrupt-4xCO2 CMIP model experiments. The parameters λ and $Q_{\mathrm{2}\times {\mathrm{CO}}_{\mathrm{2}}}$ are determined from Gregory plots, while C, C₀ and γ are determined using Geoffroy's methodology (Geoffroy et al., 2013b). We use the published values of Geoffroy et al. (2013b) for CMIP5 models. Values for CMIP6 models are given in Tables B2 and B6.

The standard deviation of white-noise forcing σ_Q is fitted for each model from the global annual mean temperature time series of either the piControl or historical experiment. This time series is linearly detrended with a rolling 55-year window. This is to isolate the T(t) response to internal variability, analogous to how Ψ is determined in the CHW18 methodology, to leave the noisy T(t) response to white noise with standard deviation σ_T. The theoretical formula is given by Williamson et al. (2018)

$$\begin{array}{}\text{(A2)}& {\mathit{\sigma }}_T^{\mathrm{2}}={\displaystyle \frac{{\mathit{\sigma }}_Q^{\mathrm{2}}}{\mathrm{2}{\mathit{\lambda }}^{\mathrm{2}}}}\left({\displaystyle \frac{a_{\mathrm{f}}^{\mathrm{2}}}{{\mathit{\tau }}_{\mathrm{f}}}}+{\displaystyle \frac{a_{\mathrm{s}}^{\mathrm{2}}}{{\mathit{\tau }}_{\mathrm{s}}}}+{\displaystyle \frac{\mathrm{4}{a}_{\mathrm{f}}{a}_{\mathrm{s}}}{{\mathit{\tau }}_{\mathrm{f}}+{\mathit{\tau }}_{\mathrm{s}}}}\right).\end{array}$$

We rearrange this relation to get σ_Q. Values of λ, a_f, a_s, τ_f and τ_s are taken from Tables B1, B2, B5 and B6. Values of σ_Q in both historical and piControl runs are also reported in Tables B5 and B6. The parameters a_f, a_s, τ_f and τ_s are complicated functions of the parameters C, C₀, λ and γ. Their exact full functional forms can be found in Geoffroy et al. (2013b) and are not given here.

Figure A1 Probability of obtaining r(Ψ,ECS) (a) and $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ (b) in the two-box CMIP emulator ensembles performing a piControl simulation of 200 years. The CMIP5 emulator ensemble histogram is given in blue, and the equivalent data for CMIP6 are shown in red. The full-complexity CMIP ensemble results performing the same experiment are shown as vertical dotted lines.

Figure A2 Probability of obtaining r(Ψ,ECS) (a) and $r(\mathrm{\Psi }/{\mathit{\sigma }}_Q,\mathrm{ECS})$ (b) in the two-box CMIP emulator ensembles performing a historical simulation of the period 1880–2005. The CMIP5 emulator ensemble histogram is given in blue, and the equivalent data for CMIP6 are shown in red. The full-complexity CMIP ensemble results performing the same experiment are shown as vertical dotted lines.

Figure A3 Probability of obtaining r(Ψ,ECS) in the two-box CMIP emulator ensembles performing a piControl (a) and historical (b) simulation if each two-box emulator is given the same value of σ_Q=0.25 W m⁻². The CMIP5 emulator ensemble histogram is given in blue, and the equivalent data for CMIP6 are shown in red. The full-complexity CMIP ensemble results performing the same experiment are shown as vertical dotted lines.

Appendix B: Parameter values

Table B1 Gregory plot determined parameters for CMIP5 models from Geoffroy et al. (2013b).

Table B2 Gregory plot determined parameters for CMIP6 models.

Table B3 H76 model parameters fitted from abrupt4xCO2 runs for CMIP5 models. The σ_Q values are calculated from detrended T(t) for either historical or piControl runs.

Table B4 H76 model parameters fitted from abrupt4xCO2 runs for CMIP6 models. The σ_Q values are calculated from detrended T(t) for either historical or piControl runs.

Table B5 Two-box parameters fitted from abrupt4xCO2 runs for CMIP5 models taken from Geoffroy et al. (2013b). The σ_Q values are calculated from detrended T(t) for either historical or piControl runs.

Table B6 Two-box parameters fitted from abrupt4xCO2 runs for CMIP6 models. The σ_Q values are calculated from detrended T(t) for either historical or piControl runs.