Articles | Volume 9, issue 2
Research article
25 Apr 2018
Research article |  | 25 Apr 2018

Earth system model simulations show different feedback strengths of the terrestrial carbon cycle under glacial and interglacial conditions

Markus Adloff, Christian H. Reick, and Martin Claussen

In simulations with the MPI Earth System Model, we study the feedback between the terrestrial carbon cycle and atmospheric CO2 concentrations under ice age and interglacial conditions. We find different sensitivities of terrestrial carbon storage to rising CO2 concentrations in the two settings. This result is obtained by comparing the transient response of the terrestrial carbon cycle to a fast and strong atmospheric CO2 concentration increase (roughly 900 ppm) in Coupled Climate Carbon Cycle Model Intercomparison Project (C4MIP)-type simulations starting from climates representing the Last Glacial Maximum (LGM) and pre-industrial times (PI). In this set-up we disentangle terrestrial contributions to the feedback from the carbon-concentration effect, acting biogeochemically via enhanced photosynthetic productivity when CO2 concentrations increase, and the carbon–climate effect, which affects the carbon cycle via greenhouse warming. We find that the carbon-concentration effect is larger under LGM than PI conditions because photosynthetic productivity is more sensitive when starting from the lower, glacial CO2 concentration and CO2 fertilization saturates later. This leads to a larger productivity increase in the LGM experiment. Concerning the carbon–climate effect, it is the PI experiment in which land carbon responds more sensitively to the warming under rising CO2 because at the already initially higher temperatures, tropical plant productivity deteriorates more strongly and extratropical carbon is respired more effectively. Consequently, land carbon losses increase faster in the PI than in the LGM case. Separating the carbon–climate and carbon-concentration effects, we find that they are almost additive for our model set-up; i.e. their synergy is small in the global sum of carbon changes. Together, the two effects result in an overall strength of the terrestrial carbon cycle feedback that is almost twice as large in the LGM experiment as in the PI experiment. For PI, ocean and land contributions to the total feedback are of similar size, while in the LGM case the terrestrial feedback is dominant.

1 Introduction

At the Last Glacial Maximum (21 000 years before present, referred to as LGM hereafter), global mean surface temperature was 4 to 5 C lower than today (Annan and Hargreaves2013). Vegetation was not only less widespread but primary productivity was also smaller (Prentice and Harrison2009). This was the consequence of the lower CO2 concentrations during that time (about 200 ppm less than today), acting physically via the resulting lower temperatures (greenhouse effect) and biogeochemically via the reduced photosynthetic activity due to less available CO2 in the atmosphere (reduced CO2 fertilization) (Prentice and Harrison2009). From measuring isotopic carbon composition in ocean sediment cores (Bird et al.1996) and the isotopic oxygen composition of air trapped in ice cores (Ciais et al.2012), it has been estimated that terrestrial carbon storage was several hundred gigatons less than today. This is consistent with less primary productivity, the effect of which on carbon storage must have been larger than the reduction in soil respiration by the lower temperatures (Prentice and Harrison2009). This describes how CO2 shaped the terrestrial carbon cycle at the LGM. But the terrestrial carbon cycle also has an effect on the atmospheric CO2 concentration. Hence, one may wonder whether the strength of this feedback was different from today in glacial times. This is what we investigate in the present paper by performing Earth system simulations for conditions of the Last Glacial Maximum and pre-industrial (PI) times. Indeed one could ask this question also for the oceanic carbon cycle component, but this paper focuses on the terrestrial component, which will be shown to dominate the difference in feedback strength between the two Earth system states.

To quantify the feedback between carbon cycle and climate, Friedlingstein et al. (2003) introduced two sensitivities characterizing the change in stored carbon (terrestrial and/or oceanic) due to different drivers: biogeochemical effects of changed atmospheric CO2 concentration, called the carbon-concentration effect measured by the β sensitivity (Pg C ppm−1), and climate change, called the carbon–climate effect measured by the γ sensitivity (Pg C K−1). For the recent climate, these sensitivities have been quantified in numerous Earth system simulations, especially within the international Coupled Climate Carbon Cycle Model Intercomparison Project (C4MIP) (Ciais et al.2013; Friedlingstein et al.2006). Attempts to quantify carbon cycle sensitivities for perturbations of climates from even earlier times are rare. The few observational studies relate reconstructions of atmospheric CO2 concentrations to reconstructions of temperature (Friedlingstein2015), but the resulting “observed” sensitivity estimates of atmospheric CO2 concentration to temperature typically involve the combined carbon-concentration and carbon–climate effect and thus neither measure β nor γ as defined by Friedlingstein et al. (2003). An exception is the study by Frank et al. (2010), who considered temperature and CO2 reconstructions for the last millennium before the industrial revolution: their estimate should be a good proxy for γ since during this period the changes in atmospheric CO2 concentration were only a few parts per million so that the carbon-concentration effect should be negligible. The resulting γ sensitivity turns out to vary in time showing values compatible with the low end of the range of values found in the C4MIP studies for the recent climate. Jungclaus et al. (2010) obtained similar values for γ from Earth system simulations of the last millennium. The compatibility of those γ values for the last millennium with those from the C4MIP for the recent climate may not be that surprising since the climates differ only moderately. On the other hand, the C4MIP values result from simulations that perturb the PI climate dramatically ( quadrupling of atmospheric CO2 concentration), while those for the last millennium are obtained from historical climate and CO2 variations (Frank et al.2010) that are rather moderate so that it is unclear what such a comparison of γ values actually means. To ensure comparability, in the present study we adopt the C4MIP methodology to determine carbon cycle sensitivities for past and recent times.

While there have been attempts to determine climate sensitivity for various climates of the deep past (PALEOSENSE2012), similar studies for carbon sensitivities are apparently missing. Nevertheless, for the climate during the LGM studied here, the underlying carbon-concentration and carbon–climate effects have been isolated in simulations to understand their separate importance for shaping the geographical distribution of vegetation as compared to today (Claussen et al.2013; Woillez et al.2011). While in these studies it was sufficient to simulate time slices for past and recent times, transient simulations are needed to determine carbon cycle sensitivities that could be compared to C4MIP values. In the present study we employ a fully coupled general circulation model including dynamic vegetation for transient simulations starting either from a climate state representing the LGM or from PI conditions and forced by a strong increase in atmospheric CO2. Letting the CO2 act either physically or biogeochemically, we isolate the individual contributions from the carbon-concentration and carbon–climate effects to changes in the terrestrial carbon budgets. Using this C4MIP-type experiment design, we quantify their contribution not only by computing β and γ for land carbon but also by performing a factor analysis following Stein and Alpert (1993) to investigate in particular the additivity of the two effects, which it is a precondition to obtain from those two sensitivities the feedback strength.

