A- A+
Alt. Display

# Power Spectrum Sensitivity Analysis of the Global Mean Surface Temperature Fluctuations Simulated in a Two-Box Stochastic Energy Balance Model

## Abstract

Climate models used in theoretical studies and the long-term projections of climate change should be able to reproduce essential features of the Earth’s climate system including natural global scale variability on timescales from years to decades. It is notable then, that models simulate a very wide range in unforced “natural” variability, for example showing a 2 ½ fold range in standard deviation of decadal surface temperature. The global mean surface temperature (GMST) temporal fluctuations used as one of the main indicators of climate variability are usefully characterized by their power spectral density, which represents the distribution of temperature variance in the frequency domain. We applied a randomly-forced two-box energy balance model (EBM) with parameters that correspond to the Coupled Model Intercomparison Project Phase 5 (CMIP5) models to estimate the influence of such crucial aspects of the climate system as feedbacks, thermal inertia and deep ocean heat uptake on the power spectra of the GMST fluctuations (climate variability). These sensitivities can provide clues to allow us to better understand the reasons for the very wide range of climate variability derived from CMIP5 models. It is found that the influence of variations (uncertainties) in the EBM parameters on power spectra of the GMST fluctuations strongly depends on periods (frequencies) of these fluctuations. In particular, it was identified that the effect of variations in the feedback parameter significantly increases with increasing periods of GMST oscillations, while the influence of uncertainties in the climate thermal inertia parameter (effective heat capacity of the “atmosphere-mixed ocean layer” system) demonstrates the opposite behaviour. Variations in the deep-ocean heat uptake parameter tangibly affect GMST fluctuations on decadal and inter-decadal time scales. Meanwhile the uncertainty in the deep-ocean heat capacity parameter is minor for GMST fluctuations over all time scales.

Keywords:
How to Cite: Soldatenko, S.A. and Colman, R.A., 2022. Power Spectrum Sensitivity Analysis of the Global Mean Surface Temperature Fluctuations Simulated in a Two-Box Stochastic Energy Balance Model. Tellus A: Dynamic Meteorology and Oceanography, 74(1), pp.68–84. DOI: http://doi.org/10.16993/tellusa.40
Published on 22 Mar 2022
Accepted on 21 Feb 2022            Submitted on 21 Feb 2022

## 1. Introduction

Climate models of various degrees of complexity represent a comprehensive, powerful and essential tool for exploring the Earth’s climate system response to both human-caused and natural external forcing (Stocker, 2011; Collins et al., 2013; Dijkstra, 2013; Kaper and Engler, 2013; McGuffie and Henderson-Sellers, 2014; Gettelman and Rood, 2016; Soldatenko and Chichkine, 2016; Jeevanjee et al., 2017). Over the past few decades, these models have successfully predicted the subsequently observed upward trend in the temperature of the troposphere and the ocean, an increase in sea level, melting of Arctic sea ice, snow cover and other changes in the climate system (Collins et al., 2013; Hayhoe et al., 2017; Hausfather et al., 2020). Although the range of climate models from different research centres unanimously project future warming under the influence of human-induced radiative forcing, indices that characterise the climate system response to a given radiative forcing (e.g. equilibrium climate sensitivity and transient climate response) vary significantly (Rahmstorf, 2008; Bindoff et al., 2013; Kaya et al., 2016; Cox et al., 2018). Accordingly, numerical modelling shows that global and regional climate change projections possess large uncertainties for a given emission scenario (Collins et al., 2013; Meehl et al., 2013). For example, the projected range of global mean surface temperature (GMST) for the late 21st century relative to the reference period of 1986–2005 likely be 0.3°C to 1.7°C, 1.1°C to 2.6°C, 1.4°C to 3.1°C and 2.6°C to 4.8°C under the Representative Concentration Pathway scenarios RCP2.6, RCP4.5, RCP6.0 and RCP8.5, respectively (IPCC, 2013). Most of this uncertainty arises from inter-model differences in the strength of radiative feedbacks, particularly those of clouds (Flato et al., 2013). Moreover, to manage complex socio-economic systems, it is vital to quantify not just the spatial-temporal patterns of change in the average values of climate variables, but also their variability over a wide range of timescales, at least from years to decades (IPCC, 2014). Climate variability occurring on different timescales is the result of natural oscillations in the climate system (Monin, 1982; Kirtman et al., 2013). Arguably measures of unforced variability show as wide a statistical range across models is seen in model response under forced climate change. For example, the decadal variance of global and hemispheric temperatures in CMIP5 (Taylor et al., 2012) models deviate from each other by a factor of more than four (Colman and Power, 2018). However, the reasons for this large inter-model spread of surface temperature variability are not obvious. Furthermore, many studies have revealed that changes in the variance of GMST simulated by CMIP5 models differ from each other significantly on inter-annual, decadal and multi-decadal timescales (Power et al., 2017; Colman and Power, 2018).

In recently published papers (Soldatenko and Colman, 2019; Colman and Soldatenko, 2020), a simple two-box energy balance model (EBM) (Gregory, 2000; Held et al., 2010) proved useful as a method of understanding the essential effects of different physical mechanisms on annual, decadal and multi-decadal variability of the GMST. Although this is a highly simplified model, it was able to represent a significant part of the range of spread of unforced variability across GCMs. It does of course ignore dynamically coupled features such as ENSO, but was able to represent, for example, features of inter-model variance differences for both interannual and decadal timescales. One novel feature of the present approach is that parameters used in the EBM are derived directly from idealised climate change experiments from CMIP5 GCMs. This includes not just the parameters of the EBM itself, but the characteristics of the radiative stochastic forcing. This means that any skill found in the depiction of power spectra shape for GCMs using these parameters would hint at deep physical links between variability and response to external forcing such as increased greenhouse gases, as the same parameters describe just such a secular response. Soldatenko and Colman (2019) derived the Fokker-Plank equation for the two-box EBM and then found the analytical expression for the variance of GMST in terms of the model parameters. The obtained formula for the GMST variance enabled one to estimate the influence of parameters describing the radiative feedback mechanisms, the climate system thermal inertia characterized by the effective heat capacities of the upper and deep ocean layers, and the ocean heat uptake on the interannual, decadal and multi-decadal variability of the GMST. Nevertheless, the study was limited to specific timescale variations in the GMST, rather than using the full spectrum of timescale variability in the models, including the relationship between variability from short to very long periods. Yet characterising the entire range of variability is necessary since impacts are important across a full range of variability timescales, and phase differences in correlations between temperature changes and top of atmosphere radiation vary across different timescales due to oceanic heat uptake (Xie et al., 2016).

In this paper, we apply a sensitivity analysis to explore the influence of crucial physical mechanisms and aspects of the climate system, including radiative feedbacks, thermal inertia and deep ocean heat uptake on climate variability over a wide range of timescales. The essence of this approach is to use sensitivity functions to estimate how variations (uncertainties) in the parameters of stochastically forced two-box EBM affect the power spectral density (PSD) of the GMST fluctuations.

Section 2 will discuss the approach and methods used in this study and, in particular, some basic background on linear time-invariant systems with random inputs, and parametric uncertainty quantification in climate models. This section also describes the stochastic two-box EBM used in calculations, derivation of analytic solutions and sensitivity analysis methods. Sections 3 will present the sensitivity analysis of the power spectrum of the GMST fluctuations with respect to the EBM parameters, discussion and comparison with CMIP5 results and Section 4 will present concluding remarks.

## 2. Approaches and methods

### 2.1. Stochastically forced two-box energy balance model

