Disentangling residence time and temperature sensitivity of microbial decomposition in a global soil carbon model

Recent studies have identified the first-order rep- resentation of microbial decomposition as a major source of uncertainty in simulations and projections of the terrestrial carbon balance. Here, we use a reduced complexity model representative of current state-of-the-art models of soil or- ganic carbon decomposition. We undertake a systematic sen- sitivity analysis to disentangle the effect of the time-invariant baseline residence time (k) and the sensitivity of micro- bial decomposition to temperature (Q10) on soil carbon dy- namics at regional and global scales. Our simulations pro- duce a range in total soil carbon at equilibrium of 592 to 2745 Pg C, which is similar to the 561 to 2938 Pg C range in pre-industrial soil carbon in models used in the fifth phase of the Coupled Model Intercomparison Project (CMIP5). This range depends primarily on the value ofk, although the impact of Q10 is not trivial at regional scales. As climate changes through the historical period, and into the future, k is primarily responsible for the magnitude of the response in soil carbon, whereas Q10 determines whether the soil re- mains a sink, or becomes a source in the future mostly by its effect on mid-latitude carbon balance. If we restrict our sim- ulations to those simulating total soil carbon stocks consis- tent with observations of current stocks, the projected range in total soil carbon change is reduced by 42 % for the histori- cal simulations and 45 % for the future projections. However, while this observation-based selection dismisses outliers, it does not increase confidence in the future sign of the soil carbon feedback. We conclude that despite this result, future estimates of soil carbon and how soil carbon responds to cli- mate change should be more constrained by available data sets of carbon stocks.