The paper is organized as follows. First we lay out the design of our simulation experiments. Next, in Sect. 3, we describe the mathematical framework used for our factor and feedback analysis. The analysis of the simulation results starts in Sect. 4 with a description of the two initial climate states representing the LGM and PI conditions (1850 AD). This is preparation for the analysis of the transient simulation in Sect. 5, which contains the main results of our investigation. By applying the factor and feedback analysis we demonstrate that the intensity of the considered feedback is very different for the Last Glacial Maximum and the recent climate and identify the underlying mechanisms explaining the observed differences in system behaviour. The paper concludes with a critical discussion of our results.

2 Experiment set-up

To quantify the feedback between carbon cycle and atmospheric CO2 concentrations, we combine the C4MIP experiment design Ciais et al. (2013) in the variant of concentration-driven simulations with a factor separation following Stein and Alpert (1993). Technically, we proceed by investigating the reaction in the climate and carbon cycle to a prescribed strong rise in atmospheric CO2. More precisely, we perform a set of four simulations called ctrl, clim, conc and full. While for the quantification of the feedbacks by the C4MIP approach, only three of these simulations are needed, by using the full set of all four simulations we are able to demonstrate that – in contrast to other models (Gregory et al.2009; Schwinger et al.2014; Zickfeld et al.2011) – the linearity assumption implicit to the C4MIP feedback analysis is indeed justified for our model. Starting from a control simulation (ctrl) performed at constant CO2 concentration, three transient simulations forced by rising CO2 concentrations are performed. In the first of those transient simulations (conc) only the carbon-concentration effect is active, which means that the rising CO2 concentration is “seen” only by the photosynthesis code of the model, while the radiation code constantly “sees” the CO2 value of the control simulation. Conversely, in the second transient simulation (clim) only the carbon–climate effect is active; i.e. only the radiation code sees the rising CO2 concentrations but not the photosynthesis model. In the third simulation (full) both effects are simultaneously active. These simulations are run once for LGM and once for PI conditions. In the following, we will use the term “experiment” to refer to one of the two cases LGM or PI. “Simulation” will refer to one of the four model runs ctrl, clim, conc or full.

Figure 1CO2 change scenarios as prescribed for the LGM and PI experiments: starting from 185 ppm (Last Glacial Maximum, green line) and starting from 285 ppm (pre-industrial, red line).


Figure 2Differences between the LGM and PI climates obtained in the respective ctrl simulations: (a) difference in global mean near-surface temperatures and (b) difference in plant water availability. Here the values in the LGM state are substracted from the values in the PI state. Land areas that are covered by ice in the LGM but not in the PI equilibrium state show soil humidity differences >0.4.


The CO2 concentrations for the ctrl simulations of the two experiments are 185 ppm (LGM) and 285 ppm (PI), which are also the initial conditions for the respective transient simulations. Experiments were performed with Earth-System Model of the Max Planck Institute (MPI-ESM; see below). In fact, we performed only the LGM experiment for this study since we could use published MPI-ESM Coupled Model Intercomparison Project Phase 5 (CMIP5) simulations (called piControl, esmFdbk1, esmFixClim1 and 1pctCo2) for our purpose that were performed for PI conditions with the same model version. The LGM simulations were initialized from restart files of the MPI-ESM CMIP5 Last Glacial Maximum spin-up experiment (1800 simulation years long), extended by another 200 years with dynamic vegetation now switched on. The PI simulations used for our study were initialized from a spin-up experiment covering more than 3000 years. For the transient simulations clim, conc and full, the same atmospheric CO2 concentration increase is imposed over a period of 150 years in both experiments (see Fig. 1), acting differently in the three simulations as explained above. The forcing for our LGM experiment is obtained by reducing the standard PI CO2 forcing by 100 ppm to account for lower glacial CO2 concentrations while preserving the rate of change. Because CO2 concentrations thereby increase by the same amount, the different reaction of the Earth system to the CO2 rise in the two experiments should mostly be attributable to the different initial conditions, i.e. the glacial–interglacial atmospheric CO2 offset and the particular characteristics of the initial climates. The distribution of ice sheets is prescribed for the appropriate LGM and PI conditions and is kept constant in all simulations.

The experiments are conducted with the MPI-ESM using the version described in Giorgetta (2013). The MPI-ESM consists of the atmosphere component ECHAM6 and the ocean component MPIOM, both including submodels for simulating the land and ocean carbon cycles. Because atmospheric CO2 concentrations are prescribed in our experiments, the oceanic and terrestrial carbon cycles are decoupled so that changes in the ocean carbon cycle are irrelevant for terrestrial carbon reservoirs that are the main interest here; nevertheless oceanic carbon fluxes play a role in calculating the overall carbon cycle feedback in our study and the physical ocean remains an important component of the climate dynamics affecting also the land carbon cycle. The land component JSBACH comprises the DYNVEG model for simulation of natural changes in the geographical distribution of vegetation controlled by competition and wind and fire disturbances (Reick et al.2013) and the BETHY model (Knorr2000) for simulation of the fast biochemical and biophysical processes of the biosphere, in particular photosynthetic production that is simulated following the Farquhar model (Farquhar et al.1980) for C3 plants and the Collatz model (Collatz et al.1992) for C4 plants. Vegetation is represented by eight plant functional types that differ in phenology and physiology and interact dynamically (Brovkin et al.2013). There is no anthropogenic land cover change considered in the experiments here. Terrestrial carbon dynamics are calculated with the CBALANCE carbon submodel of JSBACH (Reick et al.2010), representing vegetation, litter and soils by seven carbon pools, where the temperature dependence of heterotrophic respiration is modelled by a Q10 formula and turnover rates are in addition dependent on soil humidity. The oceanic biogeochemistry model HAMOCC (Ilyina et al.2013) calculates sea–air gas exchange, water column processes and sediment dynamics. CO2 exchange between sea and air is calculated with a temperature-dependent rate based on the thermodynamic disequilibrium at the interface. Carbon is then cycled as organically fixed carbon, dissolved inorganic carbon and calcium carbonate in the water column and sediments. Temperature-, nutrient- and light-dependent biological cycling of carbon within the water column is represented by an extended NPZD model (Six and Maier-Reimer1996); inorganic carbon cycling is based on Maier-Reimer and Hasselmann (1987), using updated chemical constants by Goyet and Poisson (1989).