The two-box EBM used in this study is one of the simplest and most widely used approximations of the global climate system (Gregory, 2000; Held et al., 2010), in which the atmosphere, the ocean mixed layer and the land surface are represented by the upper box, while the lower box describes the deep ocean. The state of the upper box is characterized by the GMST anomaly T, while the state of the lower box by the globally averaged temperature anomaly of the deep ocean TD. The equations of the stochastically forced two-box EBM can be written as (Colman and Soldatenko, 2020):

(1)
$C\frac{dT}{dt}=-\lambda T-\gamma \left(T-{T}_{D}\right)+{F}_{s},$
(2)
${C}_{D}\frac{d{T}_{D}}{dt}=\gamma \left(T-{T}_{D}\right),$

where C is the effective heat capacity of the upper box, λ is the proportionality coefficient between radiative response and global-mean warming that is usually called as the climate feedback parameter since the value of λ relies on feedback processes in the climate system, γ is the deep ocean heat uptake parameter describing the coupling mechanisms between the two boxes, CD is the effective heat capacity of the deep ocean, and Fs is the stochastic forcing approximated by Gaussian white noise (see Appendix A for detail).

The model (1)–(2) belongs to the class of simple climate models. Recall that EBMs were introduced almost simultaneously about 50 years ago by Budyko (1969) and Sellers (1969). The main purpose of EBMs is to obtain a better understanding of the Earth’s climate system, why the current state of the climate system is what it is, how sensitive is climate to external forcing, how feedbacks and climate system inertia affect the climate system, etc. (e.g., Ghil, 1984; Held and Suarez, 1974; North et al., 1981; Rogues et al., 2014). As shown in a number of studies (e.g. Hasselmann, 1976; Frankignoul and Hasselmann, 1977; Saravanan and McWilliams, 1998; Kleeman, 2002; Kleeman, 2011; Proistosescu et al., 2018; Cummins et al., 2020) EBMs with additive stochastic forcing parameterized by Gaussian white noise have proven to be useful for exploration of climate variability on time scales from years to decades. It is important that this approach is able to explain much of the observed climate variability on timescales varying from annual through to several decades. Despite the Gaussian white noise cannot describe the full spectrum of perturbations that affect the climate system, its use allows us partially explain natural variability in GMST (Cummins et al., 2020). In our study, the amplitude of white noise forcing was derived from simulation results using CMIP5 climate models. However, the actual picture of external forcing is very complex, and difficult, if not impossible, to be reproduced.

In addition to the feedback parameter λ, we shall consider the dimensionless feedback factor f defined by f = 1 – λ/λ0, where λ0 ≈ 3.3 W m–2 K–1 is the “Planck” feedback parameter, understood as the “black body” cooling response of the climate system with no feedbacks operating (Bony et al., 2006). In systems theory, the feedback factor represents the fraction of the system output signal fed back to the system input.

The two-box EBM has five control parameters, λ, C, CD, γ, and the magnitude of stochastic forcing qs, that have direct effects on the temporal evolution of EBM state variables T and TD. In our study, values of the parameters λ, C, CD and γ are set in accordance with the results obtained by Geoffroy et al. (2013) from the exploration of the 16 CMIP5 models under top-of-the-atmosphere radiative forcing resulting from an instantaneous quadrupling of atmospheric CO2 concentrations (Table 1). The values of the dimensionless feedback factor f listed in Table 1 are obtained using the relationship between λ and f.

Table 1

Values of EBM parameters and their multi-model mean and standard deviation given for the 16-model ensemble (all parameter values, with the exception of the climate feedback factor, are taken from Geoffroy et al., 2013).

MODEL PARAMETER

C(W yr m–2K–1) CD(W yr m–2K–1) γ(W m–2K–1) λ(W m–2K–1) f

1 BCC-CSM1-1 7.6 53 0.67 1.21 0.64

2 BNU-ESM 7.4 90 0.53 0.93 0.72

3 CanESM2 7.3 71 0.59 1.03 0.69

4 CCSM4 6.1 69 0.93 1.24 0.63

5 CNRM-CM5 8.4 99 0.50 1.11 0.67

6 CSIRO-Mk3.6.0 6.0 69 0.88 0.61 0.82

7 FGOALS-s2 7.0 127 0.76 0.88 0.74

8 GFDL-ESM2M 8.1 105 0.90 1.34 0.60

9 GISS-E2-R 4.7 126 1.16 1.70 0.49

10 HadGEM2-ES 6.5 82 0.55 0.65 0.81

11 INM-CM4 8.6 317 0.65 1.51 0.55

12 IPSL-CM5A-LR 7.7 95 0.59 0.79 0.76

13 MIROC5 8.3 145 0.76 1.58 0.53

14 MPI-ESM-LR 7.3 71 0.72 1.14 0.66

15 MRI-CGCM3 8.5 64 0.66 1.26 0.62

16 NorESM1-M 8.0 105 0.88 1.11 0.67

Mean 7.3 106 0.73 1.13 0.66

STD 1.1 62 0.18 0.31 0.09

The parameter qs that characterizes the magnitude of random radiative forcing, is estimated from the asymptotic form (Leith, 1973; Demchenko and Semenov, 2017): ${q}_{s}^{2}={\stackrel{˜}{\sigma }}_{s}^{2}{\tau }_{yr}$, where ${\stackrel{˜}{\sigma }}_{s}^{2}$ is the variance of radiative forcing averaged over the period τyr. If the duration of the averaging period is equal to one calendar year, τyr = 1 yr, and time is measured in years, as in this study, then simply ${q}_{s}^{2}={\stackrel{˜}{\sigma }}_{s}^{2}$. The estimates of the standard deviation ${\stackrel{˜}{\sigma }}_{s}$ of the global radiative forcing, calculated based on data of 16 models from the CMIP5 (Soldatenko and Colman, 2019), are in the range of 0.16–0.40 Wm–2 with a multi-model mean value of ~0.26 Wm–2. The base values of the parameters λ, C, CD and γ used in calculations represent the rounded values of the multi-model means of the CMIP5 fitted values (Table 2).

Table 2

Base values of model parameters and their ranges of change.

PARAMETER (p) MEAN VALUE (pavg) RANGE(pmin≤ p≤ pmax)

C, W yr m–2K–1 7.34 4.7 ≤ C ≤ 8.6

CD, W yr m–2K–1 105.5 53 ≤ CD ≤ 145

λ, W m–2K–1 1.13 0.61 ≤ λ ≤ 1.70

γ, W m–2K–1 0.73 0.50 ≤ γ ≤ 1.16

σs, W m–2 0.26 0.16 ≤ σs ≤ 0.40

f 0.66 0.49 ≤ f ≤ 0.82

### 2.2. Sensitivity of the GMST power spectrum to the EBM parameters

Equations (1) and (2) can be reduced to the second order differential equation that describes a randomly driven harmonic motion:

(3)
$\frac{d{T}_{D}^{2}}{d{t}^{2}}+2\beta \frac{d{T}_{D}}{dt}+{\omega }_{0}^{2}{T}_{D}=\frac{\gamma }{C{C}_{D}}{F}_{s},$

where β = [γ(C + CD) + λCD]/2CCD is the damping coefficient, and ${\omega }_{0}=\sqrt{\lambda \gamma /C{C}_{D}}$ is the natural frequency of the oscillator. The harmonic oscillator (3) is overdumped since ω0 < γ for all admissible values of its parameters listed in Table 2.

The frequency domain transfer function of the two-box EBM is derived by taking the Fourier transform of the differential equations (1) and (2) with zero initial conditions assuming that Fs(t) = δ(t):

(4)
$H\left(\omega \right)=\frac{\gamma +i\omega {C}_{D}}{C{C}_{D}\left[\left({\omega }_{0}^{2}-{\omega }^{2}\right)+i\omega 2\beta \right]}.$

Then, using (A4), we can find the PSD of GMST fluctuations:

