Journal cover Journal topic
Earth System Dynamics An interactive open-access journal of the European Geosciences Union
Journal topic
Earth Syst. Dynam., 10, 333-345, 2019
https://doi.org/10.5194/esd-10-333-2019
Earth Syst. Dynam., 10, 333-345, 2019
https://doi.org/10.5194/esd-10-333-2019

Research article 14 Jun 2019

Research article | 14 Jun 2019

# Including the efficacy of land ice changes in deriving climate sensitivity from paleodata

Efficacy of land ice changes and paleoclimate sensitivity
Lennert B. Stap, Peter Köhler, and Gerrit Lohmann Lennert B. Stap et al.
• Alfred-Wegener-Institut, Helmholtz-Zentrum für Polar- und Meeresforschung, Am Handelshafen 12, 27570 Bremerhaven, Germany
Abstract

The equilibrium climate sensitivity (ECS) of climate models is calculated as the equilibrium global mean surface air warming resulting from a simulated doubling of the atmospheric CO2 concentration. In these simulations, long-term processes in the climate system, such as land ice changes, are not incorporated. Hence, climate sensitivity derived from paleodata has to be compensated for these processes, when comparing it to the ECS of climate models. Several recent studies found that the impact these long-term processes have on global temperature cannot be quantified directly through the global radiative forcing they induce. This renders the prevailing approach of deconvoluting paleotemperatures through a partitioning based on radiative forcings inaccurate. Here, we therefore implement an efficacy factor ε[LI] that relates the impact of land ice changes on global temperature to that of CO2 changes in our calculation of climate sensitivity from paleodata. We apply our refined approach to a proxy-inferred paleoclimate dataset, using ${\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}={\mathrm{0.45}}_{-\mathrm{0.20}}^{+\mathrm{0.34}}$ based on a multi-model assemblage of simulated relative influences of land ice changes on the Last Glacial Maximum temperature anomaly. The implemented ε[LI] is smaller than unity, meaning that per unit of radiative, forcing the impact on global temperature is less strong for land ice changes than for CO2 changes. Consequently, our obtained ECS estimate of 5.8±1.3 K, where the uncertainty reflects the implemented range in ε[LI], is ∼50 % higher than when differences in efficacy are not considered.

1 Introduction

Equilibrium climate sensitivity (ECS) expresses the simulated equilibrated surface air temperature response to an instantaneous doubling of the atmospheric CO2 concentration. The simulated effect of the applied CO2 radiative forcing anomaly includes the Planck response, as well as the fast feedbacks such as those involving changes to snow, sea ice, lapse rate, clouds, and water vapour. ECS varies significantly between different state-of-the-art climate models; for instance, the CMIP5 ensemble shows a range of 1.9 to 4.4 K . Several ways have been put forward to constrain ECS, for example, through the use of paleoclimate data (e.g. Covey et al.1996; Edwards et al.2007), which is also the focus of this study. However, unlike results of models, temperature reconstructions based on paleoclimate proxy data always contain a mixed signal of all processes active in the climate system. Among these are long-term processes (or slow feedbacks) such as changes in vegetation, dust, and, arguably most importantly, land ice changes, which are kept constant in the climate model runs used to calculate ECS. Therefore, it is necessary to correct paleotemperature records for the influence of these processes, in order to make a meaningful comparison to ECS calculated by climate models.

In a coordinated community effort, the PALAEOSENS project proposed to relate the temperature response caused by these long-term processes to the globally averaged radiative forcing they induce . Consequently, the paleotemperature record can be disentangled on the basis of the separate radiative forcings of these long-term processes . If all processes are accounted for in this manner, the effect of CO2 changes and the accompanying short-term feedbacks, as described by the ECS, can be estimated. However, several studies have shown that, depending on the type of radiative forcing, the same global-average radiative forcing can lead to different global temperature changes . For instance, in a previous article , we simulated the separate and combined effects of CO2 changes and land ice changes on global surface air temperature using the intermediate-complexity climate model CLIMBER-2 and showed that the specific global temperature change per unit radiative forcing change depends on which process is involved. As a possible solution to this problem, formulated the concept of “efficacy” factors, which express the impact of radiative forcing by a certain process in comparison to the effect of radiative forcing by CO2 changes.

Based on the concept of , here we introduce an efficacy factor for radiative forcing by albedo changes due to land ice variability in our method of deriving climate sensitivity from paleodata. We first illustrate our refined approach by applying it to transient simulations over the past 5 Myr using CLIMBER-2 , obtaining a quantification of the effect on global temperature of CO2 changes and the accompanying short-term feedbacks from a simulation forced by both land ice and CO2 changes. We compare this result to a simulation where CO2 changes are the only operating long-term process. In this manner, we can assess the error resulting from using a constant efficacy factor. Thereafter, we refine a previous estimate of climate sensitivity based on a paleoclimate dataset of the past 800 kyr . In this dataset, the sole effect of CO2 is not a priori known. We therefore investigate the influence of the introduced efficacy factor on the calculated climate sensitivity. To do so, we appraise the influence of land ice changes and the associated efficacy using a range that is given by different modelling efforts of the Last Glacial Maximum (LGM; ∼21 kyr ago) (Shakun2017). The climate sensitivity resulting from applying this range provides a quantification of the consequence of the uncertain efficacy of land ice changes.

2 Material and methods

In this section, we first summarize the approach to obtaining climate sensitivity from paleodata that has been used in numerous earlier studies . We then discuss our main refinement to that approach, which is the inclusion of the efficacy of land ice changes, and a further small refinement that is meant to unify the dependent variable in cross-plots of radiative forcing and global temperature anomalies.

## 2.1 Approach to obtain climate sensitivity from paleodata

In climate model simulations used to quantify ECS, fast feedbacks, i.e. processes in the climate system with timescales of less than ∼100 years, are accounted for. However, slower processes, such as those involving changes to ice sheets, vegetation, and dust, are commonly kept constant. The resulting response is also sometimes called Charney sensitivity . Following the notation of , the ratio of the temperature change ($\mathrm{\Delta }{T}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$) to the radiative forcing due to the CO2 change ($\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$) yields Sa (in K W−1 m2, and where “a” stands for actuo):

$\begin{array}{}\text{(1)}& {S}^{\mathrm{a}}=\frac{\mathrm{\Delta }{T}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}.\end{array}$

The subscript denotes that CO2 is the only long-term process involved. Analogously, paleoclimate sensitivity (Sp) can be deduced from paleotemperature reconstructions and paleo-CO2 records as

$\begin{array}{}\text{(2)}& {S}^{\mathrm{p}}=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}.\end{array}$

In this case, the average global paleotemperature anomaly with respect to the pre-industrial period (PI) (ΔTg) is, however, also affected by the long-term processes that are typically neglected in climate simulations. Therefore, a correction to the paleotemperature record is needed to obtain $\mathrm{\Delta }{T}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ from ΔTg:

$\begin{array}{}\text{(3)}& \mathrm{\Delta }{T}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{\Delta }{T}_{\mathrm{g}}\left(\mathrm{1}-f\right),\end{array}$

or equivalently Sa from Sp:

$\begin{array}{}\text{(4)}& {S}^{\mathrm{a}}={S}^{\mathrm{p}}\left(\mathrm{1}-f\right)=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}\left(\mathrm{1}-f\right).\end{array}$

Here, f represents the effect of the slow feedbacks on paleotemperature (e.g. van de Wal et al.2011). To obtain f, proposed an approach, which has subsequently been used in numerous studies aiming to constrain climate sensitivity from paleodata and paleoclimate modelling studies . They suggested quantifying the influence of the long-term processes (X) by the radiative forcing change they induce (ΔR[X]), relative to the total radiative forcing perturbation:

$\begin{array}{}\text{(5)}& f=\frac{\mathrm{\Delta }{R}_{\left[X\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+\mathrm{\Delta }{R}_{\left[X\right]}}=\mathrm{1}-\frac{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+\mathrm{\Delta }{R}_{\left[X\right]}}.\end{array}$

Combining Eqs. (4) and (5) and following the PALAEOSENS nomenclature, we can then derive the “specific” paleoclimate sensitivity ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}$, where X represents the processes that are accounted for in the calculation of f:

$\begin{array}{ll}{S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}& =\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}\left(\mathrm{1}-\frac{\mathrm{\Delta }{R}_{\left[X\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+\mathrm{\Delta }{R}_{\left[X\right]}}\right)\\ \text{(6)}& & =\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+\mathrm{\Delta }{R}_{\left[X\right]}}=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}}.\end{array}$

If, for instance, only the most important slow feedback in the climate system, namely radiative forcing anomalies induced by albedo changes due to land ice (LI) variability, is taken into account, then one can correct Sp to derive the following specific climate sensitivity:

$\begin{array}{}\text{(7)}& {S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}}=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}}.\end{array}$

Using this approach, several studies performed a least-squares regression through scattered data from paleotemperature and radiative forcing records relating ΔTg to $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ in a time-independent manner, from which ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ could be determined. In the course of those studies, a state dependency of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ as function of background climate has been deduced for those data which are best approximated by a non-linear function. Furthermore, the quantification of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ for those state-dependent cases has been formalized in . A synthesis of estimates of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ from both colder- and warmer-than-present climates has been compiled by .

## 2.2 Refinement 1: taking the efficacy of land ice changes into account

The validity of the PALAEOSENS approach to calculating f is contingent on the notion that identical global-average radiative forcing changes lead to identical global temperature responses, regardless of the processes involved. However, it has been demonstrated that the horizontal and vertical distribution of the radiative forcing affects the resulting temperature response because, for example, different fast feedbacks are triggered depending on the location of the forcing. To address this issue, introduced the concept of efficacy factors, which we will explore further in this study. These factors (ε[X]) relate the strength of the temperature response to radiative forcing caused by a certain process X (ΔT[X]∕ΔR[X]) to a similar ratio caused by CO2 radiative forcing ($\mathrm{\Delta }{T}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}/\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$). This introduction of efficacy requires a reformulation of f as fε:

$\begin{array}{ll}{f}_{\mathit{\epsilon }}& =\frac{{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}}\\ \text{(8)}& & =\mathrm{1}-\frac{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}},\end{array}$

and hence also of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}$ as ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}^{\mathit{\epsilon }}$:

$\begin{array}{}\text{(9)}& {S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}^{\mathit{\epsilon }}=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}}.\end{array}$

In these reformulations, where in principle ε[X] can take any value, we introduce the superscript ε. This serves to clearly distinguish these newly derived sensitivities from those of the PALAEOSENS approach, in which efficacy was not taken into account, implying that identical radiative forcing of different processes leads to identical temperature changes.

To calculate ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$, we constrain the efficacy factor for radiative forcing by land ice changes (ε[LI]), using the following formulation, which is based on but slightly modified from :

$\begin{array}{}\text{(10)}& \frac{\mathrm{\Delta }{T}_{\left[\mathrm{LI}\right]}}{\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}}={\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}\frac{\mathrm{\Delta }{T}_{\mathrm{g}}-\mathrm{\Delta }{T}_{\left[\mathrm{LI}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}.\end{array}$

$\begin{array}{}\text{(11)}& {\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}=\frac{\mathit{\omega }}{\mathrm{1}-\mathit{\omega }}\frac{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}},\end{array}$

where ω represents the fractional relative influence of land ice changes on the global temperature change ($\mathit{\omega }=\mathrm{\Delta }{T}_{\left[\mathrm{LI}\right]}/\mathrm{\Delta }{T}_{\mathrm{g}}$). If ε[LI] is assumed to be constant in time (see Sect. 3.2), it can be calculated using Eq. (11) from data of any time step in the record of $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ and ΔR[LI], and consequently applied to the whole record (Fig. 1a, c). As before, with this ε[LI] a quantification of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ can be obtained by performing a least-squares regression through scattered data from paleotemperature and radiative forcing records, now relating ΔTg to $\left(\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}\right)$ in a time-independent manner.

Figure 1Time series of radiative forcing anomalies (ΔR) caused by CO2 (red) changes and land ice changes (blue) and global temperature anomalies (ΔTg) with respect to PI from (a–b) the CLIMBER-2 model dataset , with temperature data for experiment OIC in black and for experiment OC in green, and from (c–d) the proxy-inferred dataset , with solid lines for the whole dataset and dots for the data used in this study which exclude times with strong temperature–CO2 divergence (see Sect. 4.1). Note the differing axis scales.

Note that apart from the formulation based on followed here, other formulations of the efficacy factor are possible. For instance, one can define an alternative efficacy factor (ε[LI],alt) such that it relates the effect of land ice changes on global temperature directly to the radiative forcing anomaly caused by CO2 changes, leading to

$\begin{array}{}\text{(12)}& {S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right],\mathrm{alt}}^{\mathit{\epsilon }}=\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[\mathrm{LI}\right],\mathrm{alt}}\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}.\end{array}$

In this alternative case, the efficacy factor ε[LI],alt relates to our original ε[LI] as

$\begin{array}{}\text{(13)}& {\mathit{\epsilon }}_{\left[\mathrm{LI}\right],\mathrm{alt}}={\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}\frac{\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}.\end{array}$

This implies that if ε[LI] is indeed constant, any non-linearity in the relation between $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ and ΔR[LI] would demand a more complex formulation of the alternative efficacy factor ε[LI],alt (e.g. via a higher-order polynomial). Since we find such a non-linearity in our data (Fig. 2), using an F test to determine that a second-order polynomial is a significantly (p value < 0.0001) better fit to the data than a linear function, we refrain from following this alternative formulation further.

## 2.3 Refinement 2: unifying the dependent variable