3 Analysis framework

Here we introduce the mathematical framework for analysing our simulations in the next sections. First we describe how we apply the factor separation method by Stein and Alpert (1993) to separate the relative contributions of the carbon-concentration and carbon–climate effects to the overall changes in terrestrial carbon reservoirs. In the remainder of the section we describe the mathematical framework to disentangle the oceanic and atmospheric contributions to the overall carbon cycle feedback, as well as the contributions of those two effects to the feedback. This feedback framework was originally introduced by Friedlingstein et al. (2003) and further discussed by Gregory et al. (2009). We apply it here in the variant with prescribed atmospheric CO2 (Ciais et al.2013).

We apply the factor separation method of Stein and Alpert (1993) as follows. Let CL denote the total land carbon. The pure effects of the carbon-concentration and carbon–climate effects are individually quantified by the differences

(1) Δ C L,conc ( t ) := C L,conc ( t ) - C L,ctrl Δ C L,clim ( t ) := C L,clim ( t ) - C L,ctrl ,

where the indices of the right-hand side CL values refer to the simulations from which the values were obtained, while the indices of the ΔCL values at the left-hand side refer to the effect considered. The time dependence t appears only for the values from transient simulations but not for values from the control simulations which enter our calculations as mean values (indicated as a bar over the symbol). In addition, we quantify the “synergy” between the carbon-concentration and the carbon–climate effects, which is that part of the land carbon storage difference between the full and ctrl simulation that cannot be explained by a linear addition of the individual effects:

(2) Δ C L,syn ( t ) := ( C L,full ( t ) - C L,ctrl ) - ( Δ C L,conc ( t ) + Δ C L,clim ( t ) ) .

Note that in this way all separate factors sum up to the land carbon change in the full simulation:

(3) Δ C L,full ( t ) = Δ C L,conc ( t ) + Δ C L,clim ( t ) + Δ C L,syn ( t ) .

For the feedback analysis we consider the following differences in near-surface temperature and atmospheric CO2 concentration that develop in the transient simulations:

(4) Δ T clim ( t ) := T clim ( t ) - T ctrl Δ c c ( t ) := c c ( t ) - c c ctrl .

The concentration of atmospheric CO2 is denoted here by “cc” and measured in ppm CO2. Since cc(t) is the same for all transient simulations of a particular experiment, the index specifying the simulation has been omitted. With these definitions one can now introduce the two land carbon sensitivities

(5) β L ( t ) := Δ C L,conc ( t ) Δ c c ( t ) γ L ( t ) := Δ C L,clim ( t ) Δ T clim ( t ) .

βL (Pg C ppm−1) measures how strongly land carbon is affected in the conc simulation by changes in atmospheric CO2; since in the conc simulation only the carbon-concentration effect is active, βL measures the strength of this effect alone. Analogously, γL (Pg C K−1) measures how strongly land carbon is affected by temperature changes in the clim simulation; because in this simulation only the carbon–climate effect is active, it represents the strength of this effect alone. Similar sensitivities can be defined for ocean carbon but they will not be needed in this study.

In addition to βL and γL we will need the sensitivity of temperature to increasing CO2 concentrations in our simulations below, known as temperature sensitivity (K ppm−1) (Friedlingstein et al.2003):

(6) α ( t ) := Δ T clim ( t ) Δ c c ( t ) .

Note that in this framework α, βL and γL are time dependent – a point that will be further discussed below.

To introduce a measure for the feedback strength, the global carbon balance needs to be considered. Since the CO2 concentration is prescribed in our simulations, atmospheric carbon is not affected by ocean–atmosphere or land–atmosphere carbon fluxes; i.e. the global carbon budget is not closed. But one can diagnose how much external CO2 emissions into the atmosphere would be needed to close the global carbon budget. Considering our full simulation, the prescribed change in atmospheric carbon must match the imagined external carbon emissions Iext(t) minus the carbon uptake by ocean and land ΔCOL,full(t):

(7) Δ C A ( t ) = I ext ( t ) - Δ C OL,full ( t ) .

Assuming that ocean and land carbon uptake are proportional to the increase in atmospheric CO2, one can define the proportionality factor f(t) by

(8) Δ C OL,full ( t ) = : - f ( t ) Δ C A ( t ) ,

where the reason for introducing a minus sign here will become clear below. With this, one obtains from Eq. (7)

(9) Δ C A ( t ) = A ( t ) I ext ( t ) , with A ( t ) = 1 1 - f ( t ) .

A(t) is called the airborne fraction (Gregory et al.2009). If atmospheric carbon content would not be prescribed, A(t) would describe how much of the carbon Iext(t) added to the atmosphere would remain in it. Following Roe (2009), from the viewpoint of feedback analysis A(t) is the “gain” of the feedback: for A(t) larger/smaller than 1 the feedback is positive/negative; i.e. the forcing Iext(t) induces, via Eq. (8), additional carbon fluxes into/out of the atmosphere. In Eq. (9) the gain of the feedback is completely determined by the value of f(t), which – also following Roe (2009) – is called the “feedback factor”. Note that the sign in Eq. (8) is chosen such that a positive/negative feedback corresponds to a positive/negative sign of f(t).

In the present study we focus on the terrestrial contribution to the carbon cycle feedback. This contribution is obtained as follows. Splitting ΔCOL,full(t) in Eq. (8) into the separate contributions ΔCL,full(t) from land and ΔCO,full(t) from ocean, one can define individual land and ocean feedback factors

(10) Δ C O,full ( t ) = : - f O ( t ) Δ C A ( t ) Δ C L,full ( t ) = : - f L ( t ) Δ C A ( t )

so that

(11) f ( t ) = f O ( t ) + f L ( t ) .

Hence, the individual feedback factors from ocean and land contribute additively to the global feedback factor.

To disentangle the contributions of the carbon-concentration and the carbon–climate effect to fL(t), we assume that the synergy term in Eq. (3) is small compared to the others. Then one can express the carbon change in the full simulation induced by the combined action of the two effects by summing the carbon changes induced by the individual effects diagnosed in the simulations conc and clim. Using the definitions for α, βL and γL from above, and noting that atmospheric carbon content and atmospheric CO2 concentration are related via the conversion factor m=2.12Pg ppm−1 (Flato et al.2013), one thus finds