(5)
${S}_{T}\left(\omega \right)=\frac{1}{\pi }\frac{{q}_{s}^{2}}{{C}^{2}{C}_{D}^{2}}\frac{{\gamma }^{2}+{\omega }^{2}{C}_{D}^{2}}{{\left({\omega }_{0}^{2}-{\omega }^{2}\right)}^{2}+4{\beta }^{2}{\omega }^{2}}.$

The integral of ST(ω) over positive frequencies is the variance of GMST anomaly ${\sigma }_{T}^{2}$ (Yaglom, 1987):

(6)
${\sigma }_{T}^{2}=\frac{{q}_{s}^{2}}{2\lambda C}\frac{\gamma C+\lambda {C}_{D}}{\gamma C+\left(\gamma +\lambda \right){C}_{D}}.$

By successively differentiating the equation (6) with respect to the parameters λ and C, we obtain the corresponding absolute sensitivity functions ψλ and ψc characterizing the effect of changes in parameters λ and C on the power spectrum of GMST fluctuations:

(7)
${\psi }_{\lambda }=-\frac{1}{\pi }\frac{2{q}_{s}^{2}}{{C}^{2}{C}_{D}^{2}}\frac{{\gamma }^{2}+{\omega }^{2}{C}_{D}^{2}}{{\left[{\left({\omega }_{0}^{2}-{\omega }^{2}\right)}^{2}+4{\beta }^{2}{\omega }^{2}\right]}^{2}}\frac{{\omega }_{0}^{2}\left({\omega }_{0}^{2}-{\omega }^{2}\right)C+2\beta \lambda {\omega }^{2}}{\lambda C}.$
(8)
${\psi }_{C}=-\frac{1}{\pi }\frac{2{q}_{s}^{2}}{{C}^{2}{C}_{D}^{2}}\frac{{\gamma }^{2}+{\omega }^{2}{C}_{D}^{2}}{{\left[{\left({\omega }_{0}^{2}-{\omega }^{2}\right)}^{2}+4{\beta }^{2}{\omega }^{2}\right]}^{2}}\frac{{\omega }^{2}}{{C}^{2}}\left[2\beta \left(2\beta C-\lambda -\gamma \right)-\left({\omega }_{0}^{2}-{\omega }^{2}\right)C\right].$

The corresponding relative sensitivity functions ${\psi }_{\lambda }^{R}$ and ${\psi }_{C}^{R}$ are as follows:

(9)
${\psi }_{\lambda }^{R}=-2\frac{{\omega }_{0}^{2}\left({\omega }_{0}^{2}-{\omega }^{2}\right)+\left(2\lambda \beta /C\right){\omega }^{2}}{{\left({\omega }_{0}^{2}-{\omega }^{2}\right)}^{2}+4{\beta }^{2}{\omega }^{2}},$
(10)
${\psi }_{C}^{R}=-2\frac{{\left({\omega }_{0}^{2}-{\omega }^{2}\right)}^{2}\left(1-{\omega }_{0}^{2}\right)+2\beta \left[2\beta -\left(\lambda +\gamma \right)/C\right]{\omega }^{2}}{{\left({\omega }_{0}^{2}-{\omega }^{2}\right)}^{2}+4{\beta }^{2}{\omega }^{2}}.$

Similarly, expressions for sensitivity functions ψγ, ${\psi }_{\gamma }^{R}$, ${\psi }_{{C}_{D}}$ and ${\psi }_{{C}_{D}}^{R}$ can be derived that characterize the influence of variations in the parameters γ and CD on the power spectrum ST(ω).

Under the asymptotic assumption γ → 0, the two-box model (1)–(2) degenerates into a “regular” single-box EBM (Kaper and Engler, 2013; McGuffie and Henderson-Sellers, 2014) describing the evolution of the GMST anomaly under the influence of radiative (deterministic and/or stochastic) forcing, defined at the top of the atmosphere:

(11)
$C\frac{dT}{dt}=-\lambda T+{F}_{s}.$

For this model, the PSD of GMST fluctuations is given by (Klyachkin, 2010)

(12)
${S}_{T}^{\left(1\right)}\left(\omega \right)=\frac{1}{\pi }\frac{{q}_{s}^{2}/{C}^{2}}{{\left(\lambda /C\right)}^{2}+{\omega }^{2}}.$

The absolute ${\psi }_{\lambda }^{\left(1\right)}$ and ${\psi }_{C}^{\left(1\right)}$ and the relative ${\psi }_{\lambda }^{R\left(1\right)}$ and ${\psi }_{C}^{R\left(1\right)}$ sensitivity functions characterizing the effect of changes in the parameters λ and C on the power spectrum of the GMST fluctuations in the one-box EBM (11) are of the form:

(13)
(14)

For the one-box EBM the variance of the GMST anomaly ${\sigma }_{1,T}^{2}$ is given by (Soldatenko and Colman, 2019)

(15)
${\sigma }_{1,T}^{2}=\frac{{q}_{s}^{2}}{2\lambda C}.$

Comparing the two expressions (6) and (15) shows that ${\sigma }_{T}^{2}<{\sigma }_{1,T}^{2}$ for typical values of all model parameters (Table 2), which means that the deep ocean dampens fluctuations in surface temperature, reducing their amplitude.

### 2.3. Data

The two-box EBM describes the real climate very simply, because it ignores the dynamics and other fundamental features of the climate system, focusing only on global thermodynamics and energetics. Nevertheless, the model has proven very convenient and useful in examining feedback mechanisms and their effects on the climate system response to external radiative forcing. Importantly, however, it is able to represent with skill the temporal behaviour of CMIP5 GCMs under different forcing, such as from 1% compounded increases in atmospheric CO2 concentrations (Geoffroy et al., 2013). Certainly, the results from the EBM should, to some extent, be consistent with those from observations and complex global climate models. Our confidence in the PSD estimates of the GMST fluctuations from the EBM with parameters listed in Table 1, and PSD sensitivities to EBM parameters is established by comparing the power spectrum derived from the EBM with the spectra derived from historical (200 years) runs of 16 CMIP5 models, for which the EBM parameters are available from Geoffroy et al. (2013).

## 3. Results and discussion

Before proceeding to the sensitivity analysis of power spectra of the GMST temporal variability, it is important to first consider the spectrum itself derived from both one- and two-box models. It has been found from both observations and models (Pelletier, 1998; North and Kim, 2017; Zhu et al., 2019) that temporal spectra of surface temperature fluctuations ST(ν) (here ν is the ordinary frequency) follow a power-law dependence on frequency ν: ST(ν) ∝ ν–β, where the scaling exponent β varies with frequency ranges. Comprehensive analysis of power spectra estimates for the GMST fluctuations derived from general circulation models and observations (Zhu, et al., 2019) shows that the scaling exponent is 0.80 ± 0.15 over the range of frequencies corresponding to time scales greater than a few months up to several decades. Similar estimates were obtained by Rypdal and Rypdal (2014) and Fredriksen and Rypdal (2016). Analysing PSDs derived from the historical (200 years) runs of 16 CMIP5 models listed in Table 1, we found that the ensemble average power spectrum, shown in Figure 1 (a), is approximated by power law with scaling exponent 0.82 ± 0.1 at the 95% confidence level for frequencies that correspond to time scales from years to decades. Meanwhile, dividing the ensemble average power spectrum into two segments, as illustrated in Figure 1 (b), we found that in the low-frequency part of the spectrum (decadal and inter-decadal variability), the ensemble average power spectrum is approximated by power law with scaling exponent β ≈ 0.40, while the power spectrum of interannual fluctuations is approximated by β ≈ 1.53.

Figure 1