Introduction
There is a 6-fold range in the amount of carbon stored in the soil in simulations conducted as part of the fifth phase of the Coupled Model Intercomparison Project (CMIP5; Taylor et al., 2012). This 6-fold range, identified by Todd-Brown et al. (2013), is consistent with results from model intercomparison projects such as the Coupled Climate-Carbon Cycle Model Intercomparison Project (C 4 MIP; Friedlingstein et al., 2006). The analysis of carbon stores in both C 4 MIP and CMIP5 has focused on the prediction of terrestrial and soil carbon through time. In addition to demonstrating the large differences in carbon stocks (Todd-Brown et al., 2013), they have also highlighted large inter-model differences in global and regional land-atmosphere carbon (C) fluxes (e.g. Friedlingstein et al., 2006Friedlingstein et al., , 2014. This lack of agreement between simulations not only exists in fully coupled models (e.g. C 4 MIP and CMIP5) but also can be found if sources of uncertainty are narrowed by relying on one weather data set to drive multiple land models Nishina et al., 2014), or by using one land model driven by multiple climate projections (Ahlström et al., 2013).
In these previous studies, critical uncertainties have been identified in the microbial decomposition of soil organic C and the associated release of CO 2 via heterotrophic respiration (R h ). This is despite all the current state-of-the-art global soil C models relying on a similar representation of decomposition as a first-order process (see Exbrayat et al., 2013b;Nishina et al., 2014;Todd-Brown et al., 2013). This conceptualisation describes decomposition and R h as proportional to the availability of organic matter. The decay rate (or R h per unit of soil C) is modified based on an environmental scalar that intends to mimic the dynamical response of microbial biomass to soil moisture and soil temperature.
This simple model structure has recently received some criticism because it lacks explicit representation of microbial physiology (Allison et al., 2010;Todd-Brown et al., 2012;Wieder et al., 2013;Xenakis and Williams, 2014). Furthermore, the formulation of the environmental scalar is held constant in time which is not consistent with recently identified enhancing or compensatory responses of microbial communities to changes in boundary conditions (Karhu et al., 2014). Therefore, it can only explain the acclimation of decomposers to warming (Luo et al., 2001) as a result of the quick depletion of labile pools by enhanced microbial biomass (Kirschbaum, 2004;Knorr et al., 2005).
We previously identified (Exbrayat et al., 2013b) some further implications of the first-order representation of microbial decomposition. First, in climate change experiments, model pools are usually initialised using a spin-up procedure with fixed pre-industrial atmospheric CO 2 concentrations until C pool trends are removed (Xia et al., 2012). Due to the interaction with substrate availability, the decay rate simulated by the model in response to steady boundary conditions determines the size of soil C pools reached at equilibrium. Because spin-up is a long computational process, the magnitude of pool sizes is conserved during subsequent shorter simulations of climate change and, as a result, equilibrated stocks strongly explain final stocks (e.g. CMIP5 models as shown in Fig. S1 in the Supplement and in Exbrayat et al., 2014). Second, the microbial sensitivity to changing environmental conditions affects the response of the system under transient climate simulations (Falloon et al., 2011;Exbrayat et al., 2013a, b). However, because substrate availability also controls the amount of respired carbon, there is a memory control imposed by the initial conditions of this transient simulation (Exbrayat et al., 2013b and that also affects the response to perturbation in boundary conditions. The relative contribution of these two factors on soil C projections remains to be explored in detail especially since last generation models disagree on the carbon balance projected in the future (Friedlingstein et al., 2014;, making it challenging to elaborate any land-based offsetting strategy.
Here, we use a reduced complexity model representative of current state-of-the-art models of soil organic C decomposition. A systematic sensitivity analysis is performed to disentangle the effect of the time-invariant baseline residence time and the formulation of the dynamic response of microbial decomposition to climatic change on soil C dynamics at regional and global scale. Using these experiments, we seek to investigate the relative contribution of these two inter-related components that drive the absolute and relative change in soil C through time. This is a step towards understanding the origin of the disagreement between CMIP5 models' simulation of soil C and can help in reducing the uncertainty in future model intercomparisons. We also use available estimates of total soil C to assess the added value of observational data to inform the modelling procedure. We attempt to constrain the system's response to climate change by identifying model versions that simulate amounts of soil C mobilised in the active cycle that are outside the confidence intervals estimated for the observations. We argue that, due to the first-order parameterisation, such model versions are unlikely to provide reliable projections of the response of soil C pools, as they would do it for the wrong reasons. We do not aim to provide new estimates of SOC response to climate change with our reduced complexity model. Instead, we suggest that our results will help the CMIP6 community to improve the design of future intercomparisons by highlighting the need and benefits of confronting models with existing data to reduce the uncertainty.

Reduced complexity model
It is not possible to re-run each CMIP5 model or isolate the representation of soil carbon processes from each model. This would be extraordinarily computationally expensive and the associated feedbacks would make the analysis of the results problematic. A far simpler approach is required which led Todd-Brown et al. (2013 to demonstrate that the CMIP5 SOC dynamics can be successfully reproduced using a simplified model structure. In this paper we develop and then use a reduced complexity model that simulates the monthly evolution of a single soil organic carbon pool, C s , in response to input derived from net primary productivity (NPP, g C m −2 mth −1 ) and output by heterotrophic respiration (R h , g C m −2 mth −1 ). For each monthly time step, the soil carbon balance can be described as where NPP is a prescribed boundary condition in our model and R h is simulated as a first-order process dependent on the availability of substrate C s such as where k is the baseline residence time at 15 • C (Xia et al., 2013) adjusted at each time step by f T which is a function of soil temperature T s ( • C). The soil moisture (θ s ) modification function, f W , is usually expressed as a fraction of soil moisture saturation (Moyano et al., 2012). We implement a classical formulation of the soil temperature sensitivity function f T : where Q 10 is a constant factor that describes the relative increase in microbial activity for a warming of 10 • C, and T ref is the reference temperature ( • C) for which f T (T s ) = 1 (Lloyd and Taylor, 1994;Bauer et al., 2012). The chosen T ref is the commonly used 15 • C (Todd-Brown et al., 2013) so that the decomposition rate equals k −1 when moisture is non-limiting and temperature is approximately equal to the global average. We use the same formulation of f W as in the Carnegie-Ames-Stanford approach of carbon, nitrogen, and phosphorus cycles (CASA-CNP) model (Wang et al., 2010): which is a bell-shaped function that is equal to 1 for θ s = 0.55.
This first-order representation of microbial decomposition with a specified decay rate adjusted by environmental scalars is used in all 11 CMIP5 models that simulate soil carbon (Todd-Brown et al., 2013) and all 7 dynamic global vegetation models (DGVMs) used in the Inter-Sectoral Impact Model Intercomparison Project (ISI-MIP) Nishina et al., 2014). Typically, these models rely on a multi-pool architecture to represent the diversity in organic matter. Each pool has its own residence time that corresponds to a degree of resistance to decomposition (Davidson and Janssen, 2006). Usually, part of the decomposition occurring in one pool is routed to one or several other pools while the rest is emitted via R h . At the ecosystem scale, however, the same environmental scalar is applied despite the multipool architecture, and the heterotrophic respiration flux is proportional to the amount of substrate available. Therefore, our simplified model is broadly representative of the current paradigm and provides a useful framework to undertake the sensitivity analysis described hereafter.
Soil moisture also has an influence on microbial decomposition (Falloon et al., 2011, Moyano et al., 2012Exbrayat et al., 2013a, b). However, Todd-Brown et al. (2013 recently demonstrated that a one pool reduced complexity model could reproduce both total soil carbon content and its spatial distribution for most of the CMIP5 models without considering decomposition response to variations in soil moisture. We also recently showed that global features in the distribution and evolution of C s were much more related to uncertainties in f T than uncertainties in the formulation of f W (Exbrayat et al., 2013b). Therefore, in order to keep the analyses as simple as possible and isolate the effect of f T but still account for the effect of soil moisture on R h , we keep the formulation of f W constant in the experiments that follow.
We are aware that our reduced complexity model relies on simplifications such as the use of a single soil carbon pool and global values of k, Q 10 and T ref . While a multiple pool structure would provide diverging results, single pool soil carbon models similar to our design are used in 3 of the 11 CMIP5 models described by Todd-Brown et al. (2013) and two of the seven ISI-MIP models described by Nishina et al. (2014). Further, using global parameter values of k, Q 10 and T ref are consistent with these state-of-theart models (Todd-Brown et al., 2013;Nishina et al., 2014). Of course, this does not allow representing processes such as the re-mobilisation of carbon in the active cycle following permafrost thaw (Koven et al., 2011) or the probably different behaviour of biological systems in frozen conditions but these are not implemented in the land component of CMIP5 Earth system models and therefore fall beyond the scope of this paper. In summary, we fully appreciate that our reduced complexity model is a simplification of the processes that operate in various regions of the Earth system. However, we note that our study investigates the sensitivity of the firstorder parameterisation of microbial decomposition and R h processes used in current ecosystem models to its uncertain parameters (Todd-Brown et al., 2013;Nishina et al., 2014). Our approach is therefore analysing how current models behave and why current models simulate a large range in SOC. Our purpose is not to provide improved results of the response of soil carbon to climate change but rather to better understand the implications of existing approaches, used in CMIP5, to parameterisation and initial value prescription described in Sect. 2.2.

Model set-up and experiments
We configure the reduced complexity model in a spatially explicit way to represent global variations, implemented as a surrogate for the CASA-CNP biogeochemical module (Wang et al., 2010) of the CABLE land surface model . A previous simulation by CABLE coupled to the coarse-resolution CSIRO Mk3L climate model (3.2 • latitude × 5.6 • longitude; Phipps et al., 2011) and driven by CMIP5 atmospheric CO 2 data provides monthly NPP, T s and θ s to the reduced complexity model. We use both historical simulations (Exbrayat et al., 2013b) and 21st century projections using the Representative Concentration Pathway 8.5 (RCP 8.5) atmospheric concentration scenario. We perform a sensitivity analysis by running the simple model with various combinations of a Q 10 value and a baseline residence time k. We use 11 equally spaced values of Q 10 ranging from 1.5 to 2.5 (i.e. intervals of 0.1), and 31 equally spaced values of k ranging from 120 months to 480 months (i.e. intervals of 12 months). These values are based on the range of results previously obtained by Todd-Brown et al. (2013) with their own reduced complexity model. Each value of Q 10 is applied with each value of k for a total of 341 simulations. Model versions are initialised via a classical spin-up procedure (Xia et al., 2012) using input data from 1850 to 1859 for 10 000 years to ensure all soil carbon pools reach a steady state. We then continue simulations with NPP, T s and θ s data from 1850 to 2005, and continue with RCP 8.5 projections to 2100. We note that these drivers do not include the representation of land use and land cover change and their effect on NPP, T s and θ s . Therefore, SOC input are likely to be higher than in reality. However, as stated earlier we are using the reduced complexity framework to understand the behaviour of the SOC model in response to variations in its parameters, and we do not aim to provide improved estimates of global scale terrestrial carbon sinks. In each model version, both k and the sensitivity of R h to temperature (represented by Q 10 ) are constant globally, in accordance with observations (Mahecha et al., 2010) and state-of-the-art models (Todd-Brown et al., 2013;. However, the actual value of the environmental scalar f T will of course vary spatially and temporally as a function of T s . Since we keep the same formulation of f W between model versions, we can attribute differences in results to the values of Q 10 or k.

Harmonized World Soil Database
The Harmonized World Soil Database (HWSD; FAO, 2012) combines several national inventories and provides a number of chemical and physical soil properties at a 30 arcsec resolution globally. However, despite the availability of this data set, CMIP5 models exhibit a 6-fold range in their total soil carbon content (Todd-Brown et al., 2013) including values well outside the uncertainty boundaries of observational data.
We previously showed that simply using the global amount of SOC from the HWSD data set to discriminate between acceptable and unacceptable simulations resulted in a nonnegligible reduction of the uncertainty in historical net carbon uptake (Exbrayat et al., 2013b). While we do not aim to provide CMIP5-like projections of the soil carbon balance with our reduced complexity model, we investigate the value of using the HWSD to discriminate between plausible and implausible simulations. We follow the method described by Todd-Brown et al. (2013) to derive an estimate of current total soil carbon from the latest version of the HWSD. First, we re-grid the original 30 arcsec raster to a 0.5 • × 0.5 • resolution. Within each half-degree cell we select the dominant soil type. For each soil type, the database provides bulk density and organic carbon content for a top layer (0-30 cm depth) and a bottom layer (30-100 cm depth). This allows us to calculate soil C density (in kg C m −2 ) in each cell. We then multiply each grid cell by its area and sum to obtain a global estimate of ∼ 1170 Pg C. Similarly to Todd-Brown et al. (2013) we also consider the uncertainty associated with our re-gridding process as well as analytical measurements of soil properties. We therefore obtain a 95 % confidence interval (CI 95 ) of 29 % below the mean to 32 % above the mean, or ∼ 830-1550 Pg C. We provide these gridded data as supplementary material. Due to the 6-fold range of SOC simulated by CMIP5 models (Todd-Brown et al., 2013), we believe that global SOC stocks from the HWSD can already represent a strong constraint to discriminate between different simulations.  Fig. 1b) and at the end of the projections with forcing corresponding to RCP 8.5 (in 2100, Fig. 1c). Figure 1a shows that the spin-up procedure causes different model versions to equilibrate at widely varying levels of total soil carbon despite the use of the same boundary conditions of NPP and T s . Differences in residence time k contribute most of the ∼ 592 to 2745 Pg C range, with larger values of k resulting in larger pools (Fig. 1a). Variations in the Q 10 parameter of f T have a smaller influence on total soil carbon but lower values do result in lower total soil carbon. For the same value of k, simulations with Q 10 = 1.5 equilibrate with total soil carbon equal to 86 % ±0.005 % (mean ±1 standard deviation) of the amount with Q 10 = 2.5. Figure 1b shows that the distribution of total soil carbon between model versions does not vary much during historical simulations . Models with large total soil carbon pools over this period remain  Fig. 1c continues to indicate a strong control of k on the total soil carbon in 2100. The projected range narrows to ∼ 684 to 2825 Pg C throughout the 21st century. However, we note there is an inversion in the influence of Q 10 on simulated total soil carbon with lower values of Q 10 resulting in larger pools especially for longer baseline residence times k. Nevertheless, this is still minor compared to the influence of k on C s .

Total soil carbon and global balance
Although the range in simulated soil carbon remains similar through time, non-negligible changes occur. This is highlighted in Fig. 2 which shows C s , the change in total soil carbon as a function of model parameters k and Q 10 for the historical simulations (1850-2005, Fig. 2a) and RCP 8.5 projections (2006-2100, Fig. 2b). First, Fig. 2a clearly shows that all model versions act as a net carbon sink during historical simulations, accumulating between 81 and 283 Pg C. Model versions with longer residence time k tend to accumulate more carbon through time. However, models with the largest value of Q 10 tend to accumulate only 69 % ± 0.4 % (mean ±1 standard deviation) of the amount that the lowest Q 10 models do. By analysing Fig. 2b, we see that the influence of Q 10 on the total soil carbon balance grows during RCP 8.5 projections where Q 10 now determines whether the soil remains a sink or becomes a source. This change between a source or a sink for different Q 10 values follows a near-linear relationship with k (solid line on Fig. 2b). Interestingly, the −179 to 168 Pg C range in the change in total soil carbon during RCP 8.5 is mostly a function of Q 10 as both extremes are achieved with the longest residence time used here. In other words, while Q 10 decides the sign of the change, k, and hence the initial stocks of SOC after spin-up, drives the magnitude of the response.
If we consider only models that fall within the CI 95 of the HWSD for current total soil carbon (dashed contours on Fig. 2a and b) the spread in simulated total soil carbon balance is largely reduced. During the historical simulations, the range of this subset of models shrinks by 84 Pg C to between 87 and 205 Pg C. It corresponds to a reduction of about 42 % of the initial uncertainty. Similarly, the range in projected soil carbon balance is reduced by 157 to −129 to 61 Pg C, a reduction of about 45 % of the initial uncertainty. We note, however, that this restriction does not necessarily increase confidence in sign of the future soil carbon change under RCP 8.5.
Differences in the behaviour between the full set of models and this subset of observationally constrained models can be seen in the time series and probability density functions (PDFs) for the historical period, shown in Fig. 3. First, the time series from 1850 shows there is no noticeable difference between the full set of simulations (in grey) and the subset of simulations with acceptable current soil carbon (in green) until 1900. During the first half of the 20th century, stronger sinks are excluded as they lie outside the CI 95 range, which correspond to the upper tail of the distribution of C s (see PDF inset for 1950). However, the kurtosis of the distribution, or most probable change from our simulations, changes negligibly. After ∼ 1960, we observe a step change in cumulative C s that follows a strong response in NPP to the rapid increase in atmospheric CO 2 (please refer to Exbrayat et al., 2013b for a more detailed account of this behaviour). The spread between simulations grows and most of the excluded simulations based on the CI 95 range are the strongest sinks (as in Fig. 2a) while a few of the least accumulating simulations are also excluded. This does have a large impact on the most probable change in storage, reducing it from ∼ 200 to ∼ 140 Pg C.
We now examine future simulations and present time series and PDFs of change in total soil carbon during RCP 8.5 projections in Fig. 4. All simulations continue to accumulate carbon at the beginning of the 21st century and remain net carbon sinks until about 2060. At the end of the century, some model versions have simulated positive C s corresponding to a net carbon sink over the 21st century, while other ends their projections with negative C s , or a net carbon loss. However, all simulations show the same overall behaviour with first an increase in C s that peaks, and then a decrease in C s . The timing of the peak, i.e. when soil carbon starts to deplete, varies between ∼ 2035 and 2075 and is explained by the value of Q 10 (R 2 = 0.74, data not shown) with higher values leading to an earlier peak. This indicates that, in all simulations, soil has become a net source of carbon by the end of the 21st century, regardless how much carbon was accumulated since 2005, and hence since 1850. The PDFs in 2050 show that selecting only observationally consistent models results in the most heavily accumulating simulations, i.e. those that would peak later, to be dismissed. However, by 2100, both the lower and upper tails of the initial distribution are clipped, reducing the simulated range from −178 to 168 Pg C (all simulations) to −129 to 61 Pg C. In both cases, differences in the kurtosis of both distributions remain very small which indicates that our selection scheme dismisses outliers. We note that the lower bound of C s for both sets of models is the same until late in the projections (∼ 2085).

Regional differences
Although Fig. 1 indicates that the range in k can explain most of the variability in total soil carbon content at equilibrium and hence through transient simulations, Q 10 is likely to influence the local response of f T . Figure 5 shows the relative value of f T for different temperatures and values of Q 10 . Since the chosen T ref = 15 • C, all Q 10 values lead f T to be equal at this particular temperature. However, the more difference there is between the actual temperature and T ref , the more sensitive f T becomes to values of Q 10 . As our simulations are spatially explicit, this may introduce non-negligible regional differences in C pools at equilibrium and their response to transient changes in T s and NPP.
To investigate this more in detail, we present the zonal averages of soil C density for different values of Q 10 with k set to 180 months (Fig. 6). We choose this particular residence time as example because all corresponding simulations are within the CI 95 of the HWSD for 2005 regardless the value of Q 10 . Figure 6a shows that Q 10 values do introduce non-negligible differences in local equilibrated soil C density. Steady-state pools at low latitudes (30 • S to 30 • N) are larger with low values of Q 10 (blue in Fig. 6). Conversely, high latitude pools are larger with high values of Q 10 (red in Fig. 6). Overall, the range in the value of zonally averaged soil C density at equilibrium is up to 3-fold depending on the chosen value of Q 10 . This is particularly obvious in regions with high NPP including low-latitude tropical rainforests or northern taigas. As was the case with total C s , the zonal distribution soil C density and the relative position of simulations with different Q 10 do not vary much between 1850 and 2005 (Fig. 6b) although there is a slight shift towards uniformly higher densities as all model versions are net global carbon sinks (Figs. 2a and 3). The pattern of zonal soil carbon remains essentially the same at the end of RCP 8.5 projections. However, models with lower values of Q 10 now have more carbon than those with high values of Q 10 over a broader zone (40 • S-50 • N). Figure 7 shows the zonal change in soil C density for the same simulations as in Fig. 6. Figure 7a indicates that all simulations simulate a net sink almost everywhere during historical simulations, except at latitudes > 70 • N. However, the strength of this sink is strongly dependent upon the value of Q 10 , especially in low latitudes. There is an approximately 2-fold difference between the high accumulation of low Q 10 models, and the low accumulation of high Q 10 models. Differences between Q 10 values are negligible at higher latitudes. Figure 7b shows the same information for RCP 8.5 projections. Simulations with lower values of Q 10 almost always accumulate more C (except between 0 • and 10 • N). While all model versions with k = 180 months lose carbon at low latitudes (20 • S-20 • N), and gain carbon at high latitudes in the Northern Hemisphere (> 50 • N), the value of Q 10 , and hence the environmental scalar f T , decides the sign of the local soil C balance in the 21st century at mid-latitudes. Within the mid-latitudes, high values of Q 10 are more likely to simulate a net loss of soil carbon. We can therefore narrow down the dependence of the global C s on Q 10 to its affect at mid-latitudes.

Effect of k and Q 10 on soil carbon
In our simulations, the range in total soil carbon at equilibrium (∼ 592 to 2745 Pg C) depends on which value of Q 10 and especially k is used (Fig. 1a). This range captures the ∼ 561 to 2938 Pg C range in soil carbon in CMIP5 in 1860 (see Fig. S1). We note of course that CMIP5 models not only vary in their soil C component, but simulate different NPP and T s and also integrate a range of soil moisture limitations (Todd-Brown et al., 2013). The range achieved here at the end of the historical simulations (∼ 709 to 2943 Pg C) is, for example, larger than the 1090 to 2646 Pg C range in 2000 from seven DGVMs in the ISI-MIP project  which were driven by a harmonised weather data set. We can attribute this range to the first-order representation of decomposition and its response to the initialisation procedure used in most CMIP5 simulations. By spinning-up the model, the goal is to stabilise pools so that total NPP is exactly compensated by total R h over the selected period of time (here 10 years). In Eq. (2), a longer residence time k results in a lower decay rate (i.e. R h per unit of C s ). Therefore, model versions that have a slower turnover will require more substrate to simulate the same R h needed to compensate NPP. As the baseline residence time k is applied globally, it drives the global pool size (Fig. 1) much more than changing Q 10 affects f T . However, as seen in Fig. 6, when considered regionally, Q 10 plays a non-negligible role for the local response of decomposition and the definition of equilibrium soil C density. High values of Q 10 lead f T to trigger strong decay rates in warm regions (Fig. 5) that require less substrate (see low latitudes in Fig. 6a) to compensate the same NPP. Conversely, high Q 10 lead to low values of f T in cold regions. Therefore, more substrate is required to bring the pool to equilibrium as seen in high latitudes in Fig. 6a. Low values of Q 10 show an opposite regional behaviour. Regional differences compensate each other and therefore f T with different Q 10 values can only explain a small fraction of the range in equilibrated total soil carbon. Of course, if another T ref was used, the relative differences between f T with different Q 10 would be altered and the influence of Q 10 and its effect on f T on total and local C s would vary. Furthermore, the difference between f T with different Q 10 grows with the absolute value of the difference T s − T ref . Therefore, using a value of T ref that is outside the range of actual temperatures would lead f T with different Q 10 to keep the same relative position globally. It would introduce larger relative differences between these functions.
Comparing Fig. 1a, b and c suggest that the range in total C s at equilibrium is a good predictor of the current and future 7006 J.-F. Exbrayat et al.: Disentangling residence time and temperature sensitivity range in total soil carbon. Despite differences in the magnitude of the change in C s through time (Friedlingstein et al., 2014), equilibrium conditions achieved under pre-industrial conditions largely define current and future pool sizes as observed in CMIP5 models . Examining Fig. 6 confirms that this global effect can also be seen regionally, especially in low (20 • S-20 • N) and high (> 50 • N) latitudes, where carbon pools are largest. This is of concern as substrate availability also influences R h and hence its response to changes.
Changes in C s through time are nevertheless nonnegligible, and it is important to quantify the response of the system to perturbations. Our results show increasing atmospheric CO 2 concentrations enhances NPP more than the simultaneous warming enhances R h during historical simulations. This historical net carbon sink that is driven by the response of vegetation to increasing atmospheric CO 2 (and hence SOC in ) is in accordance with previous studies (Friedlingstein et al., 2006;Sarmiento et al., 2010;Zhang et al., 2011;Wania et al., 2012;Anav et al., 2013;Exbrayat et al., 2013b). Therefore, all model versions with longer residence time accumulate more C s over the same time period as a result of a slower turnover of carbon in soils, and this mirrors the state of the equilibrium stores. However, despite the dominance of the increased NPP on C s , the historical warming signal is influential. Specifically, those model versions more sensitive to changes in temperature (i.e. with high values of Q 10 ) accumulate less soil carbon during the 20th century even though they initially equilibrated with larger global pools. This is also true of local soil C density where high Q 10 values are less accumulating regardless of the initial soil C density. We however note that the value T ref used in our experiments is well within the range of actual temperatures. Therefore, the historical warming does not induce large changes in the values of f T with different Q 10 .
Projections under the strong-forcing RCP 8.5 scenario also see an increase in the influence of the value of Q 10 on C s . Figure 2b clearly shows that the capacity of soils to become carbon sources or remain sinks depends almost entirely on the Q 10 parameter, and both states can be achieved for any value of k used while remaining within range of previous studies (Friedlingstein et al., 2014;Nishina et al., 2014). Figure 7b indicates that this is clearly a result of differences in the local response of model versions in the mid-latitudes as a function of Q 10 . Such regional discrepancies leading to a change in the sign of global C s models have also been highlighted through a recent intercomparison project that used a harmonised weather data set to drive seven biome models . However, contrary to this previous study, none of our model versions accumulates soil carbon in the inter-tropical region during the 21st century. This is probably due to the fact that we use the same boundary conditions of NPP and T s for all our model versions, while models used by Nishina et al. (2014) used a prescribed weather data set but were left free to simulate their own NPP.
Overall, the globally applied model parameter k drives the steady-state response of our reduced complexity system. However, the more conditions change (i.e. steady state to historical to RCP 8.5 projections), the more the dynamic transition of the system towards a new equilibrium depends on the environmental scalar f T and the specific value of Q 10 . Although the same formulation of f T is applied globally, differences in its response to local T s sum up to determine the sign of total soil carbon balance. We also note that model versions that equilibrate as a result of longer baseline residence time k have a tendency to produce a larger absolute response of total soil carbon balance. Therefore, the size of pools to which the change is applied seems to dominate the response even when higher values of k imply a smaller relative change in the decay rate k −1 × f T × f W used in Eq. 2. This control of initial conditions obtained by spin-up on the response of the system is a critical aspect that needs to be better resolved, especially since recent intercomparison experiments all exhibit huge discrepancies in equilibrium conditions of participating models (Anav et al., 2013;Todd-Brown et al., 2013;Nishina et al., 2014).

Discriminating between model versions
Since k clearly influences the total soil carbon content at equilibrium in 1850, it is a good predictor of the current total soil carbon content. Therefore, k is the key parameter that decides how much carbon is active in the modelled system, and whether model versions fall within the CI 95 of the HWSD. Here, all simulations with baseline residence time between 150 and 250 months fulfil this requirement regardless of which Q 10 is used in f T .
If we isolate these simulations, the range in total soil carbon change shrinks by 42 and 45 % for the historical simulations and RCP 8.5 projections, respectively. However, while this selection dismisses outliers it does not increase confidence in the sign of the soil carbon change. This is because regional differences lead to similar values in total soil carbon for different values of Q 10 . These regional differences translate into heterogeneous responses under RCP 8.5 forcing, especially in mid-latitudes. They are sufficient to induce a change of sign in the global soil carbon balance.

Conclusions
We have used a reduced complexity model, broadly representative of current state-of-the-art models of soil organic C decomposition used in CMIP5 and ISI-MIP experiments, to explore the response of microbial decomposition to climate change on soil C dynamics at regional and global scale. We have shown that key parameters in the first-order representation of decomposition interact in markedly different ways depending on the nature of forcing and antecedent conditions. First, the time and space-invariant baseline residence time decides of the total soil carbon content at equilibrium after spin-up, typically the process used by CMIP5 models to initialise C pools. Next, the more boundary conditions imposed on the system move away from the equilibrium forcing, the more the environmental scalar describing the sensitivity of the system gains in importance. However, it is the size of the pool to which the change is applied that mostly controls the magnitude of the response.
Applying a constraint on total soil carbon that discriminates between acceptable simulations of total soil carbon leads to a drastic reduction of the range of simulated change. Meanwhile, most of the remaining uncertainty in 21st century projections of total soil carbon can be attributed to zonal differences in the response to change, especially at midlatitudes. These do not allow us to confidently project soil as either a global source or sink of carbon for the 21st century. However, it is clear that under RCP 8.5 tropical soils are not suited for long-term carbon storage, while some more potential exists in high latitudes.
Finally, we suggest that future estimates of terrestrial, and especially soil, carbon responses to climate change should be more constrained by available data sets of carbon stocks. This is critical as model structures describe fluxes as a fraction of the substrate pool size. So far, the process of spin-up has too many degrees of freedom that lead to model-specific amounts of active soil carbon.
The Supplement related to this article is available online at doi:10.5194/bg-11-6999-2014-supplement.