(12) f L ( t ) = - α ( t ) γ L ( t ) + β L ( t ) m .

Here the first term quantifies the contribution from the carbon–climate effect, while the second does so for the carbon-concentration effect.

Please note that the feedback considered here is different from that originally considered by Friedlingstein et al. (2003) or in the C4MIP study (Friedlingstein et al.2006): besides the fact that we focus on the feedbacks induced by terrestrial processes only, the more important difference to our study is that Friedlingstein et al. (2003) considered only the feedback induced by the carbon–climate effect (see Friedlingstein et al., 2003, Eq. 8b, Friedlingstein et al., 2006, Eq. 1, or Gregory et al., 2009, Eq. 17), while in our study, following Gregory et al. (2009), we quantify the feedback induced by the carbon–climate and carbon-concentration feedback together (see our Eq. 9). Please note also that there is confusion in the literature concerning the names gain and feedback factor; in our study we follow the naming convention of Roe (2009), who highlighted this confusion.

Figure 3Vegetation cover and carbon storage in the LGM and PI ctrl simulations. Vegetation cover is given as fraction of grid cell covered with vegetation, and carbon storage is given in kg C m−2.


4 Comparison of the simulated LGM and PI equilibrium states

Here we compare key climate and carbon variables from the LGM and PI ctrl simulations that are the initial states for the transient simulations analysed in the next sections. Globally, mean near-surface temperatures are 4.5 K colder in the LGM state than in the PI state, but locally temperatures differ by 20 K and more (see Fig. 2). Compared to PI, more water is available for vegetation growth in the LGM state, especially in the tropics and subtropics. This plant water availability is measured here in terms of the relative amount of water above the wilting point in the root zone of the soil, a value of 1 indicating optimal moisture levels and 0 indicating that photosynthesis is inhibited by water scarcity. Inland glaciers extend throughout most of North America and northern Europe in the LGM state, and the sea level is considerably lower, leading to a different geography, especially in the Bering Strait and the Malay Archipelago. On a global scale, less area is covered by vegetation and dense vegetation is restricted to the tropical zone (compare Fig. 3). In the PI state, vegetation reaches far more into the extratropics and the midlatitudes are more densely covered by vegetation. Terrestrial carbon reservoirs are larger in the PI experiment almost everywhere (see Fig. 3). Globally, terrestrial carbon reservoirs contain 1986 Pg C in the LGM and 3041 Pg C in the PI state. Our difference in carbon storage (1055 Pg C) matches the difference of 1030±625Pg C in non-permafrost land carbon obtained by Ciais et al. (2012) from combining model simulations with carbon and oxygen isotope data from sediment and ice cores; note that changes in permafrost carbon are not part of our simulations.

Figure 4Climatic changes in the full simulation (continuous lines), clim simulation (dashed lines) and the conc simulation (dotted lines) due to rising CO2 concentrations in the LGM experiment (red) and PI experiment (black). Panel (a) shows the globally averaged change in near-surface temperature and (b) that in plant water availability.


Figure 5Change in terrestrial carbon storage (Pg C) in the full simulations (black curves) and split into factors (coloured curves) as computed from Eqs. (1) and (2) for (a) the LGM experiment and (b) the PI experiment.


5 Reaction of the Earth system to rising CO2 concentration under different boundary conditions

The climate system reacts differently to rising CO2 concentrations under LGM and PI boundary conditions. Figure 4 shows changes in global mean near-surface temperature and plant water availability in the transient simulations. Due to rising CO2 concentrations, global mean near-surface temperature increases in the clim and full simulations while plant water availability decreases. Both of these changes are larger in the LGM experiment. The similarity of temperature changes in the clim and full simulations shows that the carbon-concentration and synergistic effects do not considerably affect global mean near-surface temperature. Nevertheless, also the carbon-concentration effect creates a small global warming towards the end of both experiments, as can be seen from the curves of the conc simulations. Gregory et al. (2009) explained this by less evapotranspiration under increased CO2 concentrations. The radiative effect of increased stomatal closure has been shown by previous studies, e.g. Doutriaux-Boucher et al. (2009). The influence of the carbon-concentration effect on other physical variables, however, is more important for the terrestrial carbon dynamics. For example, plant water availability, the second most important environmental constraint on most terrestrial carbon fluxes in the model, rises in the global average due to increased water use efficiency in connection with the carbon-concentration effect and decreases due to higher evapotranspiration losses under the higher temperatures as a consequence of the carbon–climate effect. Climate change dominates plant water availability changes in the full simulation, but a clear influence of the carbon-concentration effect and their synergies on plant water availability is also apparent.

Figure 6Sensitivities βL and γL to the carbon-concentration and carbon–climate effect, respectively, and temperature sensitivity α in the LGM (blue) and the PI experiment (red). Values are computed as a 20-year average around the indicated data point.


Figure 5 shows the change in terrestrial carbon storage in the transient simulations. Overall, the carbon-concentration effect increases terrestrial carbon storage in response to the rising CO2 concentration in both experiments (see the curves ΔCL,conc). This effect is stronger in the LGM than in the PI experiment. Carbon reservoir changes due to the carbon–climate effect are negative and of similar size in the two experiments (see curves ΔCL,clim). In both experiments, synergies of the two effects are small in the global integral (see curves ΔCL,syn). This shows that linear additivity of the carbon–climate and carbon-concentration effects can be assumed on the global scale for our experiments, even for the large climate perturbations considered here. This is important in the following because by this additivity one can separate the individual contributions of the two effects to the feedback strength by means of Eq. (12) (see the discussion there).

Figure 7Dependence of gross assimilation per square metre leaf area on ambient CO2 concentration at 20 C leaf temperature according to the implemented photosynthesis model (Farquhar et al.1980) for C3 plant physiology. Abbreviations stand for individual vegetation types: TET for tropical evergreen trees, TDT for tropical deciduous trees, EET for extratropical evergreen trees, EDT for extratropical deciduous trees, RGS for raingreen shrubs, DCS for deciduous shrubs and C3G for C3 grasses.


From Fig. 5 it becomes clear that the same absolute increase in atmospheric CO2 concentration triggers different reactions of terrestrial carbon storage in corresponding simulations of the LGM and PI experiments. This is also reflected in the terrestrial carbon cycle sensitivities as shown in Fig. 6 where the sensitivity values for the LGM and PI experiments are presented as a function of simulation time. In the following, before discussing the strength of the carbon cycle feedback, first the sensitivities and their temporal development will be studied separately.