Power spectra of the global mean surface temperature (GMST) fluctuations derived from historical (200 years) runs of 16 CMIP5 models: (a) the ensemble average power spectrum (red) and the characteristic 1/ν0.82 slope (blue); (b) the thin colored lines correspond to the individual ensemble members, while the thick black curve shows the ensemble average power spectrum, and the red lines show the characteristic 1/ν0.40 slope for frequencies less than 10–1 1/yr and the characteristic 1/ν1.53 for frequencies more than 10–1 1/yr.

Theoretical power spectra derived from one- and two-box EBMs in response to white noise forcing for the highest (fmax), lowest (fmin) and average (favd) values of feedback factor f found in the GCMs considered (see Table 2) are presented in Figure 2 on a log-log scale with base of 10. In these calculations, the remaining parameters, including stochastic forcing, hold their reference values. As previously established (e.g. Hall and Manabe, 1997; Fraedrich, 2001; North and Kim, 2017), the hallmark of the power spectrum derived from the one-box EBM is that at high frequencies (νλ/2πC) the power spectrum scales as ${S}_{T}^{\left(1\right)}\left(\nu \right)\propto {\nu }^{-2}$, which corresponds to so-called “red noise”. Meanwhile, at low frequencies (νλ/2πC), the power spectrum reaches a plateau region, where the PSD has a constant value that depends on the feedback strength: ${S}_{T}^{\left(1\right)}\left(\nu \right)\propto {\lambda }^{-2}$ (or ${S}_{T}^{\left(1\right)}\left(\nu \right)\propto {\left[{\lambda }_{0}\left(1-f\right)\right]}^{-2}$). The given expression shows that with the strengthening of feedbacks, the plateau constant value increases.

Figure 2

Power spectra of the global mean surface temperature fluctuations derived from the one- and two-box EBMs for different values of feedback factor f listed in Table 2. The orange dashed line shows the characteristic 1/ν2 slope.

From the characteristic frequency of the spectrum inflection, one can draw a conclusion on the feedback parameter: the lower the frequency at which the spectrum has an inflection, the stronger the feedback or, in other words, the larger the feedback factor f (see Figure 2). In a one-box EBM, the ratio λ/C defines the inflection point, the frequency ${\nu }_{c}^{\left(1\right)}$ at which the spectrum mode changes. For the reference values of model parameters we have .

A comparison between ${S}_{T}^{\left(1\right)}\left(\nu \right)$ and ST(ν) shows that for frequencies corresponding to multi-annual and longer time scales, the “two-box” PSD ST(ν), unlike the “one-box” PSD ${S}_{T}^{\left(1\right)}\left(\nu \right)$, can be described by a power law with a scaling exponent β ≈ 0.2–0.4 depending on the feedback factor: the stronger the feedback, the larger the scaling exponent. However, at high frequencies the spectrum of the two-box EBM is also approximated by ST(ν) ∝ ν–2. Since the bending point of the two-box EBM is determined by the ratio (λ + γ)/C, the frequency ${\nu }_{c}^{\left(2\right)}$, at which the spectrum regime changes occurs has the value of . Power spectra of the GMST fluctuations derived from the two-box EBM with 16 sets of parameters listed in Table 1 are shown in Figure 3. From this figure it is clearly seen that as the time scale of surface temperature oscillations increases (i.e. as the frequency decreases), the spread between the spectral density curves also increases. To characterize the ensemble spread we use the standard deviation with respect to the ensemble mean. For oscillations with periods 10, 30, 102 and 103 years, the standard deviations of ST(ν) are 2.19∙10–4, 1.14∙10–3, 2.74∙10–3 K2yr, respectively. Figure 4 shows the power spectrum calculated using ensemble averaging and the 95% confidence interval computed from model spread. At low frequencies, the “two-box” ensemble mean ST(ν) is approximated by a power law with a scaling exponent β ≈ 0.3.

Figure 3

Power spectra of the global mean surface temperature fluctuations derived from the two-box EBM with 16 parameter sets listed in Table 1. The thin colored lines correspond to the individual ensemble members.

Box-and-whisker plots shown in Figure 5 provides the descriptive statistics, such as the multi-model median, the interquartile model spread (the range between lower and upper quantiles), and the entire inter-model range, calculated from the ensemble simulations. These graphs drawn side-by-side, as displayed in Figure 5, show how the descriptive statistical measures change over frequency confirming the previously discussed fact that the spread of the multi-model spectral density ensemble increases with decreasing frequency of surface temperature fluctuations.

As one would expect (see for comparison Figures 1 and 4), the power spectra derived from the two-box EBM do not fully agree with power spectra estimated from complex climate models. The discrepancy between the spectra obtained using the EBM and CMIP5 models is mainly due to the fact that the EBM is an extremely simplified representation of the Earth’s climate system. In addition, it should be noted that the delta-correlated in time random process used in EBM is a very simple approximation of an external stochastic forcing, the patterns of which are in reality quite complex. Furthermore the parameters of the EBM are derived from climate change experiments, and not specifically selected for variability. Nevertheless, the EBM has some strength, in that it reproduces a gradient shift at around 10–1 per year, with an intermediary slope between 10–1 and 10–3 per year, suggesting its utility at longer timescales, and for exploring the ‘transition’ region between short and long timescales. This is consistent with findings that coupled dynamics are less important for global decadal variability than for interannual variability (Liu, 2012), where simple mixed layer oceans are found to have similar decadal variability to that of fully coupled ocean/atmosphere GCMs (Middlemas and Clement, 2016). Of cause, we should keep in mind that this is entirely correct for unforced GCMs, not of long, transient integrations. In addition, its simplicity permits explicit examination of the importance across the spectrum of feedbacks, as well as other factors such as stochastic forcing (see below). Therefore, the two-box EBM, despite its simplicity, remains a useful and convenient tool for exploring and understanding climate variability over a wide range of timescales (Soldatenko and Colman, 2019; Colman and Soldatenko, 2020). However, we should realize that the two-box EBM can only reproduce GCM results in an approximate way.

Figure 4

The ensemble average power spectrum of the global mean surface temperature fluctuations (red curve) derived from the two-box EBM with 16 sets of the parameters listed in Table 1. Grey shading shows the 95%confidence interval calculated from model spread. The blue line shows the characteristic 1/ν0.30 slope.

Figure 5

The box-and-whisker plots showing the temporal changes in power spectra across an ensemble of 16 EBMs.

Since the essential physical mechanisms such as radiative feedbacks, deep ocean heat uptake and climate system inertia, that control climate formation and variability, are parametrically described by the EBM, we can explore the influence on the power spectrum of GMST oscillations from the parametric ranges considered as sources of “uncertainty”. The vector of independent parameters for the two-box EBM has four components, viz λ, C, γ and CD, while the parameter vector for the one-box EBM has only two elements, namely λ and C. We deliberately excluded from consideration the stochastic forcing, since the PSDs of the one-box and two-box EBMs linearly depend on the variance of radiative forcing. In other words, the spectrum’s shape is invariant under rescaling.

To examine the influence of the range in these parameters suggested by GCMs on PSD of GMST variability we apply a sensitivity analysis using the OFAT approach described in Appendix 2. The main weakness of this approach is its inability to investigate interactions between parameters. However, since in the EBM interactions among the parameters are not considered and variations in the model parameters are sufficiently small, the use of OFAT approach is justified (Soldatenko and Colman, 2019).

Figures 6 and 7 show ASFs for one- and two-box EBMs, respectively, calculated for the highest (fmax), lowest (fmin) and average (favd) values of the feedback factor f (Table 2). We deliberately vary the feedback factor over such a wide range, since the uncertainty in feedbacks is one of the key sources of uncertainty in climate system response to external radiative forcing (Bony et al., 2006). As can be seen from these figures, for a given feedback factor, the absolute value of sensitivity functions with respect to the parameter λ for both one- and two-box EBMs increase as the frequency ν decreases. The feedback parameter determines how quickly the modulus of ψλ increases: the larger the feedback factor, the greater (in modulus) the sensitivity functions.