To calculate ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$, previous studies have used cross-plots of global temperature anomalies and radiative forcing. The latter is caused by a combination of CO2 and land ice changes, which is cumbersome if one wants to compare ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ to other specific paleoclimate sensitivities ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}^{\mathit{\epsilon }}$, where more and/or different long-term processes are considered. Here, we therefore reformulate our quantification of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ to unify the dependent variable as $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$.

$\begin{array}{ll}{S}_{\left[{\mathrm{CO}}_{\mathrm{2}},X\right]}^{\mathit{\epsilon }}& =\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}}\\ & =\frac{\mathrm{\Delta }{T}_{\mathrm{g}}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}\frac{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}}\\ \text{(14)}& & =\frac{\mathrm{\Delta }{T}_{\left[-X\right]}^{\mathit{\epsilon }}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}\end{array}$

Here, $\mathrm{\Delta }{T}_{\left[-X\right]}^{\mathit{\epsilon }}$ is the global temperature change (with respect to PI) stripped of the inferred influence of processes X, defined as

$\begin{array}{}\text{(15)}& \mathrm{\Delta }{T}_{\left[-X\right]}^{\mathit{\epsilon }}:=\mathrm{\Delta }{T}_{\mathrm{g}}\frac{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[X\right]}\mathrm{\Delta }{R}_{\left[X\right]}}.\end{array}$

Hence, for the calculation of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$, we use

$\begin{array}{}\text{(16)}& \mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}:=\mathrm{\Delta }{T}_{\mathrm{g}}\frac{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}+{\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}}.\end{array}$

Now, we quantify ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ by performing a least-squares regression (regfunc) through scattered data from $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$ and $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$. We use the precondition that no change in CO2 is related to no change in $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$, meaning the regression intersects the y axis at the origin ($\left(x,y\right)=\left(\mathrm{0},\mathrm{0}\right)$). Following , for any non-zero $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$, we calculate ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ as

$\begin{array}{}\text{(17)}& {S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}{|}_{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}=\frac{\mathrm{regfunc}}{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}{|}_{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}}.\end{array}$

If $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{0}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$, as is, for instance, the case for pre-industrial conditions, ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ is quantified as

$\begin{array}{}\text{(18)}& {S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}{|}_{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{0}}=\frac{\mathit{\delta }\left(\mathrm{regfunc}\right)}{\mathit{\delta }\left(\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}\right)}{|}_{\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{0}}.\end{array}$

Equations (17) and (18) yield a quantification of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$, which can be compared to the value obtained for ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ using the PALAEOSENS approach that does not consider efficacy differences (equivalent to using ε[LI]=1) .

To obtain Sa, one needs to multiply ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ by a conversion factor $\mathit{\varphi }=\mathrm{0.64}±\mathrm{0.07}$ (1σ uncertainty) that accounts for the influence of other long-term processes, namely vegetation, aerosol, and non-CO2 greenhouse gas changes . Note that this multiplication by ϕ ignores any possible state dependencies in ϕ and assumes unit efficacy for processes other than land ice changes. Because a comprehensive analysis of the efficacy and state dependency of these other processes is beyond the scope of this study, it is a source of uncertainty to be investigated in future research. Finally, we obtain the equivalent ECS by multiplying Sa by $\mathrm{\Delta }{R}_{\mathrm{2}×{\mathrm{CO}}_{\mathrm{2}}}=\mathrm{3.71}±\mathrm{0.37}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$ (1σ uncertainty), the radiative forcing perturbation representing a CO2 doubling .

3 Illustration of the approach using model simulations

In this section, we illustrate our refined approach, which considers efficacy differences, by applying it to transient simulations over the past 5 Myr using CLIMBER-2 . We obtain a quantification of the effect on global temperature of CO2 changes and the accompanying short-term feedbacks from a simulation forced by both land ice and CO2 changes. We compare this result to a simulation where CO2 changes are the only operating long-term process. By doing so, we assess the error resulting from using a constant efficacy factor.

## 3.1 CLIMBER-2 model simulations

Currently, long (∼105 to ∼106 years) integrations of state-of-the-art climate models, such as general circulation models and Earth system models, are not yet not feasible due to limited computer power. This gap can be filled by using models of reduced complexity . Using the intermediate-complexity climate model CLIMBER-2 , climate simulations over the past 5 Myr were performed and analysed in . CLIMBER-2 combines a 2.5-dimensional statistical–dynamical atmosphere model, with a three-basin zonally averaged ocean model and a model that calculates dynamic vegetation cover based on the temperature and precipitation . The simulations could be forced by solar insolation changes due to orbital (O) variations , by land ice (I) changes in both hemispheres (based on de Boer et al.2013), and by CO2 (C) changes (based on van de Wal et al.2011). In the reference experiment (OIC) all these factors are varied, while in other model integrations the land ice (experiment OC) or the CO2 concentration (experiment OI) is kept fixed at PI level. The synergy of land ice and CO2 changes is negligibly small, meaning their induced temperature changes add approximately linearly when both forcings are applied. Furthermore, the influence of orbital variations is also very small, so that experiment OC approximately yields the sole effect of CO2 changes on global temperature (ΔT[OC]). As in , we use the simple energy balance model of to analyse the applied radiative forcing of land ice albedo and CO2 changes and simulated global temperature changes, after averaging to 1000-year temporal resolution (Fig. 1a, b).

## 3.2 Analysis

First, we analyse experiment OC, which will serve as a target for our approach as deployed later in this section. We use a least-squares regression through scattered data of $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ and ΔT[OC] to fit a second-order polynomial (Fig. 3a). Using a higher-order polynomial rather than a linear function allows us to capture the state dependency of paleoclimate sensitivity. Fitting even higher-order polynomials leads to negligible coefficients for the higher powers and is not pursued further. From the fit, we calculate a specific paleoclimate sensitivity ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ of 0.74 K W−1 m2 for PI conditions ($\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{0}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$) using Eq. (18). Note that, in this case, ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ is equal to ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}^{\mathit{\epsilon }}$, ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$, and ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ as there are no land ice changes and therefore also no efficacy differences. The fit further shows decreasing ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ for rising $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$.

Figure 2The relation between radiative forcing anomalies caused by CO2 changes ($\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$) and land ice changes (ΔR[LI]) from the whole proxy-inferred dataset (pink dots). The red line represents a second-order polynomial least-squares regression through the scattered data.