5.1 The carbon-concentration effect

Initially βL increases in both experiments but the increase is steeper under glacial conditions. This stronger carbon-concentration effect in the LGM experiment is mostly due to the lower CO2 concentrations: in both experiments, photosynthesis is initially carboxylation rate limited. In other words, in both experiments the fraction of available radiative energy that the plants are able to use to build up organic matter is initially limited by low atmospheric CO2 concentrations. This initial CO2 limitation is lifted by increasing CO2 concentrations, which leads to increasing primary productivity that allows for the extension of vegetation and increasing terrestrial carbon storage. This mechanism becomes obvious from Fig. 7, which shows the dependence of the primary production rate on CO2 concentration calculated directly from the equations for C3 photosynthesis, which dominates global natural productivity, implemented in JSBACH. At low ambient CO2 concentrations, productivity increases steeply with rising CO2, but its sensitivity gets smaller at higher CO2 concentrations due to the convex nature of the underlying functionality. In our experiments the carbon-concentration effect on productivity differs most substantially in the tropics, where temperatures are similar but the lower LGM ambient CO2 concentration makes productivity more sensitive to CO2 increases in the glacial setting. Additionally, vegetation has more room to expand and can generally grow denser in the glacial tropics than under the drier pre-industrial conditions where tropical forests are more regularly perturbed by wild fires.

Figure 6a shows that, after 30 to 40 years, the increase of βL slows down and its values eventually start to decrease. Arora et al. (2013) attribute this behaviour to the different response time of primary production and biomass decomposition. While productivity increases almost instantaneously with rising CO2 concentration, biomass decomposition initially remains unchanged and increases only when after a temporal delay of the order of the lifetime of plants the additional carbon from higher plant productivity reaches the litter and soil carbon reservoirs. Additionally the carbon-concentration effect becomes less effective at high productivity levels because the carbon density of living vegetation is reaching upper limits. In fact, the amount of carbon allocatable to biomass carbon reservoirs is limited in JSBACH to account for a down-regulation of productivity in mature vegetation. But the sensitvity of productivity to ambient CO2 also changes: Fig. 7 shows a transition point from high to low dependence on CO2 changes. Below the transition point photosynthesis is carboxylation rate limited, while beyond the transition point it is limited by a lack of radiation (see any textbook on photosynthesis). Accordingly, as long as CO2 availability stays to be the main limitation for productivity, the carbon-concentration effect of rising CO2 concentration leads to large increases in productivity. In our experiments, the prescribed CO2 concentration rise is, however, large enough to reach a point where insolation becomes more limiting to productivity than CO2 availability. From that transition point on, the effectivity of the carbon-concentration effect saturates. In the PI experiment ambient CO2 concentration reaches that point of saturation earlier than in the LGM experiment, leading to a shorter period in the PI experiment where primary productivity is limited by CO2 availability and thus highly sensitive to rising CO2 concentrations.

5.2 The carbon–climate effect

The sensitivity γL grows increasingly negative in both experiments (see Fig. 6b) and increasingly larger in absolute value in the PI experiment than in the LGM experiment. Although γL values are clearly different in the two experiments, the overall terrestrial carbon reservoir changes in the clim simulations, from which the γL values are computed (see Eq. 5), are almost similar (compare Fig. 5). The reason for this is that the temperature sensitivity α also varies between the two experiments. Throughout the simulation α is larger in the LGM experiment. The higher temperature sensitivity and the lower carbon cycle sensitivity γL partially compensate for differences between the PI and LGM cases as is seen from Fig. 6d where the product αγL is plotted; it is this combination of sensitivities that determines the strength of the carbon–climate effect (compare Eq. 12). Thereby, the carbon–climate effect differs much less between the LGM and PI case than the carbon-concentration effect discussed above.

Table 1Sensitivity of net primary productivity (NPP) and soil respiration, Rh, to the carbon–climate effect. These sensitivities (ΔNPP∕ΔT and ΔRh∕ΔT) are computed from the clim simulation by first integrating NPP and Rh over the particular region (tropics, extratropics) and over the full simulation period and then dividing by the temperature change in this region. ΔRh∕ΔNPP is the quotient of the two sensitivities. “Tropics” refers here to the latitudinal belt between 30 S and 30 N and “extratropics” to the remaining part of the globe. Here, ΔNPP and ΔRh are considered positive for plant carbon uptake and soil carbon loss, respectively.

Download Print Version | Download XLSX

To understand the processes behind the different γL sensitivity in the two experiments, it is useful to analyse first how climate change induces carbon losses differently in the tropics and extratropics. Table 1 lists the change in soil respiration ΔRh and net primary productivity ΔNPP per degree temperature change as well as their ratio separately for tropics and extratropics in the two clim simulations. In both simulations this ratio is smaller than 1 in the tropics (carbon fluxes into land reservoirs change more than fluxes into the atmosphere) but larger than 1 in the extratropics (carbon fluxes into the atmosphere change more than fluxes into land carbon reservoirs), indicating a very different reaction of the carbon cycle under climate change in these two regions. In the tropics, net primary productivity and soil respiration decrease (see Table 1), indicating that living conditions deteriorate. This has two reasons: firstly, it gets drier so that plant productivity and also soil decomposition decrease. Secondly, the already hot tropical climate is getting even hotter so that physiological limitations are reached more frequently, deteriorating plant productivity by damaging the photosynthetic apparatus (implemented as “heat inhibition” in JSBACH). The reduction in NPP is much larger than the reduction in soil respiration; hence, in the tropics land carbon losses are mostly driven by reduced plant productivity. In the extratropics the situation is different: values of NPP and soil respiration (see Table 1) both rise under the warming climate because physiological processes speed up. But since ultimately soil respiration is fed from NPP, the considerably larger increase in soil respiration cannot be a result of the enhanced carbon input. The explanation, instead, is the enhanced decomposition of soil carbon that had accumulated in those vast cold boreal areas already in the control simulation from which the transient simulations are initialized. Hence, in the extratropics land carbon losses are mostly driven by the enhanced soil respiration of “old” carbon.