Figure 6

Absolute sensitivity functions ${\psi }_{\lambda }^{\left(1\right)}$ and ${\psi }_{C}^{\left(1\right)}$ for the one-box EBM power spectral density with respect to parameters λ and C calculated for the highest (fmax), lowest (fmin) and average (favd) values of feedback factor f listed in Table 2.

Figure 7

Absolute sensitivity functions ψλ, ψC, ψγ and ψCD for the two-box EBM power spectral density with respect to parameters λ, C, γ and CD, respectively, calculated for the highest (fmax), lowest (fmin) and average (favd) values of feedback factor f listed in Table 2.

For any particular frequency ν, the impact of the “uncertainty” (here interpreted as being a measure of the GCM derived range) in some parameter α on the uncertainty in power spectrum ST can be quantified as follows: δ(ST) ≈ ±ST |δα|, where the parameter variation δα characterises the absolute uncertainty in the parameter α, while δ(ST) is the variation (uncertainty) in the PSD caused by δα. In addition to ASFs, we also computed RSFs (see Figures 8 and 9), which allow us, first, to quantitatively estimate the fractional (percentage) uncertainty in the power spectrum due to the uncertainty in the corresponding parameter and, second, to rank the relative importance of the parameters according to their degree of influence on power spectrum. Recall that the fractional uncertainty in the power spectrum due to the uncertainty in the parameter α is given by $\delta {S}_{T}/\left|{S}_{T}\right|\approx ±{\psi }_{\alpha }^{R}\left|\delta \alpha /{\alpha }_{0}\right|×100%$.

Figure 8

Relative sensitivity functions ${\psi }_{\lambda }^{R\left(1\right)}$ and ${\psi }_{C}^{R\left(1\right)}$ for the one-box EBM power spectral density with respect to parameters λ and C calculated for the highest (fmax), lowest (fmin) and average (favd) values of feedback factor f listed in Table 2.

Figure 9

Relative sensitivity functions ${\psi }_{\lambda }^{R}$, ${\psi }_{C}^{R}$, ${\psi }_{\gamma }^{R}$ and ${\psi }_{{C}_{D}}^{R}$ for the two-box EBM power spectral density with respect to parameters λ, C, γ and CD calculated for the highest (fmax), lowest (fmin) and average (favd) values of feedback factor f listed in see Table 2.

As an example, assume the one-sigma absolute uncertainty in all parameters (see Table 1), i.e. δλ = 0.31 W m–2 K–1, δC = 1.1 W yr m–2 K–1, δγ = 0.18 W m–2 K–1 and δCD = 62 W yr m–2 K–1. In other words, we suppose that the parameter variations are one-sigma deviations from the GCM means. Table 3 shows how the uncertainties in model parameters affect both the absolute and relative uncertainty in the power spectrum of the GMST fluctuations with periods 2, 10, 30, 102 and 103 years. The corresponding ASFs and RSFs are also presented in Table 3.

Table 3

The modulus of absolute and relative sensitivity functions with respect to the two-box EMB parameters, and the corresponding absolute δ(ST) (K2yr) and relative [δ(ST)/ST] (%) uncertainties in power spectrum caused by one-sigma uncertainty in model parameters.

PARAMETER λ (W m–2K–1) C(W yr m–2K–1) CD(W yr m–2K–1) γ (W m–2K–1)

One-sigma parameter uncertainties from GCMs ±0.31 ±1.10 ±62.60 ±0.18

Period of GMST fluctuations T = 2 yr

|ψα| 2.79 × 10–5 1.09 × 10–5 5.28 × 10–11 2.95 × 10–7

$|{\psi }_{\alpha }^{R}|$ 0.008 1.99 1.38 × 10–4 0.005

δ(ST) (K2yr) ±8.66 × 10–8 ±1.20 × 10–5 ±3.30 × 10–9 ±5.30 × 10–8

[δ(ST)/ST] × 100% ±0.2 ±29.8 ±0.008 ±0.1

Period of GMST fluctuations T = 10 yr

|ψα| 1.13 × 10–4 2.03 × 10–4 2.44 × 10–8 1.37 × 10–4

$|{\psi }_{\alpha }^{R}|$ 0.17 1.72 0.003 0.12

δ(ST) (K2yr) ±4.03 × 10–5 ±2.23 × 10–4 ±1.53 × 10–6 ±2.47 × 10–5

[δ(ST)/ST] × 100% ±4.7 ±25.7 ±0.2 ±2.8

Period of GMST fluctuations T = 30 yr

|ψα| 2.30 × 10–3 4.05 × 10–4 4.06 × 10–7 2.42 × 10–3

$|{\psi }_{\alpha }^{R}|$ 0.71 0.81 0.01 0.48

δ(ST) (K2yr) ±7.14 × 10–4 ±4.46 × 10–4 ±2.54 × 10–5 ±4.36 × 10–4

[δ(ST)/ST] × 100% ±19.6 ±12.2 ±0.7 ±11.9

Period of GMST fluctuations T = 102 yr

|ψα| 5.76 × 10–3 1.06 × 10–4 2.92 × 10–7 5.92 × 10–3

$|{\psi }_{\alpha }^{R}|$ 1.12 0.13 0.005 0.75

δ(ST) (K2yr) ±1.793 × 10–3 ±1.17 × 10–4 ±1.83 × 10–5 ±1.06 × 10–3

[δ(ST)/ST] × 100% ±30.9 ±2.0 ±0.3 ±18.4

Period of GMST fluctuations T = 103 yr

|ψα| 1.19 × 10–2 2.09 × 10–5 4.19 × 10–5 9.83 × 10–4

$|{\psi }_{\alpha }^{R}|$ 1.44 0.02 0.47 0.08

δ(ST) (K2yr) ±3.69 × 10–3 2.30 × 10–5 ±2.59 × 10–3 1.77 × 10–4

[δ(ST)/ST] × 100% ±39.4 ±0.2 ±28.0 ±1.9

Analysis of the results in Table 3 and Figures 69 shows the following:

1. The ASFs change significantly with frequency and all ASFs are nonlinear.
2. ASFs with respect to parameters λ, γ and C are quantities of the same order of magnitude for fluctuations with period of 10 yr. While the oscillation period increases, only the ASF with respect to the feedback parameter λ increases in magnitude, whilst the remaining ASFs demonstrate more complex behaviour: starting out with a high-frequency domain of the spectrum, the magnitudes of ASFs with respect to parameters C, CD and γ gradually increase in modules with an increase in the period of fluctuations; once the oscillation frequency exceeds a certain (critical) value (e.g. for the parameter C the critical value is ν* = (λ + γ)/2πC), the magnitudes of ASFs begin to reduce in absolute values. Physically this is consistent with the mixed layer ocean only influencing variability over ‘intermediate’ frequencies, at long periods the mixed layer is not important. An important finding of the current research is that radiative feedbacks operate strongly across the full range of frequencies, and continue to strengthen at very long timescales.
3. Importantly too the effect of all parameter ranges is greatly strengthened for strong feedback compared with weak. Again, this is physically understandable as weak feedback results in a strongly radiatively damped response, whereas stronger (positive) feedback produces stronger oscillations at all frequencies. Moreover, these differences are more marked at long timescales when the ‘full effect’ of radiative amplification is seen dynamically.
4. The effect of uncertainties in the EBM parameters δα, where α = (λ, C, CD, γ), on the uncertainty in power spectrum δST caused by δα is linearly dependent on the ASFs (recall that δ(ST) ≈ ±ST |δα|). Therefore the behaviour of δST is similar to the behaviour of the ASFs.
5. The percentage uncertainty in the power spectrum (δST/ST) × 100% due to the uncertainty in the feedback parameter λ increases with the decrease in fluctuation frequency: from ±5% at ν =1/(10 yr) to ±40% at ν = 1/103 yr. Meanwhile, the effect of the uncertainty in the parameter C is diametrically opposed: the relative uncertainty in the PSD decreases from ±26% at ν =1/10 yr to negligible value at ν = 1/103 yr. Again this is consistent with the negligible role of mixed layer variations at very low frequencies, but the stronger role of feedbacks at longer timescales.
6. The uncertainty in the parameter γ that characterises the deep ocean heat uptake, affects the relative uncertainty in the PSD for decadal (Θ = 10 yr) and inter-decadal (Θ = 30 yr) oscillations. Obviously, on short time scales, the role of the deep ocean in the formation of surface temperature variability is insignificant. On longer time scales (centennial time scales and more), the effect of the deep ocean on surface temperature fluctuations is also insignificant, since the climate system approaches a state of equilibrium, at which there are no surface temperature fluctuations at all.
7. The uncertainty in the parameter CD that characterises the deep ocean heat capacity is tangible only for very low frequency oscillations, which is consistent with its role mainly at very long timescales. Here, we need to clarify the difference between the effect of γ and CD at very long timescales. At extremely long timescales the two-box EBM approaches the approximation of a single, very deep ocean layer. Therefore, the parameter γ becomes unimportant, in contrast to the parameter CD, which keeps its significant role in modelling very low frequency oscillations.
8. The RSFs, which allows us to rank the impact of EBM parameters on the power spectrum of the GMST fluctuations, show that there is intrinsic nonlinearity in all RSFs and these functions themselves vary significantly over the frequency domain. Whereas RSFs with respect parameters λ and C exhibit monotonic behavior, the behavior of RSFs with respect to parameters γ and CD is clearly non-monotonic and resembles an inverted bell. The extrema of RSFs with respect to parameters γ and CD are mainly determined by the feedback parameter and the heat capacity of the mixed ocean layer, which, in turn, determine the time scale of surface temperature fluctuations.

Figure 10 gives a visual representation on the rank of the EBM parameters in accordance with their degree of influence on PSD at different frequencies. Figure 9 (a) shows that there is a certain critical frequency ν* at which the relative roles of the parameters λ and C are reversed. For the one-box EBM we can find that ν* = λ/2πC, while for the two-box EMB ν* = (λ + γ)/2πC. We already saw that at these frequency points the one- and two-box ASFs with respect to feedback parameter have extreme values. As discussed above, the frequency ν* defines the bending point at which the spectrum regime changes. Figure 10 (b) allows us to evaluate the relative importance of the two-box EBM parameters at different frequencies. For example, for oscillations with period of 10 yr the parameter C has the largest rank. Feedback parameter λ and heat exchange parameter γ rank second and third, respectively, followed by deep ocean effective heat capacity CD. For fluctuations with period of 100 yr the feedback parameter ranks first, followed by heat exchange parameter and then upper-box heat capacity. The relative role of the parameter CD remains negligible except at extremely long timescales. It is important to emphasize that varying the parameters within their ranges (see Table 2) does not qualitatively change the character of RSFs behavior and, consequently, the rank of the model parameters.

Figure 10

Relative sensitivity functions of the power spectral density of the global mean surface temperature fluctuations with respect to (a) one- and (b) two-box EBMs parameters calculated for the average value of feedback factor favd.

## 4. Concluding remarks

Climate models, ranging from simple conceptual models to highly complex coupled models of the Earth’s climate system, serve as one of the main tools for the theoretical study and projection of climate dynamics, climate change and variability. Simple climate models, such as the stochastically-forced EBM used in this paper are useful for studying new ideas and various problems arising in climate research when coupled climate models, due to their complexity and computational expense, are impractical, inefficient or difficult to understand. However, the results of simulations obtained using simple models should be at least qualitatively consistent with observations and results from complex coupled climate models. Recent studies (Soldatenko and Colman, 2019; Colman and Soldatenko, 2020) have analysed the ability of a two-box stochastic EBM to simulate climate variability, and also estimated the influence of EBM parameters on annual, decadal, and multi-decadal variability of GMST fluctuations. Based on this analysis, it was concluded that the adequacy of the two-box stochastic EBM is acceptable enough for obtaining “first-order” estimates results on variability and to study the sensitivity of the model to parameter changes corresponding to large scale features in climate models, such as feedbacks, mixed layer depths, stochastic forcing and deep ocean mixing. In this article, which is a logical continuation of our previous research, we, firstly, investigate the power spectrum of GMST fluctuations produced by both two-box and one-box stochastic EBM and, secondly, analyze the sensitivity of power spectra with respect to parameters derived from GCMs. This is important for two main reasons: (1) GMST and its fluctuations are used as one of the key indicators of climate change and variability, and (2) the power spectrum reflects the amplitude of the GMST fluctuations at different frequencies.

We found that the PSD of GMST fluctuations derived from the two-box EBM can be approximated by a power law with a scaling exponent of 0.82 ± 0.1 at the 95% confidence level for frequencies that correspond to time scales from years to decades, while the scaling exponent for power spectrum estimates of GMST fluctuations derived from complex climate models and observations (Zhu et al., 2019) is 0.80 ± 0.15 over the range of frequencies corresponding to time scales greater than a few months up to several decades. This suggests that the EBM fairly realistically reproduces the spectrum of surface temperature fluctuations on inter-annual and inter-decadal time scales. The computational results also show the acceptable agreement between the EBM’s power spectrum and PSDs derived from historical (200 years) runs of 16 CMIP5 models.

Since the parameters of the EBM are uncertain, we used the absolute and relative sensitivity functions to explore how the spread in the EMB parameters derived from GCMs affect the power spectrum of the GMST fluctuations. We have obtained analytic expressions for the absolute and relative sensitivity functions via the transfer function of the EBM assuming the external random forcing is represented by Gaussian white noise.

Using sensitivity functions, we found that the influence of variations (uncertainties) in the EBM parameters on power spectrum of the GMST fluctuations strongly depends on periods (frequencies) of these fluctuations. In particular, it was identified that the effect of variations in the feedback parameter significantly increases with increasing periods of GMST oscillations, while the influence of uncertainties in the climate ‘mixed layer’ thermal inertia parameter (effective heat capacity of the “atmosphere-mixing ocean layer” system) exhibits the opposite behaviour. Variations in the deep-ocean heat uptake parameter tangibly affect GMST fluctuations on decadal and inter-decadal time scales. Meanwhile the uncertainty in the deep-ocean heat capacity parameter is minor for GMST fluctuations over all time scales. The advantage of this is it provides a clear, physically understandable, picture of the major factors expected to affect variability across the time spectrum, and of the key ‘parameters’ that determine them. This combined with the actual range in parameters derived from GCMs under simplified climate change scenarios highlight their expected importance at a given timeframe. For example, the role of radiative feedbacks is seen to become more important as timescales increase suggesting better ‘analogues’ between variability and climate change at longer timescales. Furthermore it can help suggest avenues for exploration of the source of ranges in GCM variability that are currently poorly understood, such as the large range in decadal variability found across GCMs (Colman and Power, 2018).

The results obtained should be interpreted with caution primarily due to the simplicity of the models used in this study. We associate our further research with the development of a more complex and realistic climate model, the analysis on its basis of parametric sensitivity and uncertainties in climate variability simulations, as well as the study of various stochastic processes used to describe the external random forcing.

Appendix A

A linear time-invariant system with random input. DOI: https://doi.org/10.16993/tellusa.40.s1

Appendix B

Quantification of parametric uncertainty. DOI: https://doi.org/10.16993/tellusa.40.s2

## Acknowledgements

