2. Data and Methods
2.1. Last Millennium Global Mean Surface Temperature Data
We derived GMST records for the last millennium from a collection of 18 GCMs (Supplemental Table 2). Our ensemble spans multiple model generations because of the relatively low number of modeling centers that have published past1000/past2k simulations as a part of Paleoclimate Model Intercomparison Projects 3 and 4 (PMIP3, Braconnot et al., 2011; Braconnot et al., 2012 and PMIP4, Jungclaus et al., 2017; Kageyama et al., 2018). Ten of the GCMs assessed here come from CMIP5 (Taylor et al., 2012), seven come from CMIP6 (Eyring et al., 2016), and another stems from non-CMIP activities (CESM1-LME, Otto-Bliesner et al., 2016). For each GCM, GMST from the past1000 scenario (850-1849) was stitched with that of the historical scenario (1850-1999) to form a 1150-year timeseries. Model-based GMST anomalies were computed separately relative to either the 1800-1849 baseline (for simulation years before 1850) or the 1850-1899 baseline (for simulation years after 1850) to reduce any discontinuities between past1000 and historical experiments. For one model (MPI-ESM1-2-LR, Mauritsen et al., 2019), a past2k experiment was run instead of past1000. However, for consistency, we only consider the 850-1849 from the past2k experiment. For another model (FGOALS-gl; Zhou et al., 2018), we only consider the years that were simulated, 1000-1999.
To assess temperature variability in the GCMs over the past millennium, we leveraged an ensemble of Common Era GMST reconstructions published by the PAGES 2k consortium (PAGES 2k Consortium, 2019) and based on a global collection of temperature-sensitive proxies covering all or parts of the Common Era (PAGES 2k Consortium, 2017). Seven statistical methods were used to reconstruct GMST, such as the principal component regression (PCR) method shown in Figure 2, and 1000 realizations are available for each method to sample the uncertainty space. From these reconstructions, we sub-selected the years 850-1999 to match model output. We used 1800-1899 as a climatological baseline to compute the reconstruction based GMST anomaly so that the record would be easily compared to model output.
From 850-1849, the most dramatic external forcing stems from volcanic eruptions, which have a substantial cooling effect on GMST due to their abundant stratospheric aerosol emissions (PAGES 2k Consortium, 2019). These eruptions are also notoriously poorly represented in last millennium climate simulations, leading to potential systemic biases which can distort the emergent relationship (Eyring et al., 2019). To reduce the impact of major volcanic eruptions on GMST variability, we regressed the timeseries generated in Section 2.1 against the volcanic forcing profiles corresponding to each past1000 experiment (Crowley, 2000; Crowley et al., 2008; Gao et al., 2008; Toohey & Sigl, 2017). We then used the resulting residuals to compute temperature variability. This intervention substantially improved the strength of the emergent relationship with ECS, increasing the correlation coefficient from r=0.51 to r=0.62 for \(\psi\) and from r=0.42 to r=0.59 for\(\sigma_{b}\) (Supplemental Figure 1), so we used these residuals in our analysis (Figure 2b, 3, 4) instead of the original GMST timeseries (Figure 2a).
Historical external forcing is dominated by anthropogenic sources, which have a persistent, strong warming effect. Anthropogenic forcing has been shown to potentially impact the emergent relationship between interannual temperature variability (\(\psi\)) and ECS (Brown et al., 2018; Po-Chedley et al., 2018; Rypdal et al. 2018). To test the impact of anthropogenic forcing on our analysis, we modified records to exclude the historical (1850-1999) experiment and recomputed the emergent relationship. Doing so had minimal impact on the emergent relationship (Supplemental Figure 2). Therefore, we chose to include historical records in our overall analysis to reduce sampling uncertainty in the variability metric.
2.2. Calculating Temperature Variability Metrics and Equilibrium Climate Sensitivity
We computed interannual (\(\psi\), 2) and decadal trend (\(\sigma_{b}\), 3) temperature variability for each timeseries generated in Section 2.1. Observational uncertainty in each case was approximated as the standard deviation across the seven available reconstruction methods, while the observation mean was taken as the average of these methods (Supplemental Figure 3bd). Temperature variability metrics computed using the PAGES 2k reconstructions showed close agreement with those computed using the HadCRUT4 dataset (Morice et al., 2012) over the 1850-1999 period (Supplemental Figure 3ac). We removed volcanic forcing from reconstruction output via linear regression against the evolv2k version 3 volcanic forcing profile before computing both metrics (Toohey & Sigl, 2017).
The ECS values associated with each GCM were taken from Zelinka et al. (2020), which uses model output from the abrupt-4xCO2 experiment and the Gregory method (Gregory et al., 2004) to estimate ECS. In the case that Zelinka et al. (2020) did not calculate ECS for a GCM in our ensemble, we took ECS directly from publications pertaining to that model (Meehl et al., 2013; Phipps et al., 2012; Randall et al., 2007; Zhang et al., 2013). For each GCM, these values and their sources are listed in Supplemental Table 2.
2.3. Constraining Equilibrium Climate Sensitivity:
To infer possible values of ECS suggested by the emergent constraint, we employed the Hierarchical Emergent Constraint statistical framework (Bowman et al., 2018). This method offers a more comprehensive approach to constraining ECS compared to traditional techniques, as it accounts for the emergent constraint’s reliance on the signal-to-noise ratio and the correlation strength.