Having identified the major drivers for carbon losses in the tropics and extratropics, one can now understand why the sensitivity γL is larger in the PI than in the LGM experiment. In the tropics reduced plant productivity is the major driver, and productivity is more sensitive in the PI than the LGM experiment (see Table 1) because growth conditions deteriorate from already initially drier and hotter levels. In the extratropics the enhancement of soil respiration was found to be the major driver, and soil respiration reacts more sensitively in the PI than in the LGM experiment (see Table 1) because vegetation extends much farther north under the warmer conditions and in the absence of ice sheets, going along with vastly more extratropical old soil carbon. Hence, both in the tropics and in the extratropics the land carbon cycle is more sensitive to climate change in the PI experiment.

While our model set-up allows us to study the reaction of active carbon reservoirs to perturbations, it does not include inert carbon reservoirs which could be activated under a strong forcing (i.e. permafrost soils). This might be particularly important for the comparison of γL between the LGM and the PI state since Ciais et al. (2012) estimate that there was a considerably larger amount of inert carbon stored on land at the LGM than in the Holocene. Therefore, it has to be stressed that the sensitivities found in this study do only consider active carbon reservoirs.

Figure 8Feedback strength fL computed from Eq. (12) for the terrestrial carbon cycle in the LGM and the PI experiment.


5.3 Feedback strength of the terrestrial carbon cycle

The carbon–climate and the carbon-concentration effect cause a feedback of the terrestrial carbon cycle to rising atmospheric CO2 concentrations. The constantly negative values of the strength fL of this feedback (see Fig. 8) demonstrate that it dampens the effect of the forcing so that less carbon is left in the atmosphere than emitted. Accordingly, the feedback is negative in both experiments. From the beginning of the simulations, the feedback strength grows increasingly negative in both experiments, a trend that reverses later on with an earlier minimum in the PI experiment. This reflects the different development of βL that dominates the feedback strength for both PI and LGM (compare values of βL and αγL in Fig. 6). The dominance of βL is particularly visible towards the end of the simulations, where the timing of the reversal of the trends in fL match those in βL (compare Fig. 6). The constantly higher absolute values of fL in the LGM setting show that the feedback is much stronger under LGM conditions, especially towards the end of the simulations. Because βL is dominating fL, the stronger terrestrial LGM feedback is also explained by the mechanisms identified in Sect. 5.1 to cause the higher LGM sensitivity to the carbon-concentration effect.

6 Discussion and conclusion

In the present study we investigated in simulations how the terrestrial carbon cycle feedback differs between pre-industrial (PI) times and during the Last Glacial Maximum (LGM). This was done by separating the contributions from the carbon-concentration and carbon–climate effects that induce this feedback in C4MIP-type simulations. These simulations starting either at PI or LGM conditions are rather artificial, since the CO2 forcing scenario used to probe the feedbacks neither resembles the atmospheric CO2 changes during the Holocene nor is it realistic for recent times (compare Fig. 1). But they are not meant to be historically realistic. Instead, such artificial scenarios have been introduced to facilitate the comparison of the carbon cycle feedback across different models (Gregory et al.2009). In our study we adopted this approach for a comparison of this feedback between different climate states.

An important question for the applicability of the C4MIP-type feedback analysis is the additivity of the two effects for global land carbon storage because only then can the feedback strength be properly split into separate contributions from the two effects (Gregory et al.2009). Our factor separation analysis (Stein and Alpert1993) revealed that their synergy is rather small for both the PI and LGM case, meaning that we can indeed consider the two effects independently to understand the simulated feedback behaviour. Concerning this additivity models seem to behave differently: Gregory et al. (2009) reported significant deviations from additivity for the HadCM3LC model.

Table 2Terrestrial carbon sensitivities βL and γL, associated feedback factor fL, as well as the global feedback factor f that includes the oceanic feedback (see Eq. 11) from our simulations for PI and LGM, as well as their published CMIP5 model range for PI. Our values (columns PI and LGM) are taken as their value after 140 years of simulation. The CMIP5 model range is taken from Arora et al. (2013), considering only models without nitrogen cycle. The CMIP5 ranges for fL and f have been computed using the published CMIP5 sensitivities in Eq. (12) and its ocean analogue together with Eq. (11). Because the intermodel range for α is not given in Arora et al. (2013), αs were calculated from the gain g^E provided in Fig. 9 in Arora et al. (2013).

Download Print Version | Download XLSX

Generally, the values of the carbon sensitivities βL and γL are time dependent (compare Fig. 6), but for easier comparison they are usually reported taking their values at the end of the simulation period (Ciais et al.2013). The respective values from our simulations are given in Table 2, together with their CMIP5 intermodel range. Since we used the same data from MPI-ESM for the analysis of PI conditions that entered the CMIP5 study by Arora et al. (2013), one should expect that published values for βL and γL should be similar. This is indeed true for βL, for which we find 1.42 PgC ppm−1 while Arora et al. (2013) find 1.46 PgC ppm−1. But for γL we calculate 68.6 Pg C K−1 while they report 83.2 Pg C K−1. We attribute this apparent inconsistency to differences in the way we and Arora et al. (2013) compute sensitivities: while we use as a reference mean values from the control simulation (see Eqs. 1 and 4), we guess that Arora et al. (2013) use as a reference the value from the first year of the respective transient simulation. Thereby, the resulting sensitivity values are not only sensitive to random climate variations at the end of the simulations (which are typically smaller than changes from the strong forcing) but also sensitive to such variations at their beginning. For the considered sensitivities, this effect should be the largest for temperature that varies at much shorter timescales than carbon stocks. Accordingly, βL values should be less sensitive to the way they are computed and this may explain why our βL values are similar but γL values differ.

For our further considerations it is interesting to see how our LGM carbon sensitivities relate to published PI values. In view of the technical complications just mentioned, such a comparison makes sense only for βL. We see from Table 2 that our LGM βL is considerably larger than the PI value of any CMIP5 model. This may be taken as an indication that our result for differences in βL between PI and LGM is even robust against uncertainties in representing climate and carbon cycle in models. Since, as we discussed in Sect. 5.3, the terrestrial feedback strength, as measured by the feedback factor fL, is dominated by the contribution from βL (compare also Eq. 12), it is clear that for LGM and PI the feedback is dominated by the carbon-concentration effect. Hence, also the much larger LGM feedback factor fL – almost twice the PI value – should be a robust result from our study.

So far, we have concentrated our study on the terrestrial part of the Earth system, but it is interesting to consider for a moment also the oceanic contributions to the feedback to discuss the relevance of our results for the carbon cycle feedback in the Earth system as a whole. Our simulations have been performed also with the ocean carbon cycle being active. Accordingly, one can calculate from our simulations also the ocean feedback factor fO (see Eq. 10). A basic property of the global feedback strength is that ocean and land contributions to the overall feedback factor f are additive (compare Eq. 11). Obtaining in this way the global feedback strength, one sees from the values in Table 2 that in our simulations the terrestrial component dominates the global feedback in the LGM case, while both contributions are of approximately equal size for pre-industrial climate.