We sincerely thank Scott Power, Meelis Zidikheri and two anonymous reviewers for critically reading the manuscript and suggesting substantial improvements. We acknowledge the World Climate Research Programme’s Working Group on Coupled Modelling, which is responsible for CMIP, and we thank the climate modelling groups of the Australian Bureau of Meteorology for producing and making available their model output.

## Competing Interests

The authors have no competing interests to declare.

## References

1. Bindoff, NL, Stott, PA, AchutaRao, KM, Allen, MR, Gillett, N, Gutzler, D, Hansingo, K, Hegerl, G, Hu, Y, Jain, S, Mokhov, II, Overland, J, Perlwitz, J, Sebbari, R and Zhang, X. 2013. Detection and attribution of climate change: From global to regional. In: Stocker, TF, Qin, D, Plattner, G-K, Tignor, M, Allen, SK, Doschung, J, Nauels, A, Xia, Y, Bex, V and Midgley, PM (eds.), Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press. pp. 867–952. DOI: https://doi.org/10.1017/CBO9781107415324.022

2. Bony, S, Colman, R, Kattsov, VM, Allan, RP, Bretherton, CS, Dufresne, JL, Hall, A, Hallegatte, S, Holland, MM, Ingram, W and Randall, DA. 2006. How well do we understand and evaluate climate change feedback processes? J. Clim., 19(15): 3445–3482. DOI: https://doi.org/10.1175/JCLI3819.1

3. Budyko, MI. 1969. The effect of solar radiation variations on the climate of the Earth. Tellus, 21: 611–619. DOI: https://doi.org/10.3402/tellusa.v21i5.10109

4. Collins, M, Knutti, R, Arblaster, J, Dufresne, J-L, Fichefet, T, Friedlingstein, P, Gao, X, Gutowski, WJ, Johns, T, Krinner, G, Shongwe, M, Tebaldi, C, Weaver, AJ and Wehner, M. 2013. Long-term Climate Change: Projections, Commitments and Irreversibility. In: Stocker, TF, Qin, D, Plattner, G-K, Tignor, M, Allen, SK, Doschung, J, Nauels, A, Xia, Y, Bex, V and Midgley, PM (eds.), Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press. pp. 1029–1136.

5. Colman, RA and Power, SB. 2018. What can decadal variability tell us about climate feedbacks and sensitivity? Clim. Dyn., 51: 3815–3828. DOI: https://doi.org/10.1007/s00382-018-4113-7

6. Colman, R and Soldatenko, S. 2020. Understanding the links between climate feedbacks, variability and change using a two-layer energy balance model. Clim. Dyn., 54: 3441–3459. DOI: https://doi.org/10.1007/s00382-020-05189-3

7. Cox, PM, Huntingford, C and Williamson, MS. 2018. Emergent constraint on equilibrium climate sensitivity from global temperature variability. Nature, 553: 319–322. DOI: https://doi.org/10.1038/nature25450

8. Cummins, DP, Stephenson, DB and Stott, PA. 2020. Optimal estimation of stochastic energy balance model parameters. J. Clim., 33: 7909–7926. DOI: https://doi.org/10.1175/JCLI-D-19-0589.1

9. Demchenko, PF and Semenov, VA. 2017. Estimation of uncertainty in surface air temperature climatic trends related to the internal dynamics of the atmosphere. Dokl. Earth Sc., 476: 1105–1108. DOI: https://doi.org/10.1134/S1028334X17090239

10. Dijkstra, HA. 2013. Nonlinear Climate Dynamics. Cambridge: Cambridge University Press, p. 367.

11. Flato, G, Marotzke, J, Abiodun, B, Braconnot, P, Chou, SC, Collins, W, Cox, P, Driouech, F, Emori, S, Eyring, V, Forest, C, Gleckler, P, Guilyardi, E, Jakob, C, Kattsov, V, Reason, C and Rummukainen, M. 2013. Evaluation of Climate Models. In: Stocker, TF, Qin, D, Plattner, G-K, Tignor, M, Allen, SK, Doschung, J, Nauels, A, Xia, Y, Bex, V and Midgley, PM (eds.), Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press. pp. 741–866.

12. Fraedrich, K. 2001. Simple climate models. In: Imkeller, P and von Storch, JS (eds.), Stochastic Climate Models. Progress in Probability, 49: 65–100. Basel: Birkhäuser. DOI: https://doi.org/10.1007/978-3-0348-8287-3_2

13. Frankignoul, C and Hasselmann, K. 1977. Stochastic climate models, part II. Application to sea surface temperature anomalies and thermocline variability/ Tellus, 29: 289–305. DOI: https://doi.org/10.3402/tellusa.v29i4.11362

14. Fredriksen, H and Rypdal, K. 2016. Spectral characteristics of instrumental and climate model surface temperatures. J. Clim., 29: 1253–1268. DOI: https://doi.org/10.1175/JCLI-D-15-0457.1

15. Geoffroy, O, Saint-Martin, D, Olivie, DJL, Voldoire, A, Bellon, G and Tytéca, S. 2013. Transient climate response in a two-layer energy-balance model. Part I: Analytical solution and parameter calibration using representation of the efficacy of deep-ocean heat uptake and validation for CMIP5 AOGCMs. J. Clim., 26: 1841–1857. DOI: https://doi.org/10.1175/JCLI-D-12-00195.1

16. Gettelman, A and Rood, RB. 2016. Demystifying climate models: A user’s guide to Earth System Models. Berlin-Heidelberg: Springer-Verlag, p. 247. DOI: https://doi.org/10.1007/978-3-662-48959-8

17. Ghil, M. 1984. Climate sensitivity, energy balance models, and oscillatory climate models. J. Geophys. Res., 89: 1280–1284. DOI: https://doi.org/10.1029/JD089iD01p01280

18. Gregory, JM. 2000. Vertical heat transport in the ocean and their effect on time-dependent climate change. Clim. Dyn., 16: 501–515. DOI: https://doi.org/10.1007/s003820000059

19. Hall, A and Manabe, S. 1997. Can local linear stochastic theory explain sea surface temperature and salinity variability? Clim. Dyn., 13(3): 167–180. DOI: https://doi.org/10.1007/s003820050158

20. Hasselmann, K. 1976. Stochastic climate models. Part I. Theory. Tellus, 28: 473–485. DOI: https://doi.org/10.3402/tellusa.v28i6.11316

21. Hayhoe, K, Edmonds, J, Kopp, RE, LeGrande, AN, Sanderson, BM, Wehner, MF and Wuebbles, DJ. 2017. Climate models, scenarios, and projections. In: Wuebbles, DJ, Fahey, DW, Hibbard, KA, Dokken, DJ, Stewart, BC and Maycock, TK (eds.), Climate Science Special Report: Fourth National Climate Assessment, I. Washington, DC, USA: U.S. Global Change Research Program. pp. 133–160. DOI: https://doi.org/10.1029/2019GL085378

22. Hausfather, Z, Drake, HF, Abbott, T and Schmidt, GA. 2020. Evaluating the performance of past climate model projections. Geophys. Res. Lett., 47: e2019GL085378. DOI: https://doi.org/10.3402/tellusa.v26i6.9870

23. Held, IM and Suarez, M. 1974. Simple albedo feedback models of the icecaps. Tellus, 36: 613–629. DOI: https://doi.org/10.1175/2009JCLI3466.1

24. Held, IM, Winton, M, Takahashi, K, Delworth, T, Zeng, F and Vallis, GK. 2010. Probing the fast and slow components of global warming by returning abruptly to preindustrial forcing. J. Clim., 23: 2418–2427. DOI: https://doi.org/10.1175/2009JCLI3466.1

25. IPCC. 2013. Summary for Policymakers. In: Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change [Stocker, TF, Qin, D, Plattner, G-K, Tignor, M, Allen, SK, Boschung, J, Nauels, A, Xia, Y, Bex, V and Midgley, PM (eds.)], Climate Change 2013: The Physical Science Basis. Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press. p. 1535.