Now, we apply our approach to the results of experiment OIC, in which both CO2 and land ice cover vary over time, with the aim of deducing the sole effect of CO2 changes on global temperature. We calculate the efficacy of land ice changes for the LGM (21 kyr ago) from experiment OI, in which the CO2 concentration is kept constant. We obtain $\mathit{\omega }=\mathrm{\Delta }{T}_{\left[\mathrm{LI}\right]}/\mathrm{\Delta }{T}_{\mathrm{g}}=\mathrm{\Delta }{T}_{\left[\mathrm{OI}\right]}/\mathrm{\Delta }{T}_{\left[\mathrm{OIC}\right]}=\mathrm{0.54}$. Consequently, we find ε[LI]=0.58 from Eq. (11) and apply this value to the whole record of $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ and ΔR[LI]. In this manner, we calculate $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$ using Eq. (16). In principle, ε[LI] can be obtained using data from any time step of the record, preferably when the radiative forcing anomalies are large to prevent outliers resulting from divisions by small numbers. For example, using the results from all glacial marine isotope stages of the past 810 kyr (Marine Isotope Stage (MIS) 2, 6, 8, 10, 12, 14, 16, 18, and 20), instead of only the LGM, leads to a mean (±1σ) ε[LI] of 0.56±0.09.

We then fit a second-order polynomial to the scattered data of the thusly obtained $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$ from the results of experiment OIC and $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ (Fig. 3b, c). Between $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=-\mathrm{0.5}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$ and $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{0.5}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$, outliers resulted from division by small numbers (not shown in Fig. 3b). To remove these outliers, we first calculate the root mean square error (RMSE) between the fit and the data in the remainder of the domain. We then exclude all values from the range $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=-\mathrm{0.5}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$ to $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{0.5}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$, where the fit differs from the data by more than 3 times the RMSE, and perform the regression again. This yields an ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ of 0.72 K W−1 m2 for PI (Fig. 3b) in the LGM-only case and ${\mathrm{0.73}}_{-\mathrm{0.05}}^{+\mathrm{0.06}}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ in the case where all glacial periods are used (Fig. 3c). This supports our approach since it is only slightly lower than the ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ of 0.74 K W−1 m2 obtained from experiment OC, which it should approximate. The relationship between $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$ and $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ (Fig. 3b) is less linear than that between ΔT[OC] and $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ (Fig. 3a); hence, the state dependency of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ is enhanced. However, the difference between the ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ obtained from both experiments remains smaller than 0.07 K W−1 m2 throughout the entire 5 Myr interval in the LGM-only case, indicating that a constant efficacy is an acceptable assumption which only introduces a negligible additional uncertainty. However, the possible time dependency of efficacy could be investigated more rigorously in future research using more sophisticated climate models.

Figure 3Temperature anomalies with respect to PI over the last 5 Myr from CLIMBER-2 against imposed radiative forcing of CO2. (a) Simulation with fixed PI land ice distribution (experiment OC) (ΔT[OC]). (b) Calculated global temperature perturbations from experiment OIC stripped of the inferred influence of land ice ($\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$) using Eq. (16) with ε[LI]=0.58. Here, ε[LI] is obtained from matching climate sensitivity with the target value at the LGM. (c) Same as in panel (b) but using ε[LI]=0.47 (cyan dots), ε[LI]=0.56 (pink dots), and ε[LI]=0.65 (yellow dots), Here, ε[LI] is obtained from the mean (±1σ) of matching climate sensitivity with the target value at all glacial marine isotope stages of the past 810 kyr (MIS 2, 6, 8, 10, 12, 14, 16, 18, and 20). (d) Same as in panel (b) but using ε[LI]=1, which is equivalent to the PALAEOSENS approach where efficacy differences were not considered. The red lines – and in panel (c) also the orange and blue lines – represent second-order polynomial least-squares regressions through the scattered data.

The PALAEOSENS approach that does not consider efficacy differences (ε[LI]=1) yields a PI ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ of 0.54 K W−1 m2 (Fig. 3d). This is clearly much more off-target than the results of our approach, signifying the importance of considering efficacy.

4 Application to proxy-inferred paleoclimate data

In this section, we compare our refined approach to calculating ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ incorporating efficacy to our previous quantification of ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ by reanalysing the same paleoclimate dataset (introduced in Köhler et al.2015a). Other than for climate model simulations, in proxy-based datasets the influence of land ice changes on global temperature perturbations cannot be directly obtained and hence is a priori unknown. We therefore base the value of ε[LI] we implement here on a multi-model assemblage of simulated relative influences of land ice changes on the LGM temperature anomaly (Shakun2017).

## 4.1 Proxy-inferred paleoclimate dataset

The dataset to be investigated contains reconstructions of ΔTg, $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$, and ΔR[LI] for the past 800 kyr. Although the dataset covers the past 5 Myr, here we focus only on the past 800 kyr (Fig. 1c, d) because over this period $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ is constrained by high-fidelity measurements of CO2 within ice cores, whereas Pliocene and Early Pleistocene CO2 levels are still heavily debated . Radiative forcing by CO2 is obtained from Antarctic ice core data compiled by , using $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=\mathrm{5.35}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}\cdot \mathrm{ln}\left({\mathrm{CO}}_{\mathrm{2}}/\left(\mathrm{278}\phantom{\rule{0.125em}{0ex}}\mathrm{ppm}\right)\right)$ . The revised formulation for $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ from leads to very similar results with less than 0.01 W m−2 differences between the approaches for typical late Pleistocene CO2 values . Radiative forcing caused by land ice albedo changes, as well as the global surface air temperature record (ΔTg), are based on results of the 3-D ice-sheet model ANICE forced by northern hemispheric temperature anomalies with respect to a reference PI climate. The ANICE results are here considered to be proxy-inferred because, unlike climate models, ANICE is not constrained by climatic boundary conditions such as insolation and greenhouse gases. The temperature anomalies follow directly from a benthic δ18O stack using an inverse technique. Nevertheless, the results are model-dependent and therefore subject to uncertainty. ANICE provided geographically specific land ice distributions and hence radiative forcing due to albedo changes with respect to PI in both hemispheres. In , the northern hemispheric (NH) temperature anomalies (ΔTNH) are translated into global temperature perturbations (ΔTg1 in ) using polar amplification factors (${f}_{\mathrm{PA}}=\mathrm{\Delta }{T}_{\mathrm{NH}}/\mathrm{\Delta }{T}_{\mathrm{g}}$) as follows: at the LGM, fPA=2.7 is taken from the average of PMIP3 model data , while at the mid-Pliocene Warm Period (mPWP, about 3.2 Myr ago), fPA=1.6 is calculated from the average of PlioMIP results . At all other times, fPA is linearly varied as a function of the NH temperature. In Appendix A, we investigate the influence of the chosen polar amplification factor on our results. The temporal resolution of the dataset is 2000 years.

Analysing this dataset, found a temperature–CO2 divergence appearing mainly during, or in connection with, periods of decreasing obliquity related to land ice growth or sea level fall. For these periods, a significantly different ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ was obtained than for the remainder of the time frame. However, in the future, we expect sea level to rise; hence, these intervals of strong temperature–CO2 divergence should not be considered for the interpretation of paleodata in the context of future warming, e.g. by using paleodata to constrain ECS. In the following analysis, we therefore exclude these times with strong temperature–CO2 divergence, leaving 217 data points as indicated in Fig. 1c, d.