As discussed in Sect. 5.1 and 5.2, the difference in carbon sensitivities between the LGM and the PI experiments comes mostly from the different initial conditions of these experiments. But there is also a strong dependence on the strength of the CO2 forcing. For example, the difference in βL depends largely on whether the CO2 reaches values high enough to produce a switch from carboxylation-limited assimilation to radiation-limited assimilation. Additionally, bioclimatic limits of vegetation, model-specific maximum productivity rates, the choice of the global value for the wilting point and the assumed maximum vegetation density introduce limitations to the system that shape the behaviour of terrestrial carbon storage in the model. Such limitations should also exist in reality but are hard to quantify.

Besides the dependence on the forcing scenario, the calculated sensitivity parameters are also time dependent. This is due to the fact that the Earth system's response to the imposed forcing is not entirely instantaneous. Many physical and biogeochemical processes react on longer timescales (e.g. plant and ecosystem growth and inertia in heat and carbon reservoirs), which also interact and thereby complicate the system's response. This simultaneous dependence of the α, β and γ values on system state and forcing is well known (Arora et al.2013; Gregory et al.2009). Accordingly, these sensitivity metrics do not characterize an Earth system state as such but only a combination of initial state and forcing scenario. Hence, to isolate their state dependence one must consider simulations with similar forcing. This is the reason why in our study we subjected the LGM and PI state to the same increase in CO2.

To conclude, the present study has demonstrated that C4MIP-type simulations can be used to understand why the Earth system may react differently to rising CO2 concentrations under LGM and PI conditions. In the two experiments performed here for LGM and PI conditions, the terrestrial biosphere and associated land carbon dynamics show a clear, climate-state-dependent transient reaction to increasing CO2 concentrations. More precisely, under conditions of the Last Glacial Maximum, the terrestrial carbon flux balance is more sensitive to the carbon-concentration effect than under pre-industrial conditions. This is due to the lower CO2 concentration in the LGM initial state that allows for a larger productivity increase under CO2 concentration rise. The carbon–climate effect, in contrast, is larger under PI conditions, which is caused by higher initial temperatures and larger amounts of extratropical terrestrial carbon in the PI initial state. As a consequence of this behaviour, the terrestrial feedback is stronger for LGM than PI conditions.

Code availability

The model code is publicly available after registration at (MPI, 2018).

Data availability

Simulation data are available on request from the authors.

Author contributions

The study was led by MA who also performed the simulations and data analysis. All authors contributed to the design of the study and the manuscript.

Competing interests

The authors declare that they have no conflict of interest.


We would like to thank Irina Fast and the DKRZ team for their technical support and Gitta Lasslop for her careful comments on the final version of our draft. Additionally, we would like to thank our editor Christoph Heinze and the three anonymous reviewers as well as Peter Rayner for their constructive comments and suggestions. Their critical assessments helped us to significantly improve our manuscript, especially regarding the theoretical framework of the feedback analysis.

The article processing charges for this open-access
publication were covered by the Max Planck Society.

Edited by: Christoph Heinze
Reviewed by: Peter Rayner and three anonymous referees


Annan, J. D. and Hargreaves, J. C.: A new global reconstruction of temperature changes at the Last Glacial Maximum, Clim. Past, 9, 367–376,, 2013. a

Arora, V. K., Boer, G. J., Friedlingstein, P., Eby, M., Jones, C. D., Christian, J. R., Bonan, G., Bopp, L., Brovkin, V., Cadule, P., Hajima, T., Ilyina, T., Lindsay, K., Tjiputra, J. F., and Wu, T.: Carbon-Concentration and Carbon-Climate Feedbacks in CMIP5 Earth System Models, J. Climate, 26, 5289–5314, 2013. a, b, c, d, e, f, g, h, i

Bird, M., Lloyd, J., and Farquhar, G.: Terrestrial carbon storage from the last glacial maximum to the present, Chemosphere, 33, 1675–1685, 1996. a

Brovkin, V., Boysen, L., Raddatz, T., Gayler, V., Loew, A., and Claussen, M.: Evaluation of vegetation cover and land-surface albedo in MPI-ESM CMIP5 simulations, J. Adv. Model. Earth Sy., 5, 48–57,, 2013. a

Ciais, P., Tagliabue, A., Cuntz, M., Bopp, L., Scholze, M., Hoffmann, G., Lourantou, A., Harrison, S. P., Prentice, I., Kelley, D., Koven, C., and Piao, S.: Large inert carbon pool in the terrestrial biosphere during the Last Glacial Maximum, Nat. Geosci., 5, 74–79,, 2012. a, b, c

Ciais, P., Sabine, C., Bala, G., Bopp, L., Brovkin, V., Canadell, J., Chhabra, A., DeFries, R., Galloway, J., Heimann, M., Jones, C., Le Quéré, C., Myneni, R. B., Piao, S., and Thornton, P.: Carbon and Other Biogeochemical Cycles, in: Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Stocker, T., Qin, D., Plattner, G.-K., Tignor, M., Allen, S., Boschung, J., Nauels, A., Xia, Y., Bex, V., and Midgley, P., Cambridge University Press, Cambridge, UK, and New York, NY, USA, 465–570, 2013. a, b, c, d

Claussen, M., Selent, K., Brovkin, V., Raddatz, T., and Gayler, V.: Impact of CO2 and climate on Last Glacial maximum vegetation – a factor separation, Biogeosciences, 10, 3593–3604,, 2013. a

Collatz, G. J., Ribas-Carbo, M., and Berry, J. A.: Coupled photosynthesis-stomatal conductance model for leaves of C4 plants, Aust. J. Plant Physiol., 19, 519–538, 1992. a

Doutriaux-Boucher, M., Webb, M., Gregory, J. M., and Boucher, O.: Carbon dioxide induced stomatal closure increases radiative forcing via a rapid reduction in low cloud, Geophys. Res. Lett., 36, L02703,, 2009. a

Farquhar, G. D., Caemmerer, S., and Berry, J.: A biochemical model of photosynthesis in leaves of C3 species, Planta, 149, 78–90, 1980. a, b