26. IPCC. 2014. Climate Change 2014: Synthesis Report. In: Core Writing Team, Pachauri, RK and Meyer, LA (eds.), Contribution of Working Groups I, II and III to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. Geneva, Switzerland: IPCC. p. 222.

27. Jeevanjee, N, Hassanzadeh, P, Hill, S and Sheshardi, A. 2017. A perspective on climate model hierarchies. J. Adv. Modelling Earth systems, 9: 1760–1771. DOI: https://doi.org/10.1002/2017MS001038

28. Kaper, H and Engler, H. 2013. Mathematics and climate. SIAM, Philadelphia, USA, p. 295. DOI: https://doi.org/10.1137/1.9781611972610

29. Kaya, Y, Yamaguchi, M and Akimoto, K. 2016. The uncertainty of climate sensitivity and its implication for the Paris negotiation. Sustainability Science, 11: 515–518. DOI: https://doi.org/10.1007/s11625-015-0339-z

30. Kirtman, B, Power, SB, Adedoyin, JA, Boer, GJ, Bojariu, R, Camilloni, I, Doblas-Reyes, FJ, Fiore, AM, Kimoto, M, Meehl, GA, Prather, M, Sarr, A, Schär, C, Sutton, R, van Oldenborgh, GJ, Vecchi, G and Wang, HJ. 2013. Near-term Climate Change: Projections and Predictability. In: Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change [Stocker, TF, Qin, D, Plattner, G-K, Tignor, M, Allen, SK, Doschung, J, Nauels, A, Xia, Y, Bex, V and Midgley, PM (eds.)], Climate Change 2013: The Physical Science Basis. Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press. pp. 741–866.

31. Kleeman, R. 2002. Measuring dynamical prediction utility using relative entropy. J. Atmos. Sci., 59: 2057–2072. DOI: https://doi.org/10.1175/1520-0469(2002)059<2057:MDPUUR>2.0.CO;2

32. Kleeman, R. 2011. Spectral analysis of multi-dimensional stochastic geophysical models with an application to decadal ENSO variability. J. Atmos. Sci., 68: 13–25. DOI: https://doi.org/10.1175/2010JAS3546.1

33. Klyachkin, VI. 2010. Lectures on Dynamics of Stochastic Systems. Amsterdam: Elsevier, p. 410.

34. Leith, CE. 1973. The standard error of time-average estimates of climate means. J. Appl. Meteorol., 12: 1066–1069. DOI: https://doi.org/10.1175/1520-0450(1973)012<1066:TSEOTA>2.0.CO;2

35. Middlemas, E and Clement, A. 2016. Spatial patterns and frequency of unforced decadal scale changes in global mean surface temperature in climate models. Journal of Climate, 29: 6245–6257. DOI: https://doi.org/10.1175/JCLI-D-15-0609.1

36. Monin, AS. 1982. Introduction to the theory of climate. Leningrad, USSR: Hydrometeorological Publ., p. 246.

37. McGuffie, K and Henderson-Sellers, A. 2014. The climate modelling primer. 4th ed. Chichester, UK: Wiley-Blackwell. p. 456.

38. North, GR, Cahalan, RF and Coakley, JA. 1981. Energy balance climate models. Rev. Geophys. Space Phys., 19: 91–121. DOI: https://doi.org/10.1029/RG019i001p00091

39. North, GR and Kim, K-Y. 2017. Energy balance climate models. Weinheim, Germany: Wiley-VCH, p. 392. ISBN 9783527411320. DOI: https://doi.org/10.1002/9783527698844

40. Power, SB, Surral, R, Chung, C, Colman, RA, Kharin, V, Boer, G, Gergis, J, Henley, B, McGregor, S, Arblaster, JM, Holbrook, N and Liguori, G. 2017. Towards the prediction of multi-year to decadal climate variability in the Southern Hemisphere. PAGES Magazine, 25: 32–40. DOI: https://doi.org/10.22498/pages.25.1.32

41. Pelletier, JD. 1998. The power spectral density of atmospheric temperature from time scales of 10–2 to 106 yr. Earth and Planetary Sci. Lett., 158: 157–164. DOI: https://doi.org/10.1016/S0012-821X(98)00051-X

42. Proistosescu, C, Donohoe, A, Armour, KC, Roe, GH, Stuecker, MF and Bitz, CM. 2018. Radiative feedbacks from stochastic variability in surface temperature and radiative imbalance. Geophys. Res. Lett., 45: 5082–5094. DOI: https://doi.org/10.1029/2018GL077678

43. Rahmstorf, S. 2008. Anthropogenic climate change: Revisiting the facts. In: Zedillo, E (ed.), Global warming: looking beyond Kyoto. Washington, DC: Brookings Institution Press. p. 34–53.

44. Rogues, L, Chekroun, MD, Cristofol, M, Soubeyrand, S and Ghil, M. 2014. Parameter estimation for energy balance models with memory. Proc. R. Soc., A 470: 20140349. DOI: https://doi.org/10.1098/rspa.2014.0349

45. Rypdal, M and Rypdal, K. 2014. Long-memory effects in linear response models of Earth’s temperature and implications for future global warming. J. Clim., 27: 5240–5258. DOI: https://doi.org/10.1175/JCLI-D-13-00296.1

46. Saravanan, R and McWilliams, J. 1998. Advective ocean–atmosphere interaction: An analytical stochastic model with implications for decadal variability. J. Clim., 11: 165–188. DOI: https://doi.org/10.1175/1520-0442(1998)011<0165:AOAIAA>2.0.CO;2

47. Sellers, WD. 1969. A global climate model based on the energy balance of the Earth atmosphere system. J. Appl. Meteorol., 21: 391–400. DOI: https://doi.org/10.1175/1520-0450(1969)008<0392:AGCMBO>2.0.CO;2

48. Soldatenko, SA and Chichkine, D. 2016. Climate model sensitivity with respect to parameters and external forcing. In: Hromadka, T and Rao, P (eds.), Topics in Climate Modelling. Rijeka, Croatia: InTech Publ. pp. 105–135. DOI: https://doi.org/10.5772/64710

49. Soldatenko, S and Colman, R. 2019. Climate variability from annual to multi-decadal timescales in a two-layer stochastic energy balance model: analytic solutions and implications for general circulation models. Tellus A, 71: 1554421. DOI: https://doi.org/10.1080/16000870.2018.1554421

50. Stocker, T. 2011. Introduction to climate modelling. Berlin-Heidelberg: Springer-Verlag, p. 182. DOI: https://doi.org/10.1007/978-3-642-00773-6

51. Taylor, KE, Stouffer, RJ and Meehl, GA. 2012. An overview of the CMIP5 and the experimental design. Bull. American Meteorol. Soc., 93: 485–498. DOI: https://doi.org/10.1175/BAMS-D-11-00094.1

52. Xie, S-P, Kosaka, Y and Okumura, YM. 2016. Distinct energy budgets for anthropogenic and natural changes during global warming hiatus. Nature Geosci., 9: 29–34. DOI: https://doi.org/10.1038/ngeo2581

53. Yaglom, AM. 1987. Correlation Theory of Stationary and Related Random Functions, Volume: Basic results. New York: Springer-Verlag. p. 526. DOI: https://doi.org/10.1007/978-1-4612-4628-2

54. Zhu, F, Emile-Geay, J, McKay, NP, Hakim, GJ, Khider, D, Ault, TR, Steig, EJ, Dee, S and Kirchner, JW. 2019. Climate models can correctly simulate the continuum of global-average temperature variability. PNAS, 116: 8728–8733. DOI: https://doi.org/10.1073/pnas.1809959116