Figure 4The global temperature perturbations stripped of the inferred influence of land ice ($\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$) calculated using Eq. (16) against $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ from the proxy-inferred paleoclimate dataset , using (a) ε[LI]=0.79 (maroon dots), ε[LI]=0.45 (cyan dots), and ε[LI]=0.25 (green dots). Here, ε[LI] is obtained by converting the multi-model assemblage of simulated relative influences of land ice changes on the LGM temperature anomaly (0.46±0.14) (Shakun2017). (b) Same as in panel (a) but using ε[LI]=1 (grey dots), which is equivalent to the PALAEOSENS approach. The brown, blue, dark green (a), and black lines (b) represent second-order polynomial least-squares regressions through the data.

## 4.2 Analysis

compiled model-based estimates of the relative impact of land ice changes on the LGM temperature anomaly (ω in Eq. 11) using an ensemble of 12 climate models and estimated ω to be 0.46±0.14 (mean ±1σ; full range 0.200.68). Applying these values, in combination with the LGM values (taken here as the mean of the data at 20 and 22 kyr ago) $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=-\mathrm{2.04}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$ and $\mathrm{\Delta }{R}_{\left[\mathrm{LI}\right]}=-\mathrm{3.88}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$, yields ${\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}={\mathrm{0.45}}_{-\mathrm{0.20}}^{+\mathrm{0.34}}$. Implementing this range for ε[LI] in Eq. (16), we calculate $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$ over the whole 800 kyr period. Fitting second-order polynomials by least-squares regression to the scattered data of $\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$ and $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$, we infer a PI ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ of ${\mathrm{2.45}}_{-\mathrm{0.56}}^{+\mathrm{0.53}}\phantom{\rule{0.125em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ (Fig. 4a). The substantial uncertainty given here only reflects the 1σ uncertainty in ε[LI]. Similar to , we also detect a state dependency with decreasing ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ towards colder climates for this dataset, more strongly so in the case of lower ε[LI]. This state dependency is opposite to the one found in the CLIMBER-2 results (Sect. 3). The difference may be related either to the fact that fast climate feedbacks are too linear or that some slow feedbacks are underestimated in intermediate-complexity climate models like CLIMBER-2 (see Köhler et al.2018, for a detailed discussion). At $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}=-\mathrm{2.04}\phantom{\rule{0.125em}{0ex}}\mathrm{W}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{-\mathrm{2}}$, the LGM value, ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ is only ${\mathrm{1.45}}_{-\mathrm{0.37}}^{+\mathrm{0.33}}\phantom{\rule{0.125em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$. The PALAEOSENS approach, which does not consider efficacy and is therefore equivalent to our approach using ε[LI]=1, yields ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}=\mathrm{1.66}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ for PI and ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}=\mathrm{0.93}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ for the LGM (Fig. 4b). The specific paleoclimate sensitivities we find using the refined approach are hence generally larger than those obtained when neglecting efficacy differences. This is because, for the range of the impact of land ice changes on the LGM temperature anomaly implemented ($\mathit{\omega }=\mathrm{0.46}±\mathrm{0.14}$), the efficacy factor ε[LI] is smaller than unity. In other words, these land ice changes contribute comparatively less per unit radiative forcing to the global temperature anomalies than the CO2 changes.

Our inferred PI ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ is equivalent to an Sa of ${\mathrm{1.6}}_{-\mathrm{0.4}}^{+\mathrm{0.3}}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ when only considering the uncertainty caused by the implemented range in ε[LI] and to an Sa of ${\mathrm{1.6}}_{-\mathrm{0.2}}^{+\mathrm{0.1}}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ when only considering the uncertainty in the conversion factor ϕ. The equivalent ECS is 5.8±1.3 K per CO2 doubling when only considering the uncertainty caused by the implemented range in ε[LI] and 5.8±0.6 K per CO2 doubling when only considering the uncertainty in the conversion factor $\mathrm{\Delta }{R}_{\mathrm{2}×{\mathrm{CO}}_{\mathrm{2}}}$. The ECS we find is thus at the high end of the results of other approaches to obtain ECS , e.g. the 2.0 to 4.3 K 95 % confidence range from a large model ensemble and the 2.2 to 3.4 K 66 % confidence range from an emerging constraint from global temperature variability and CMIP5 . Hence, the low end of our ECS estimate is in the best agreement with these other estimates. This could mean that the relative influence of land ice changes on the LGM temperature anomaly is on the high side of or possibly higher than the 0.46±0.14 range we consider here. Alternatively, the conversion factor $\mathit{\varphi }=\mathrm{0.64}±\mathrm{0.07}$ we use to convert ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ to Sa is an overestimation, which could be caused by a larger-than-unity efficacy of long-term processes besides CO2 and land ice changes. We have focused primarily on the effect of ε[LI] on ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ in this analysis, and therefore we have, for simplicity, ignored uncertainties in the investigated proxy-inferred records themselves. A comprehensive description of these uncertainties and their influence on the calculated climate sensitivity can be found in .

5 Conclusions

We have incorporated the concept of a constant efficacy factor , which interrelates the global temperature responses to radiative forcing caused by land ice changes and CO2 changes, into our framework of calculating specific paleoclimate sensitivity ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$. The aim of this effort has been to overcome the problem that land ice and CO2 changes can lead to significantly different global temperature responses, even when they induce the same global-average radiative forcing. Firstly, we have assessed the usefulness of considering efficacy differences by applying our refined approach to results of 5 Myr CLIMBER-2 simulations , where the separate effects of land ice changes and CO2 changes can be isolated. In the results of these simulations, the error from assuming the efficacy factor to be constant in time is negligible. Thereafter, we have used our approach to reanalyse an 800 kyr proxy-inferred paleoclimate dataset . We have inferred a range in the land ice change efficacy factor ε[LI] from the relative impact of land ice changes on the LGM temperature anomaly simulated by a 12-member climate model ensemble (Shakun2017). The thusly obtained efficacy factor ${\mathit{\epsilon }}_{\left[\mathrm{LI}\right]}={\mathrm{0.45}}_{-\mathrm{0.20}}^{+\mathrm{0.34}}$ is smaller than unity, implying that the impact on global temperature per unit of radiative forcing is less strong for land ice changes than for CO2 changes. Consequently, our derived PI ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ of ${\mathrm{2.45}}_{-\mathrm{0.56}}^{+\mathrm{0.53}}\phantom{\rule{0.125em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ is ∼50 % larger than when efficacy differences are neglected. The equivalent Sa and ECS corresponding to this ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ are ${\mathrm{1.6}}_{-\mathrm{0.4}}^{+\mathrm{0.3}}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ and 5.8±1.3 K per CO2 doubling, respectively. The uncertainty in these estimates is only caused by the implemented range in ε[LI].

Data availability
Data availability.

Appendix A: Influence of the polar amplification factor

In the analysis performed in Sect. 4.2, we used a global temperature record that was obtained from northern high-latitude temperature anomalies using a polar amplification factor fPA that varies from 2.7 in the coldest to 1.6 in the warmest conditions (Sect. 4.1). However, recent climate model simulations of the Pliocene using updated paleogeographic boundary conditions show that in warmer times polar amplification could have been nearly the same as in colder times . We therefore repeat the analysis using the same range in ε[LI] and the same dataset but with an applied constant fPA=2.7 over the entire past 800 kyr to generate ΔTg (ΔTg2 in ).

Figure A1The global temperature perturbations stripped of the inferred influence of land ice ($\mathrm{\Delta }{T}_{\left[-\mathrm{LI}\right]}^{\mathit{\epsilon }}$) calculated using Eq. (16) against $\mathrm{\Delta }{R}_{\left[{\mathrm{CO}}_{\mathrm{2}}\right]}$ from the proxy-inferred paleoclimate dataset , using (a) ε[LI]=0.79 (maroon dots), ε[LI]=0.45 (cyan dots), and ε[LI]=0.25 (green dots). Here, ε[LI] is obtained from converting the multi-model assemblage of simulated relative influences of land ice changes on the LGM temperature anomaly (0.46±0.14) (Shakun2017). (b) Same as in panel (a) but using ε[LI]=1 (grey dots), which is equivalent to the PALAEOSENS approach. The brown, blue, dark green (a), and black lines (b) represent second-order polynomial least-squares regressions through the data. Here, the global temperature anomalies are derived from the northern high-latitude temperature anomaly reconstruction assuming a constant polar amplification factor (fPA) of 2.7, as opposed to the variable fPA used in Fig. 4.

The constant polar amplification used here counteracts increasing state dependency towards low temperatures, as the temperature differences are no longer amplified by changing polar amplification. Hence, ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ is smaller at PI (${\mathrm{1.96}}_{-\mathrm{0.44}}^{+\mathrm{0.42}}\phantom{\rule{0.125em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ compared to ${\mathrm{2.45}}_{-\mathrm{0.56}}^{+\mathrm{0.53}}\phantom{\rule{0.125em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ using the variable fPA) but diminishes less strongly towards colder conditions (Fig. A1a, cf. Fig. 4a). As before, the PALAEOSENS approach (equivalent to our approach using ε[LI]=1) yields a lower PI ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}$ of 1.34 K W−1 m2 (Fig. A1b). The PI ${S}_{\left[{\mathrm{CO}}_{\mathrm{2}},\mathrm{LI}\right]}^{\mathit{\epsilon }}$ inferred here using our refined approach corresponds to an Sa of ${\mathrm{1.3}}_{-\mathrm{0.3}}^{+\mathrm{0.2}}\phantom{\rule{0.25em}{0ex}}\mathrm{K}\phantom{\rule{0.125em}{0ex}}{\mathrm{W}}^{-\mathrm{1}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}$ and an ECS of ${\mathrm{4.6}}_{-\mathrm{1.3}}^{+\mathrm{1.0}}\phantom{\rule{0.125em}{0ex}}\mathrm{K}$ per CO2 doubling.

Author contributions
Author contributions.

LBS designed the research. LBS and PK performed the analysis. LBS drafted the paper, with input from all co-authors.

Competing interests
Competing interests.

The authors declare that they have no conflict of interest.

Acknowledgements
Acknowledgements.

This work is institutionally funded at AWI via the research program PACES-II of the Helmholtz Association. We thank Roderik van de Wal for commenting on an earlier draft of the paper and two anonymous referees for their constructive comments, which have helped to improve the quality of the paper.

Financial support
Financial support.

The article processing charges for this open-access publication were covered by a Research Centre of the Helmholtz Association.

Review statement
Review statement.

This paper was edited by Daniel Kirk-Davidoff and reviewed by two anonymous referees.

References

Badger, M. P. S., Lear, C. H., Pancost, R. D., Foster, G. L., Bailey, T. R., Leng, M. J., and Abels, H. A.: CO2 drawdown following the middle Miocene expansion of the Antarctic Ice Sheet, Paleoceanography, 28, 42–53, https://doi.org/10.1002/palo.20015, 2013. a

Bereiter, B., Eggleston, S., Schmitt, J., Nehrbass-Ahles, C., Stocker, T. F., Fischer, H., Kipfstuhl, S., and Chappellaz, J.: Revision of the EPICA Dome C CO2 record from 800 to 600 kyr before present, Geophys. Res. Lett., 42, 542–549, https://doi.org/10.1002/2014GL061957, 2015. a

Braconnot, P., Harrison, S. P., Kageyama, M., Bartlein, P. J., Masson-Delmotte, V., Abe-Ouchi, A., Otto-Bliesner, B., and Zhao, Y.: Evaluation of climate models using palaeoclimatic data, Nat. Clim. Change, 2, 417–424, https://doi.org/10.1038/nclimate1456, 2012. a

Brovkin, V., Ganopolski, A., and Svirezhev, Y.: A continuous climate-vegetation classification for use in climate-biosphere studies, Ecol. Model., 101, 251–261, https://doi.org/10.1016/S0304-3800(97)00049-5, 1997. a

Chalk, T. B., Hain, M. P, Foster, G. L., Rohling, E. J., Sexton, P. F., Badger, M. P. S., Cherry, S. G., Hasenfratz, A. P., Haug, G. H., Jaccard, S. L., Martínez-García, A., Pälike, H., Pancost, R. D., and Wilson, P. A.: Causes of ice age intensification across the Mid-Pleistocene Transition, P. Natl. Acad. Sci. USA, 114, 13114–13119, https://doi.org/10.1073/pnas.1702143114, 2017. a

Chandan, D. and Peltier, W. R.: Regional and global climate for the mid-Pliocene using the University of Toronto version of CCSM4 and PlioMIP2 boundary conditions, Clim. Past, 13, 919–942, https://doi.org/10.5194/cp-13-919-2017, 2017. a

Chandan, D. and Peltier, W. R.: On the mechanisms of warming the mid-Pliocene and the inference of a hierarchy of climate sensitivities with relevance to the understanding of climate futures, Clim. Past, 14, 825–856, https://doi.org/10.5194/cp-14-825-2018, 2018. a

Charney, J. G., Arakawa, A., Baker, D. J., Bolin, B., Dickinson, R. E., Goody, R. M., Leith, C. E., Stommel, H. M., and Wunsch, C. I.: Carbon dioxide and climate: a scientific assessment, National Academy of Sciences, Washington, DC, 1979. a

Claussen, M., Mysak, L., Weaver, A., Crucifix, M., Fichefet, T., Loutre, M.-F. ,Weber, S., Alcamo, J., Alexeev, V., Berger, A., Calov, R., Ganopolski, A., Goosse, H., Lohmann, G., Lunkeit, F., Mokhov, I., Petoukhov, V., Stone, P., and Wang, Z.: Earth system models of intermediate complexity: closing the gap in the spectrum of climate system models, Clim. Dynam., 18, 579–586, https://doi.org/10.1007/s00382-001-0200-1, 2002. a

Covey, C., Sloan, L. C., and Hoffert, M. I.: Paleoclimate data constraints on climate sensitivity: the paleocalibration method, Climatic Change, 32, 165–184, https://doi.org/10.1007/BF00143708, 1996. a

Cox, P. M., Huntingford, C., and Williamson, M. S.: Emergent constraint on equilibrium climate sensitivity from global temperature variability, Nature, 553, 319, https://doi.org/10.1038/nature25450, 2018. a

de Boer, B., van de Wal, R. S. W., Lourens, L. J., Bintanja, R., and Reerink, T. J.: A continuous simulation of global ice volume over the past 1 million years with 3-D ice-sheet models, Clim. Dynam., 41, 1365–1384, https://doi.org/10.1007/s00382-012-1562-2, 2013. a

de Boer, B., Lourens, L. J., and van de Wal, R. S. W.: Persistent 400,000-year variability of Antarctic ice volume and the carbon cycle is revealed throughout the Plio-Pleistocene, Nat. Commun., 5, 2999, https://doi.org/10.1038/ncomms3999, 2014. a

Dyez, K. A., Hönisch, B., and Schmidt, G. A.: Early Pleistocene obliquity-scale pCO2 variability at ∼1.5 million years ago, Paleoceanography and Paleoclimatology, 33, 1270–1291, https://doi.org/10.1029/2018PA003349, 2018. a

Edwards, T. L., Crucifix, M., and Harrison, S. P.: Using the past to constrain the future: how the palaeorecord can improve estimates of global warming, Prog. Phys. Geog., 31, 481–500, https://doi.org/10.1177/0309133307083295, 2007. a

Etminan, M., Myhre, G., Highwood, E. J., and Shine, K. P.: Radiative forcing of carbon dioxide, methane, and nitrous oxide: A significant revision of the methane radiative forcing, Geophys. Res. Lett., 43, 542–549, https://doi.org/10.1002/2016GL071930, 2016. a

Friedrich, T., Timmermann, A., Tigchelaar, M., Timm, O. E., and Ganopolski, A.: Non-linear climate sensitivity and its implications for future greenhouse warming, Science Advances, 2, e1501923, https://doi.org/10.1126/sciadv.1501923, 2016. a, b, c, d, e

Ganopolski, A., Petoukhov, V., Rahmstorf, S., Brovkin, V., Claussen, M., Eliseev, A., and Kubatzki, C.: CLIMBER-2: a climate system model of intermediate complexity. Part II: model sensitivity, Clim. Dynam., 17, 735–751, https://doi.org/10.1007/s003820000144, 2001. a

Goodwin, P., Katavouta, A., Roussenov, V. M., Foster, G. L., Rohling, E. J., and Williams, R. G.: Pathways to 1.5 C and 2 C warming based on observational and geological constraints, Nat. Geosci., 11, 102–107, https://doi.org/10.1038/s41561-017-0054-8, 2018. a

Hansen, J., Sato, M., Ruedy, R., Nazarenko, L., Lacis, A., Schmidt, G. A., Russell, G., Aleinov, I., Bauer, M., Bauer, S., Bell, N., Cairns, B., Canuto, V., Chandler, M., Cheng, Y., Del Genio, A., Faluvegi, G., Fleming, E., Friend, A., Hall, T., Jackman, C., Kelley, M., Kiang, N., Koch, D., Lean, J., Lerner, J., Lo, K., Menon, S., Miller, R., Minnis, P., Novakov, T., Oinas, V., Perlwitz, Ja., Perlwitz, Ju., Rind, D., Romanou, A., Shindell, D., Stone, P., Sun, S., Tausnev, N., Thresher, D., Wielicki, B., Wong, T., Yao, M., and Zhang, S.: Efficacy of climate forcings, J. Geophys. Res.-Atmos., 110, D18104, https://doi.org/10.1029/2005JD005776, 2005. a, b, c, d, e, f, g, h

Haywood, A. M., Hill, D. J., Dolan, A. M., Otto-Bliesner, B. L., Bragg, F., Chan, W.-L., Chandler, M. A., Contoux, C., Dowsett, H. J., Jost, A., Kamae, Y., Lohmann, G., Lunt, D. J., Abe-Ouchi, A., Pickering, S. J., Ramstein, G., Rosenbloom, N. A., Salzmann, U., Sohl, L., Stepanek, C., Ueda, H., Yan, Q., and Zhang, Z.: Large-scale features of Pliocene climate: results from the Pliocene Model Intercomparison Project, Clim. Past, 9, 191–209, https://doi.org/10.5194/cp-9-191-2013, 2013. a

Kamae, Y., Yoshida, K., and Ueda, H.: Sensitivity of Pliocene climate simulations in MRI-CGCM2.3 to respective boundary conditions, Clim. Past, 12, 1619–1634, https://doi.org/10.5194/cp-12-1619-2016, 2016. a

Knutti, R., Rugenstein, M. A. A., and Hegerl, G. C.: Beyond equilibrium climate sensitivity, Nat. Geosci., 10, 727, https://doi.org/10.1038/ngeo3017, 2017. a

Köhler, P., Bintanja, R., Fischer, H., Joos, F., Knutti, R., Lohmann, G., and Masson-Delmotte, V.: What caused Earth's temperature variations during the last 800,000 years? Data-based evidences on radiative forcing and constraints on climate sensitivity, Quaternary Sci. Rev., 29, 129–145, https://doi.org/10.1016/j.quascirev.2009.09.026, 2010. a

Köhler, P., de Boer, B., von der Heydt, A. S., Stap, L. B., and van de Wal, R. S. W.: On the state dependency of the equilibrium climate sensitivity during the last 5 million years, Clim. Past, 11, 1801–1823, https://doi.org/10.5194/cp-11-1801-2015, 2015a. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o

Köhler, P., de Boer, B., von der Heydt, A. S., Stap, L. B., and van de Wal, R. S. W.: Model-based changes in global annual mean surface temperature change (Delta T_g) and radiative forcing due to land ice albedo changes (Delta R_[LI]) over the last 5 Myr, supplementary material, data set, PANGAEA, https://doi.org/10.1594/PANGAEA.855449, 2015b. a

Köhler, P., Nehrbass-Ahles, C., Schmitt, J., Stocker, T. F., and Fischer, H.: A 156 kyr smoothed history of the atmospheric greenhouse gases CO2, CH4, and N2O and their radiative forcing, Earth Syst. Sci. Data, 9, 363–387, https://doi.org/10.5194/essd-9-363-2017, 2017a. a

Köhler, P., Stap, L. B., von der Heydt, A. S., de Boer, B., van de Wal, R. S. W., and Bloch-Johnson, J.: A state-dependent quantification of climate sensitivity based on paleo data of the last 2.1 million years, Paleoceanography, 32, 1102–1114, https://doi.org/10.1002/2017PA003190, 2017b. a, b, c, d, e, f

Köhler, P., Knorr, G., Stap, L. B., Ganopolski, A., de Boer, B., van de Wal, R. S. W., Barker, S., and Rüpke, L. H.: The effect of obliquity-driven changes on paleoclimate sensitivity during the late Pleistocene, Geophys. Res. Lett., 45, 6661–6671, https://doi.org/10.1029/2018GL077717, 2018. a, b, c, d, e, f, g, h, i, j

Laskar, J., Robutel, P., Joutel, F., Gastineau, M., Correia, A. C. M., and Levrard, B.: A long-term numerical solution for the insolation quantities of the Earth, Astron. Astrophys., 428, 261–285, https://doi.org/10.1051/0004-6361:20041335, 2004. a

Lisiecki, L. E. and Raymo, M. E.: A Pliocene-Pleistocene stack of 57 globally distributed benthic δ18O records, Paleoceanography, 20, PA1003, https://doi.org/10.1029/2004PA001071, 2005. a

Martínez-Botí, M. A., Foster, G. L., Chalk, T. B., Rohling, E. J., Sexton, P. F., Lunt, D. J., Pancost, R. D., Badger, M. P. S., and Schmidt, D. N.: Plio-Pleistocene climate sensitivity evaluated using high-resolution CO2 records, Nature, 518, 49–54, https://doi.org/10.1038/nature14145, 2015. a, b, c, d, e

Myhre, G., Highwood, E. J., Shine, K. P., and Stordal, F.: New estimates of radiative forcing due to well mixed greenhouse gases, Geophys. Res. Lett., 25, 2715–2718, https://doi.org/10.1029/98GL01908, 1998. a, b

PALAEOSENS Project Members: Making sense of palaeoclimate sensitivity, Nature, 491, 683–691, https://doi.org/10.1038/nature11574, 2012. a, b, c, d, e, f

Petoukhov, V., Ganopolski, A., Brovkin, V., Claussen, M., Eliseev, A., Kubatzki, C., and Rahmstorf, S.: CLIMBER-2: a climate system model of intermediate complexity. Part I: model description and performance for present climate, Clim. Dynam., 16, 1–17, https://doi.org/10.1007/PL00007919, 2000. a

Shakun, J. D.: Modest global-scale cooling despite extensive early Pleistocene ice sheets, Quaternary Sci. Rev., 165, 25–30, https://doi.org/10.1016/j.quascirev.2017.04.010, 2017. a, b, c, d, e, f

Stap, L. B., de Boer, B., Ziegler, M., Bintanja, R., Lourens, L. J., and van de Wal, R. S. W.: CO2 over the past 5 million years: Continuous simulation and new δ11B-based proxy data, Earth Planet. Sc. Lett., 439, 1–10, https://doi.org/10.1016/j.epsl.2016.01.022, 2016. a

Stap, L. B., van de Wal, R. S. W., de Boer, B., Bintanja, R., and Lourens, L. J.: The influence of ice sheets on temperature during the past 38 million years inferred from a one-dimensional ice sheet–climate model, Clim. Past, 13, 1243–1257, https://doi.org/10.5194/cp-13-1243-2017, 2017. a, b

Stap, L. B., Van de Wal, R. S. W., de Boer, B., Köhler, P., Hoencamp, J. H., Lohmann, G., Tuenter, E., and Lourens, L. J.: Modeled influence of land ice and CO2 on polar amplification and paleoclimate sensitivity during the past 5 million years, Paleoceanography and Paleoclimatology, 33, 381–394, https://doi.org/10.1002/2017PA003313, 2018a. a, b, c, d, e, f, g, h, i

Stap, L. B., van de Wal, R. S. W., de Boer, B., Köhler, P., Hoencamp, J. H., Lohmann, G., Tuenter, E., and Lourens, L. J.: Simulation of northern hemispheric, southern hemispheric and global temperature over the past 5 million years, data set, PANGAEA, https://doi.org/10.1594/PANGAEA.887427, 2018b. a

Stocker, T. F., Mysak, L. A., and Wright, D. G.: A zonally averaged, coupled ocean-atmosphere model for paleoclimate studies, J. Climate, 5, 773–797, https://doi.org/10.1175/1520-0442(1992)005<0773:AZACOA>2.0.CO;2, 1992. a

Stuber, N., Ponater, M., and Sausen, R.: Why radiative forcing might fail as a predictor of climate change, Clim. Dynam., 24, 497–510, https://doi.org/10.1007/s00382-004-0497-7, 2005. a, b

van de Wal, R. S. W., de Boer, B., Lourens, L. J., Köhler, P., and Bintanja, R.: Reconstruction of a continuous high-resolution CO2 record over the past 20 million years, Clim. Past, 7, 1459–1469, https://doi.org/10.5194/cp-7-1459-2011, 2011. a, b

Vial, J., Dufresne, J.-L., and Bony, S.: On the interpretation of inter-model spread in CMIP5 climate sensitivity estimates, Clim. Dynam., 41, 3339–3362, https://doi.org/10.1007/s00382-013-1725-9, 2013. a

von der Heydt, A. S., Köhler, P., van de Wal, R. S. W., and Dijkstra, H. A.: On the state dependency of fast feedback processes in (paleo) climate sensitivity, Geophys. Res. Lett., 41, 6484–6492, https://doi.org/10.1002/2014GL061121, 2014. a, b, c

von der Heydt, A. S., Dijkstra, H. A., van de Wal, R. S. W., Caballero, R., Crucifix, M., Foster, G. L., Huber, M., Köhler, P., Rohling, E., Valdes, P. J., Ashwin, P., Bathiany, S., Berends, C. J., van Bree, L. G. J., Ditlevsen, P., Ghil, M., Haywood, A. M., Katzav, J., Lohmann, G., Lohmann, J., Lucarini, V., Marzocchi, A., Pälike, H., Ruvalcaba Baroni, I., Simon, D., Sluijs, A., Stap, L. B., Tantet, A., Viebahn, J., and Ziegler, M.: Lessons on climate sensitivity from past climate changes, Current Climate Change Reports, 2, 148–158, https://doi.org/10.1007/s40641-016-0049-3, 2016.  a

Willeit, M., Ganopolski, A., Calov, R., Robinson, A., and Maslin, M.: The role of CO2 decline for the onset of Northern Hemisphere glaciation, Quaternary Sci. Rev., 119, 22–34, https://doi.org/10.1016/j.quascirev.2015.04.015, 2015. a

Yoshimori, M., Hargreaves, J. C., Annan, J. D., Yokohata, T., and Abe-Ouchi, A.: Dependency of feedbacks on forcing and climate state in physics parameter ensembles, J. Climate, 24, 6440–6455, https://doi.org/10.1175/2011JCLI3954.1, 2011. a, b