Flato, G., Marotzke, J., Abiodun, B., Braconnot, P., Chou, S., 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.: Evaluation of Climate Models, in: Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Stocker, T., Qin, D., Plattner, G.-K., Tignor, M., Allen, S., Boschung, J., Nauels, A., Xia, Y., Bex, V., and Midgley, P., Cambridge University Press, Cambridge, UK, and New York, NY, USA, 741–866, 2013. a

Frank, D. C., Esper, J., Raible, C. C., Buntgen, U., Trouet, V., Stocker, B., and Fortunat, J.: Ensemble reconstruction constraints on the global carbon cycle sensitivity to climate, Nature, 463, 527–530, 2010. a, b

Friedlingstein, P.: Carbon cycle feedbacks and future climate change, Phil. Trans. R. Soc. A, 373, 2014042,, 2015. a

Friedlingstein, P., Dufresne, J.-L., Cox, P. M., and Rayner, P.: How positive is the feedback between climate change and the carbon cycle?, Tellus, 55B, 692–700, 2003. a, b, c, d, e, f, g, h

Friedlingstein, P., Cox, P., Betts, R., Bopp, L., von Bloh, W., Brovkin, V., Cadule, P., Doney, S., Eby, M., Fung, I., Bala, G., John, J., Jones, C. D., Joos, F., Kato, T., Kawamiya, M., Knorr, W., Lindsay, K., Metthews, H. D., Raddatz, T., Rayner, P., Reick, C., Roeckner, E., Schnitzler, K.-G., Schnur, R., Strassmann, K., Weaver, A. J., Yoshikawa, C., and Zeng, N.: Climate-Carbon Cycle Feedback Analysis: Results from the C4MIP Model Intercomparison, J. Climate, 19, 3337–3353, 2006. a, b, c, d

Giorgetta, M. A.: Climate and carbon cycle changes from 1850 to 2100 in MPI-ESM simulations for the Coupled Model Intercomparison Project phase 5, J. Adv. Model. Earth Sy., 5, 572–597, 2013. a

Goyet, C. and Poisson, A.: New determination of carbonic acid dissociation constants in seawater as a function of temperature and salinity, Deep-Sea Res. Pt I, 36, 1635–1654, 1989. a

Gregory, J. M., Jones, C. C., Cadule, P., and Friedlingstein, P.: Quantifying Carbon Cycle Feedbacks, J. Climate, 22, 5232–5250, 2009. a, b, c, d, e, f, g, h, i, j, k

Ilyina, T., Six, K. D., Segschneider, J., Maier-Reimer, E., Li, H., and Núñez-Riboni, I.: Global ocean biogeochemistry model HAMOCC: Model architecture and performance as component of the MPI-Earth system model in different CMIP5 experimental realizations, J. Adv. Model. Earth Sy., 5, 287–315, 2013. a

Jungclaus, J. H., Lorenz, S. J., Timmreck, C., Reick, C. H., Brovkin, V., Six, K., Segschneider, J., Giorgetta, M. A., Crowley, T. J., Pongratz, J., Krivova, N. A., Vieira, L. E., Solanki, S. K., Klocke, D., Botzet, M., Esch, M., Gayler, V., Haak, H., Raddatz, T. J., Roeckner, E., Schnur, R., Widmann, H., Claussen, M., Stevens, B., and Marotzke, J.: Climate and carbon-cycle variability over the last millennium, Clim. Past, 6, 723–737,, 2010.  a, b, c

Knorr, W.: Annual and interannual CO2 exchanges of the terrestrial biosphere: Process-based simulations and uncertainties, Global Ecol. Biogeogr., 9, 225–252, 2000. a

Maier-Reimer, E. and Hasselmann, K.: Transport and storage of CO2 in the ocean – an inorganic ocean-circulation carbon cycle model, Clim. Dynam., 2, 63–90, 1987. a

MPI (Max-Planck-Institut für Meteorologie): Model distribution, available at:, last access: 17 April 2018. 

PALEOSENSE: Making sense of palaeoclimate sensitivity, Nature, 491, 683–691, 2012. a

Prentice, I. C. and Harrison, S. P.: Ecosystem effects of CO2 concentration: evidence from past climates, Clim. Past, 5, 297–307,, 2009. a, b, c

Reick, C., Raddatz, T., Pongratz, J., and Claussen, M.: Contribution of anthropogenic land cover change emissions to pre-industrial atmospheric CO2, Tellus B, 62, 329–336, 2010. a

Reick, C. H., Raddatz, T., Brovkin, V., and Gayler, V.: The representation of natural and anthropogenic land cover change in MPI-ESM, J. Geophys. Res., 5, 1–24, 2013. a

Roe, G.: Feedbacks, timescales, and seeing red, Annu. Rev. Earth Pl. Sc., 37, 93–115, 2009. a, b, c

Schwinger, J., Tjiputra, J. F., Heinze, C., Bopp, L., Christian, J. R., Gehlen, M., Ilyina, T., Jones, C. D., Salas-Mélia, D., Segschneider, J., Séférian, R., and Totterdell, I.: Nonlinearity of Ocean Carbon Cycle Feedbacks in CMIP5 Earth System Models, J. Climate, 27, 3869–3888, 2014. a

Six, K. D. and Maier-Reimer, E.: Effects of plankton dynamics on seasonal carbon fluxes in an ocean general circulation model, Global Biogeochem. Cy., 10, 559–583, 1996. a

Stein, U. and Alpert, P.: Factor Separation in Numerical Simulations, J. Atmos. Sci., 50, 2107–2115, 1993. a, b, c, d, e

Woillez, M.-N., Kageyama, M., Krinner, G., de Noblet-Ducoudré, N., Viovy, N., and Mancip, M.: Impact of CO2 and climate on the Last Glacial Maximum vegetation: results from the ORCHIDEE/IPSL models, Clim. Past, 7, 557–577,, 2011. a

Zickfeld, K., Eby, M., Matthews, H. D., Schmittner, A., and Weaver, A. J.: Nonlinearity of carbon cycle feedbacks, J. Climate, 24, 4255–4275, 2011. a

Short summary
Computer simulations show that during an ice age a strong atmospheric CO2 increase would have resulted in stronger carbon uptake of the continents than today. Causes are the larger potential of glacial vegetation to increase its photosynthetic efficiency under increasing CO2 and the smaller amount of carbon in extratropical soils during an ice age that can be released under greenhouse warming. Hence, for different climates the Earth system is differently sensitive to carbon cycle perturbations.
Final-revised paper