Oscillatory Behavior of Two Nonlinear Microbial Models of Soil Carbon Decomposition

A number of nonlinear models have recently been proposed for simulating soil carbon decomposition. Their predictions of soil carbon responses to fresh litter input and warming differ significantly from conventional linear models. Using both stability analysis and numerical simulations, we showed that two of those nonlinear models (a two-pool model and a three-pool model) exhibit damped oscillatory responses to small perturbations. Stability analysis showed the frequency of oscillation is proportional to ε −1 − 1 K s /V s in the two-pool model, and to ε −1 − 1 K l /V l in the three-pool model, where ε is microbial growth efficiency, K s and K l are the half saturation constants of soil and litter carbon, respectively, and V s and V l are the maximal rates of carbon decomposition per unit of microbial biomass for soil and litter carbon, respectively. For both models, the oscillation has a period of between 5 and 15 years depending on other parameter values, and has smaller amplitude at soil temperatures between 0 and 15 • C. In addition, the equilibrium pool sizes of litter or soil carbon are insensitive to carbon inputs in the nonlinear model, but are proportional to carbon input in the conventional linear model. Under warming, the microbial biomass and litter carbon pools simulated by the nonlinear models can increase or decrease, depending whether ε varies with temperature. In contrast, the conventional linear models always simulate a decrease in both microbial and litter carbon pools with warming. Based on the evidence available, we concluded that the oscillatory behavior and insensitivity of soil carbon to carbon input are notable features in these nonlinear models that are somewhat unrealistic. We recommend that a better model for capturing the soil carbon dynamics over decadal to centennial timescales would combine the sensitivity of the conventional models to carbon influx with the flexible response to warming of the nonlinear model.


Introduction
A number of soil and litter carbon decomposition models based on Michaelis-Menton kinetics, or nonlinear soil carbon models, have recently been developed (Schimel and Weintraub, 2003;Allison et al., 2010;German et al., 2012).These models can simulate the priming of existing soil carbon stocks with additional litter inputs (Kuzyakov et al., 2000;Wutzler and Reichstein, 2013), can replicate the acclimatory response of soil carbon decomposition to warming Published by Copernicus Publications on behalf of the European Geosciences Union.(Allison et al., 2010;Bradford et al., 2008), and can represent the spatial variation of global soil carbon better than models based on conventional linear soil carbon dynamics (Wieder et al., 2013).The key to this success has been to account explicitly for the effects of both microbial biomass and substrate concentration on the rate of soil carbon decomposition in the nonlinear models, compared to the linear models for which it is assumed that only the amount of substrate is limiting to the rate of soil carbon decomposition.Consequently, the conventional linear models cannot simulate the effect of priming on soil carbon decomposition well without significant modifications (see Wutzler and Reichstein, 2008).
However, the conventional linear models have been successfully used to capture the soil carbon dynamics over interannual to decadal timescales (see Parton et al., 1993).A global synthesis of litter decomposition (Zhang et al., 2008) and a decade-long litter decomposition study in diverse climates (Adair et al., 2008;Bonan et al., 2013) both indicate the monotonic decay of a quantity of litter over time in diverse ecosystems.Long-term soil incubation data show the similar responses (Li et al., 2013).These dynamics can be predicted well by conventional linear soil carbon models provided they use multiple carbon pools (Bolker et al., 1998).However, the interannual to decadal response of soil carbon to inputs has not been analyzed in detail for the nonlinear models of soil decomposition to verify whether they can also capture these dynamics.Such analyses are needed because nonlinear models can potentially exhibit a much richer range of behaviors than linear models, not all of which may be desired or intended.For example, it is well known that a system of nonlinear ordinary differential equations, such as a nonlinear soil model, can become unstable in response to a small perturbation to its initial pool sizes (Raupach, 2007) or inputs and can switch between different equilibrium states in response to climatic variation (Manzoni and Porporato, 2007), although there is presently no evidence that soil carbon dynamics exhibits such characteristics over interannual to decadal timescales.In this paper, we analyze the stability of two recently published nonlinear soil carbon models in response to nonperiodic change in carbon input: the two-pool model developed by German et al. (2012) and the three-pool model simplified from Wieder et al. (2013), and pay particular attention to the time course of the responses to perturbation over decadal timescales.We focus on the intrinsic oscillatory responses of the modeled system to perturbations and, therefore, do not consider the forced responses of the system to oscillations in external factors, such as through diurnal variation in soil temperature or seasonal variation in carbon input.We address the following questions.(1) What are the responses of these two models' stable to small perturbations in initial pool sizes?(2) What determines the stability of these two models?(3) How different are the simulated responses of soil carbon to climate-induced changes in litter inputs and warming between the nonlinear and conventional linear models?We conclude by discussing which of the pre- dictions from two types of models are more consistent with empirical evidence from field and laboratory studies.

Two nonlinear soil carbon models
The two nonlinear soil carbon models developed by German et al. (2012) and Wieder et al. (2013) were based on the earlier work by Schimel and Weintraub (2003) and Allison et al. (2010).Both models assume that the decomposition rate of litter or soil carbon is proportional to the biomass of decomposers (soil microbes), and varies with substrate concentration (litter or soil carbon) following Michaelis-Menten kinetics.Growth rate of soil microbes is proportional to the rate of carbon decomposition, and mortality of soil microbes is proportional to their biomass (see Fig. 1).
The equations for the two-pool soil model developed by German et al. (2012) are and where C b and C s are the pool sizes of soil microbial biomass and soil organic matter (in g C m −2 ), ε is microbial growth efficiency or fraction of assimilated carbon that is converted into microbial biomass (unit-less), µ b is the turnover rate of microbial biomass per year, F npp is carbon influx into soil (in g C m −2 year −1 ), and V s and K s represent the maximum rate of soil carbon assimilation per unit microbial biomass per year and the half-saturation constant for soil carbon assimilation by microbial biomass (in g C m −2 ), respectively.
We also analyze a simplified version of a three-pool nonlinear soil carbon model developed by Wieder et al. (2013).Compared with the two-pool model, the three-pool model includes an additional litter carbon pool (C l ).Only a fraction of the carbon influx goes to the soil organic matter pool and the rest goes to the litter pool.The dynamics of these three carbon pools are and where α is the fraction of carbon influx that directly enters the soil organic matter pool, and usually is less than 0.05, V l and K l represent the maximum rate of litter carbon assimilation per unit microbial biomass per year and the halfsaturation constant for litter carbon assimilation by microbial biomass (in g C m −2 ).To simplify mathematical analysis, we assume that α = 0, and will discuss the case with α > 0 later.The equilibrium microbial pool sizes are identical for both models: where * denotes equilibrium values.However the equilibrium soil carbon pool sizes are different between the two models.For the two-pool model it is whereas it is less in the three-pool model: because part of the nonmicrobial equilibrium biomass is held in litter, with when α > 0, the equilibrium pool sizes of litter or soil carbon pools depend on α, and are given by

Parameter values
To ensure positive pool sizes, the following constraints are applicable: > 1, for the two-pool soil model and > 1, for the three-pool soil model.The parameter values given by German et al. (2012) are for a soil column with finite volume.Here we use both models to represent the dynamics of the top 1 m of soil.Based on the results of German et al. (2012), the temperature dependence of the model parameters evaluated here is given by where T s is mean temperature of the top 1 m soil (in • C) and parameters x vl , x vs , x kl and x ks are tunable parameters that we use to scale the rate of litter and soil decomposition, with x vl , x vs ∈ [1, 10], and x kl , x ks ∈ [0.1, 1] .V in Eq. ( 12) is multiplied by a factor 24 × 365 to convert from per hour to per year, and K is multiplied by 1000 to convert from milligrams of carbon per cubic centimeter to grams of carbon per square meter (i.e., mg C cm −3 to g C m −2 ) for the top 1m soil depth.
The value used in this study for µ b is 4.38 year −1 , analogous to rates in German et al. (2012) and Wieder et al. (2013).We impose a maximum value of 0.6 for ε, based on the work of Sinsabaugh et al. (2013) and a lower limit of 0.001.Both numerical simulations and analytical techniques were used to study the stability of the two models.For all numerical simulations, unless specified otherwise, we used a constant carbon influx of 345 g C m −2 year −1 representing the mean net primary production of the global land biosphere under present climate conditions (Field et al., 1998), varied soil temperature from −10 to 35 • C, and modified x vl , x vs , x kl and x ks within their respective ranges.The default values are T s =15 • C, x vl = x kl =8, and x kl = x ks =0.2.At steady state, the simulated carbon pool sizes using the default values are 50.Changes in microbial biomass (in g C m −2 ; upper panel) or soil organic carbon (in g C m −2 ; middle panel) (both in g C m −2 , 0-100 cm soil depth) in time or soil organic carbon against microbial biomass (lower panel) following a 10 % reduction in initial pool sizes at time t = 0 from their respective equilibrium values for the two-pool soil model.Soil temperature was constant at 15 • C for this simulation.

Numerical simulations
To study the stability of the two nonlinear models numerically, we initialized the models with 90 % of the microbial biomass and soil organic carbon equilibrium pool sizes, and 100 % of the litter carbon equilibrium pool, and ran both models for 500 years.
Figure 2 shows that the responses of both microbial biomass and soil organic carbon pools in the two-pool model to a small perturbation are oscillatory, and the amplitude of oscillation decreases with time.During the first 50 years, the microbial biomass pool varies between 12 and 114 g C m −2 while soil organic carbon varies between 11508 and 13412 g C m −2 .Following perturbation, both microbial biomass and soil organic carbon spiral towards their equilibrium states with the width of spiral decreasing exponentially with time (see Fig. 2).
For the three-pool model, the response of the soil organic carbon is monotonic while the responses of litter carbon and microbial biomass are oscillatory (see Fig. 3).At a soil tem- Changes in litter carbon (upper panels), microbial biomass (middle panels) or soil organic carbon (lower panels) following a 10 % reduction in the initial pool sizes in microbial biomass and soil organic carbon at time t = 0 from their respective equilibrium pool sizes for the three-pool soil carbon model.The three panels on the left are for soil temperatures of 15 • C and the right panels for soil temperatures of 35 • C. The unit (g C m −2 ) is the same for all pool sizes (0-100 cm soil depth).perature (T s ) of 15 • C, the amplitude of the oscillation decreases to less than 5 % of its initial value after 20 years.Thus, both models have oscillatory responses to a perturbation to their initial values; however, the oscillatory response in the three-pool model has a smaller amplitude than the two-pool model at the same T s .At higher temperatures (T s = 35 • C) the oscillatory responses of litter carbon and microbial biomass are much stronger than at T s = 15 • C for the three-pool model (see Fig. 3).
The oscillatory responses of both models are a result of the interaction between substrate availability (litter or soil carbon) and decomposer biomass (soil microbial biomass).In both models, the rate of carbon decomposition is proportional to both the biomass of decomposers and the carbon concentration of the substrate (Eqs.1-5).For the three-pool model, when soil microbial biomass is reduced below its equilibrium value at t =0, the rate of litter carbon decomposition and the growth rate of soil microbial biomass are lower than their respective values at equilibrium.As a consequence, litter carbon increases and microbial biomass decreases.This leads to more litter carbon being available for soil microbes, and the growth rate of soil microbial biomass and rate of litter decomposition rate consequently increase, soil microbial biomass increases and litter carbon decreases.These changes in substrate and microbial biomass carbon pools will result in oscillatory behavior (see Fig. 3).The amplitude of oscillation decreases exponentially with time until both oscillatory pools reach their respective equilibrium states.

Mathematical analysis
Numerical simulations show that the responses of both models to a small perturbation are oscillatory and stable.We therefore conducted linear stability analyses to evaluate the stability of the two models more generally.This technique has been used in many studies of ecological models, biogeochemical models and human-carbon cycle interactions (see Manzoni et al., 2004;Raupach, 2007).
The linear stability of a system dy/dt = f (y, t), where y is a state vector, and f is a function -to small perturbations can be determined by the Jacobian matrix (J = df /dy) of the system.We linearized the nonlinear systems by assuming that changes of carbon pool sizes near their equilibrium values are proportional to the difference between the size of each pool and its equilibrium value, i.e., dy/dt = J (y − y * ), where y * is the equilibrium value of y (see Appendix A for further details).If the response is stable, then the system will return to its equilibrium state, otherwise, the system states (or carbon pools in our case) will depart from their initial values indefinitely.The stability of a linearized system is fully characterized by its eigenvalues and corresponding eigenvectors of matrix J .A stable system always has nonpositive eigenvalues or nonpositive real parts, in case of eigenvalues being complex numbers, and an equilibrium is locally stable if all eigenvalues of the Jacobian are negative or have negative real parts, respectively.The response to a small perturbation is oscillatory if the eigenvalues are complex, or monotonic if the eigenvalues are real (Drazin, 1992).
If an eigenvalue related to the response of a carbon pool, λ j , is complex (λ j = a j + b j i, where i 2 = −1), the change of that carbon pool size with time close to equilibrium is oscillatory with a frequency of oscillation of b j or period (p) of 2π /b j (Roe, 2009).The amplitude of oscillation changes with time at a rate of exp(a j t).The half-life time (t 0.5 ) of the amplitude decrease, or the time taken for the amplitude to be reduced to half the initial value, can be calculated as −ln(2)/a j if a j < 0 (Drazin, 1992).
As shown in Appendix A, the eigenvalues (λ 1 , λ 2 ) of the linearization of the two-pool model are where whereas the corresponding eigenvalues for the three-pool model (λ 1 , λ 2 and λ 3 ) are where > 0 and Both models have two complex eigenvalues with negative real parts, while the third eigenvalue of the three-pool model is a negative real number.This indicates that small perturbations to the equilibrium carbon return to those equilibrium states through a series of damped oscillations, confirming what was observed in the numerical simulations.
We used Eqs.( 18) and ( 20) to quantify the properties of the oscillatory responses using p (period) and t 0.5 (half-life).For the two-pool model these are and Similarly, for the three-pool model: and while the third pool (soil organic carbon) decays monotonically towards equilibrium according to exp(λ 3 t) and λ 3 ∝ − V s (ε −1 −1)K s .Thus, the properties (period and half-life) of the oscillatory approach to equilibrium in the microbial pools in both models, and the soil pool in the two-pool model or litter pool in the three-pool model, depend on the ratio of the half saturation constant to the maximum decomposition rate.In addition, the response of organic soil carbon to a small perturbation in the three-pool model is totally decoupled from the responses of litter and microbial carbon (see Appendix A). x ks at different values of x vs (1, 3, 5, 7, 9), (b) 1/t 0.5 with x vs at different x ks (0.1, 0.3, 0.5, 0.7, 0.9); (c) period of oscillation (p) with x ks at different values of x vs (1, 3, 5, 7, 9) and (d) frequency of oscillation (1/p) with x vs at different x ks (1,3,5,7,9).(1,3,5,7,9), (b) 1/t 0.5 with x vs at different x ks (0.1, 0.3, 0.5, 0.7, 0.9); (c) period of oscillation (p) with x ks at different values of x vs (1, 3, 5, 7, 9) and (d) frequency of oscillation (1/p) with x vs at different x ks (1, 3, 5, 7, 9).
To illustrate the numerical consequences of these analytic results, we plot the values of t 0.5 and p, or their reciprocals, t −1 0.5 and 1/p (frequency) with different tuning parameters x ks or x vs , respectively (see Eqs. 13, 14 and 16) in Fig. 4. As expected, half-life (t 0.5 ) increases linearly with x ks for a given x vs , or 1/t 0.5 decreases linearly with x vs for a given x ks .The period increases with x ks or the frequency of oscillation (1/p) increases with x vs , and both increases are nonlinear.
The oscillatory response with longer period and half-life is likely to reach the equilibrium steady state much slower.Both the period and half-life increase with x ks /x vs , therefore the system with larger value of x ks /x vs will reach steady state much slower after some perturbation.Since the oscillatory responses of both models are related to the complex eigenvalues, and the dependence of two complex eigenvalues on x ks /x vs in the two-pool model is quantitatively similar to the dependence on x kl /x vl in the three-pool model, it might seem that the oscillatory response of two models should be quite similar for a given carbon influx (F npp ), with x vl = x vs , and x kl = x ks .However, the simulated oscillatory responses of two-pool model are much stronger than those of the threepool model.
This can be explained by studying the solution to the linearized system.As shown in Appendix B, the maximal value of the microbial biomass pool during the first cycle of oscillation is proportional to (C s0 -C b0 a 11 )/b 11 for the two-pool model (C b max 2 ), and to C l0 −C b0 a 11 +C s0 a 32 a 11 −C s0 a 31 b 11 for the three-pool model (C b max 3 ), where a 11 and b 11 are the values of the real and imaginary parts of the first element of the eigenvector for the microbial biomass carbon and are the same for both models; C s0 , C b0 , and C l0 are the initial pool sizes of soil organic carbon, microbial biomass and litter carbon; and a 31 and a 32 are the first two elements of the third eigenvector corresponding to λ 3 .Given F npp , ε and µ b , the equilibrium pool sizes of microbial biomass in the two models are equal, but the soil carbon pool is at least one order of magnitude greater than litter carbon, the maximal amplitude of the oscillation cycle of microbial biomass in the two-pool model is therefore much greater than that of the three-pool model (see Appendix B).The half-life and period are similar because the value of K s /V s in the two-pool model is equal to K l /V l in the three-pool model if x ks /x vs = x kl /x vl .

The temperature dependence of system stability and soil carbon decomposition
The responses and stability of the linearized systems vary with soil temperature (see Figs. 3 and 4).In this study, the values of V s , K s , V l or K l at a reference temperature (T s =15 • C) are varied by a factor of x vs , x ks , x vl or x kl (see Eqs. 13, 14 and 16), respectively.To quantify the temperature responses of both models, we calculated t 0.5 and p for both models at different T s and different combinations of x ks and x vs for the two-pool model or x kl and x vl for the three-pool model using the analytic solutions (Eqs.22-25).
Figure 5 shows that t 0.5 and p increase with x ks /x vs for the two-pool model.Results for x kl /x vl from the three-pool model are quantitatively similar because the temperature dependence is quantitatively the same for V l and V s or K l and K s , and the differences in t 0.5 and p between the two models result from the terms in the square brackets in Eq. ( 19) for A and Eq. ( 21) for B1.Both those terms have a value between 0.3 and 0.98 within the soil temperature range from −10 to 35 • C. Figure 5 shows that the period of oscillation (p) increases with x ks /x vs for the two-pool model.For a given value of x ks /x vs , t 0.5 decreases with T s when T s < 0 • C, or when T s > 20 • C; and p decreases with an increase in T s when T s < 10 • C, and increases with T s when T s > 30 • C. The temperature responses of the t 0.5 and p for the threepool model are very similar to those for the two-pool model, and are therefore not shown here.The results here suggest that both model systems will recover much faster after some perturbation at a soil temperature between 0-20 • C than at other temperatures, depending on the value of x ks /x vs for the two-pool model or x kl /x vl for the three-pool model.
These responses can be understood by analyzing the temperature dependence of different terms of t 0.5 in Eq. ( 22) for the two-pool model and Eq. ( 24) for the three-pool model.For a given F npp , there are three terms dependent on T s in t 0.5 for both models: the first term is ε −1 − 1, common to t 0.5 for both models; the second term is K s /V s in the twopool model or K l /V l in the three-pool model; and the third oscillation with x ks /x vs and T s (lower panel) in the two-pool model.The unit is year for both t 0.5 and 794 period.795 Fig. 5. Variation of t 0.5 with x ks /x vs and T s (upper panel) and period of oscillation with x ks /x vs and T s (lower panel) in the two-pool model.The unit is year for both t 0.5 and period.
term is that in square brackets with an exponent of −2 in the Eq. ( 22) for the two-pool model or Eq. ( 24) for the threepool model.As shown in Fig. 6 for the three-pool model, the first term (ε −1 -1) increases with T s , the second term K l /V l decreases exponentially with T s , and the third term decreases with T s but only varies within a small range.As a result, variation of t 0.5 with soil temperature is dominated by the responses of the first and second terms in Eq. ( 22) for the two-pool model or in Eq. ( 24) for the three-pool model, or t 0.5 ∝ (ε −1 − 1)K s /V s in the two-pool model or ∝ (ε −1 −1)K l /V l in the three-pool model.For the parameter values used in our study, t 0.5 is at a minimum when T s = 2 • C, and increases in warmer and colder conditions.However, if ε does not vary with T s , t 0.5 and the rate of the oscillatory pools approaching their steady state after some disturbance (exp(a 1 t)) increase exponentially with T s , with an equivalent Q 10 of 1.8 (note a 1 = − ln(2)t 0.5 ).
The response of p to soil temperature is quite similar to t 0.5 for both models, and can be calculated as the product of three terms (see Fig. 6).Because the third term only varies between 0.38 and 0.58 as soil temperature changes from −10 to 35 • C, variation of p with soil temperature is proportional to (ε −1 − 1)K s /V s in the two-pool model, or (ε −1 − 1)K l /V l in the three-pool model.The temperature dependence of microbial growth efficiency also affects the turnover rate of soil organic carbon.If ε varies with T s , the turnover rate of soil organic carbon in the three-pool model (−λ 3 ) is maximum around T s = 2 • C (see Fig. 7).However if ε does not vary with T s , the turnover rate of soil organic carbon increases exponentially with T s with an equivalent Q 10 of 1.8.warming if microbial growth efficiency does not vary with soil temperature; or microbial biomass decreases and litter carbon increases with soil warming if microbial growth efficiency decreases with an increase in soil temperature for the three-pool model (see Fig. 8).However soil carbon will decrease with soil warming, independent of the temperature dependence of microbial growth efficiency with an equivalent Q 10 of 1.8.This is explained later.The equilibrium pool sizes of litter or soil carbon do not change with carbon input (F npp ), but soil microbial biomass pool size at equilibrium will increase proportionally with an increase in F npp .This is easily seen by examining Eqs.(8-10) for the three-pool model with α = 0.For a conventional linear model of soil carbon decomposition, its pool sizes (litter, microbial biomass or soil carbon) will increase proportionally with an increase F npp , and decrease with soil warming (Xia et al., 2013).Therefore the responses to soil warming can be quite different for the litter and soil microbial biomass in the nonlinear models, and the response to an increase F npp can also be quite different for litter or soil carbon pools between the nonlinear and conventional linear models.
The above analysis is done by assuming α =0 although the results will be similar if α is small (< 0.05 for most field sites).To help explain the differences in the responses of the equilibrium carbon pools to change in F npp between the linear and nonlinear models, we developed a three-pool linear model that has same equilibrium pool sizes as the nonlinear model (see Appendix C).It should be noted that the transient responses of the three-pool linear model are different from the three-pool nonlinear model.We find that both the turnover rates of the litter and soil carbon pools at equilibrium are proportional to F npp and therefore their equilibrium pool sizes are independent of F npp .When α = 0, both the turnover rate of soil carbon and input of carbon to soil car- − 1), therefore the response of soil carbon to warming is independent of the temperature dependence of microbial growth efficiency (see Fig. 8).When α = 0 and is small, the response of soil carbon to warming will also depend on the temperature response of microbial growth efficiency, but the dependence is rather weak.
All of the analyses above relate to the properties of the carbon pools at equilibrium but do not tell us about how warming or changes in F npp affect the transient responses.
As shown in Fig. 9, the simulated responses of litter carbon or microbial biomass to a 5 • C warming by the nonlinear model are oscillatory, and converge to their equilibrium pool sizes along quite different trajectories from those by the conventional linear model.Furthermore the simulated litter carbon in response to warming by the nonlinear model decreases initially below the initial value, and then increases to above the initial value after 200 years if the microbial growth efficiency is reduced from 0.39 to 0.31 after 5 • C warming.This response cannot be reproduced by the linear model, as the response of the linear model to warming will always approach the new equilibrium exponentially, that is , where C l0 is the initial pool size, µ l is the turnover rate of litter carbon and C * l is the new equilibrium pool size.Therefore if the final equilibrium pool size is greater than the initial pool size, the pool size after time t = 0 is always greater than the initial pool size for the conventional linear models.
The simulated response of soil carbon to a 5 • C warming by the nonlinear model decreases with time much faster than that by the linear model for the first 50 years, then slower than the linear model after 100 years (see Fig. 9c and f).The simulated response of soil carbon by the nonlinear model cannot be accurately reproduced using one soil carbon pool, this is because the decomposition of soil carbon in the nonlinear model depends on microbial biomass linearly and soil carbon nonlinearly.However the simulated response of the nonlinear model can be accurately approximated using two soil carbon pools with first-order kinetics, as shown by the fitted pink curves in Fig. 9.

Discussion
In this study we show both analytically and numerically that the two nonlinear models of soil carbon decomposition have an oscillatory response to a small perturbation (Figs. 2, 3), and provide quantitative understanding of how model parameter values affect those oscillations (Figs. 4,5).
Our analysis shows that responses of microbial biomass and soil carbon in the two-pool model, or litter carbon and soil microbial biomass in the three-pool models, converge to new states by a damped oscillation.We quantify the oscillatory responses using two parameters, t 0.5 and p of oscillation for both models (Fig. 3).Both t 0.5 and p increase with (ε −1 -1)K s /V s in the two-pool model or (ε −1 -1)K l /V l in the three-pool model (Fig. 4).Because K s /V s or K l /V l decreases exponentially with an increase in soil temperature, and (ε −1 -1) increases with soil temperature, the optimal temperature at which t 0.5 and p are smallest is between 0 and 15 • C, depending on other parameter values.Therefore the two nonlinear models approach their respective equilibrium faster at a soil temperature between 0 and 15 • C than at other higher or lower temperatures.
Oscillatory responses of soil microbial biomass to nutrient or carbon inputs in the rhizosphere have been observed (Se-  (a-c) for an instant acclimation to warming with a microbial growth efficiency of 0.39 as calculated for a soil temperature of 15 • C, and the panels (d-f) for no acclimation to warming with a microbial growth efficiency of 0.31 as calculated for a soil temperature of 20 • C. α = 0.02 for all simulations here.The pink curves on panels (c) and (f) are the best fitted regressions to the simulation by the nonlinear model (dashed black curves).menov et al., 1999), in which the oscillation of soil microbial biomass carbon lasted for over a month, and was predominately driven by microbe-substrate interaction only in the first week or so (Zelenev et al., 2006), which is much shorter than the strong oscillation of multiple decades simulated by the nonlinear models here.The response of soil organic carbon to a step change in carbon inputs is oscillatory only in the two-pool nonlinear model, but is monotonic in the three-pool nonlinear model.Furthermore, oscillatory responses of soil carbon to perturbations under relatively constant conditions are yet to be observed.Instead, field observations support the simulated monotonic responses of soil carbon to a perturbation by the three-pool nonlinear model or conventional linear model, and are inconsistent with the simulated response by the two-pool microbial model (Adair et al., 2008;Zhang et al., 2008;Yang et al., 2011;Li et al., 2013).Therefore the two-pool nonlinear model may be limited to predicting microbial dynamics over relatively short timescales (< 1 year), and is probably not suitable for applications at regional or global scales over multiple years or longer, and therefore will not be discussed further.
When no carbon input enters the soil carbon pool directly (α = 0), the response of soil carbon to warming is independent of microbial growth efficiency.When α is small (< 0.05), the response of soil carbon to warming with temperature-dependent ε is weaker than that with temperature-independent ε.This explains why, in the studies of Allison et al. (2010) and Wieder et al. (2013), who used a value of 0.02 for α, the soil with temperature acclimation of ε was predicted to lose less carbon than the soil without temperature acclimation.
Mathematically there are two fundamental differences between the nonlinear microbial models and conventional linear models: (1) the rate of decay of the substrate (litter or soil carbon) is proportional to microbial biomass and varies with substrate concentration following the Michaelis-Menten kinetics in the nonlinear model, and is proportional to substrate concentration only in the conventional linear model; and (2) the microbial growth efficiency can vary with soil temperature in the nonlinear model, but typically does not vary with soil temperature in the conventional linear model.The first difference results in oscillatory responses to a small perturbation in the nonlinear model and monotonic responses in the conventional linear model, and the second difference results in different sensitivity of soil carbon to warming between the two types of models.
The nonlinear and linear models represent two paradigms of our diverging understanding of soil carbon dynamics at present.In the following we will discuss which one is better supported by the evidence available from the observed responses of soil carbon to carbon input or warming at timescales of years to decades.
By mathematical coincidence, the formulation in the nonlinear models gives very different sensitivities of equilibrium pool sizes to carbon influx from the conventional linear model.The equilibrium pool sizes of litter or soil carbon are insensitive to a change in carbon input in the nonlinear three-pool model, and are proportional to change in carbon input in the conventional linear model.Although the equilibrium microbial carbon pool size is linearly proportional to carbon influx, which agrees with observations (Fierer et al., 2009), soil microbial biomass only accounts for about 1-4 % of total soil carbon (Sparling, 1992;Serna-Chavez et al., 2013), therefore the sensitivity of total soil carbon (microbial biomass and soil organic carbon) to carbon influx in the nonlinear models is close to zero.
This difference in their sensitivity to carbon influx of the two types of models can be assessed using measurements from different ecosystems with different net primary production within a region, such as forests and grasslands, and the same ecosystems under different ambient [CO 2 ] treatments.Measurements of soil carbon change before and after conversion of forests to grasslands generally support the predictions by the conventional linear model (Emanuel et al., 1984;Wang et al., 1997;Murty et al., 2002).Plant invasion usually increases carbon input into ecosystems, leading to increased soil carbon storage (Liao et al., 2008).Several metaanalyses also showed that elevated CO 2 increased photosyn-thesis, plant production, litter mass, and soil carbon, particularly when nitrogen input to the ecosystem is high (Hungate et al., 2009).In contrast, leaf litter manipulation studies demonstrate that augmenting litter inputs does not necessarily increase soil C storage (Sulzman et al., 2005;Nadelhoffer et al., 2004;Leff et al., 2012).Elsewhere, a recent study found that mycorrhiza-mediated competition between plants and microbes, rather than NPP (net primary productivity), is the major driver of soil carbon storage (Averill et al., 2014).Therefore more studies are needed to determine the sensitivity of soil carbon storage to carbon input.
The simulated responses to soil warming by the nonlinear and conventional linear models also are quite different.This has attracted considerable attention recently (Allison et al., 2010;German et al., 2012;Wieder et al., 2013).The nonlinear model simulated a decrease in microbial biomass carbon and an increase in litter carbon when ε decreases with warming or no change in microbial biomass and a decrease in litter carbon if ε does not vary with warming; whereas the conventional linear models always simulate a decrease in both microbial biomass and litter carbon pools with warming, given everything else is independent of warming.These differences can be assessed by comparing changes in soil and microbial biomass carbon from various experiments and observations.Zhang et al. (2005) found that warming did not significantly change the soil microbial biomass amount but the species composition of the microbial community.As a result, microbial activities acclimated to the warming and soil respiration did not increase at the rate predicted by conventional linear models (Luo et al., 2001).The shifted microbial community toward more fungi and less bacteria under warming than in the control was found to be the dominant driver of warming acclimatization of soil respiration (Zhang et al., 2005).Furthermore, climate warming depletes labile carbon (Xu et al., 2012) and responses, such as changes in nutrient dynamics and use efficiency, ecohydrology and plant phenology.(see Luo, 2007;Melillio et al., 2011).It is still a challenge to realistically simulate this cascade of mechanisms underlying acclimation of soil carbon decomposition.
Theoretically microbial growth efficiency also varies with substrate quality and the composition of the soil microbial community (Frey et al., 2013;Sinsabaugh et al., 2013).The nonlinear model seems to be quite flexible to capture a range of the observed responses, and its simulated responses to warming also encompass the simulated responses by the conventional linear model if appropriate temperature sensitivity of ε is used in the nonlinear model.However the nonlinear model has an unrealistic sensitivity to carbon influx.Therefore a better model would have the sensitivity to carbon influx of the conventional linear model with the flexible response of the nonlinear model to warming.
In theory, the nonlinear soil carbon models represent the carbon decomposition process more realistically than the conventional model, but are yet to be tested against a wide range of field observations.Some key processes are still missing from the nonlinear soil carbon models, such as the dependence of key model parameters, such as V l , V s , K l , K s ε and µ, on soil moisture, substrate quality, and soil texture (see Parton et al. 1993) and interactions of carbon cycles with nitrogen and phosphorous cycles (see Wang et al., 2010).Although initial steps to address these needs are being made (e.g., Moorhead et al., 2012;Wang et al., 2013;Wieder et al., 2014), developing the appropriate scaling relationships between parameters and environmental drivers will be critical to constrain the response of nonlinear soil models to perturbations.Improving the models to better predict global soil carbon dynamics at decadal or centennial timescales might therefore require the incorporation of (1) some degree of sensitivity of soil carbon pool sizes to carbon input; (2) responses of microbial growth efficiency to soil temperature, warming, and dependence of key parameters on soil temperature, moisture, texture and substrate quality; (3) representing different microbial communities, such as fungi and bacteria; and (4) different sensitivities of different soil organic carbon pools to microbial attack.We then need to confront different model structures with empirical data focusing on (1) the spatial variation of soil carbon stocks at regional to global scales, and (2) response of soil carbon dynamics to disturbance (e.g., fire, land use change, warming, CO 2 enrichment) from longterm (> decade) field experiments (Guo and Gifford, 2002;Bradford et al., 2008;Hungate et al. 2009).Advanced modeldata fusion techniques can be used to identify better models for global change studies by comparing model predictions against empirical data (Wang et al., 2009;Smith et al., 2013).
In summary, our study has quantified how soil carbon dynamics is affected by the nonlinear interactions between substrate and decomposers, and what determines the stability of the soil carbon.We have also identified some key differences in the simulated responses of soil carbon to warming and carbon input between the nonlinear models and conventional linear models.These analyses will help us to develop better models of soil carbon decomposition and to eventually lead to more accurate predictions of global soil carbon under future climate conditions.
for the two-pool model The approximation in Eq. ( B9) is made because the size of soil organic carbon is at least two-orders of magnitude greater than microbial biomass carbon, and a 11 <0, a 31 <0 and a 32 <0.Numerical calculations also show that 0 < a 11 a 32 − a 31 < 1, therefore the amplitude of the oscillation of microbial biomass in the two-pool model is greater than that in the three-pool model.

Fig. 1 .
Fig. 1.Schematic diagrams showing the structure of two-pool and three-pool models.Carbon input enters the soil organic carbon (SOC) pool directly in the two-pool model, and litter carbon pool (LIT) in the three-pool model.When carbon from litter or soil organic carbon is consumed by soil microbes, a fraction of the consumed carbon, 1 − ε, is lost as CO 2 , where ε is microbial growth efficiency.See main text for the definitions of all other symbols.
4 g C m −2 for C b and 12650 g C m −2 for C s for the two-pool model, and 688 g C m −2 for C l , 50.4 g C m −2 for C b and 13170 g C m −2 for C s for the three-pool soil model.These values represent mean pool sizes within the global land biosphere.

Figure 2 .
Figure 2. Changes in microbial biomass in g C m -2 (upper panel) or soil organic carbon in g C m -2 770 (middle panel) (both in g C m -2 , 0 -100 cm soil depth) in time or soil organic carbon against 771 microbial biomass (lower panel) following a 10% reduction in initial pool sizes at time t=0 from their 772 respective equilibrium values for the two-pool soil model.Soil temperature was constant at 15 o C for 773 this simulation.774

Figure 3 .Fig
Figure3.Changes in litter carbon (upper panels), microbial biomass (middle panels) or soil org 777 carbon (lower panels) following 10% reduction in the initial pool sizes in microbial biomass an 778 organic carbon at time t=0 from their respective equilibrium pool sizes for the three-pool soil ca 779 model.The left three panels were for soil temperature of 15 o C and the right three panels for soi 780 temperature of 35 o C. The unit is g C m -2 for all pool sizes (0 -100 cm soil depth).781 782 Fig.3.Changes in litter carbon (upper panels), microbial biomass (middle panels) or soil organic carbon (lower panels) following a 10 % reduction in the initial pool sizes in microbial biomass and soil organic carbon at time t = 0 from their respective equilibrium pool sizes for the three-pool soil carbon model.The three panels on the left are for soil temperatures of 15 • C and the right panels for soil temperatures of 35 • C. The unit (g C m −2 ) is the same for all pool sizes (0-100 cm soil depth).

Figure 4 .
Figure 4. Variation in the two-pool soil model at a soil temperature of 15 o C of (a) half-time (t 0.5 ) with

Figure 6 . 800 Fig. 6 .
Figure 6.Temperature dependence of different terms of the half-life time (t0.5) or period for the two 797 oscillatory pools in the three-pool model.The solid and dashed curves represent the temperature 798 dependences with temperature-sensitive and temperature-insensitive microbial growth efficiency, 799 respectively.800 Fig. 6.Temperature dependence of different terms of t 0.5 or period for the two oscillatory pools in the three-pool model.The solid and dashed curves represent the temperature dependence with temperature-sensitive and temperature-insensitive microbial growth efficiency, respectively.

4. 4 Figure 7 .
Figure 7. Temperature dependence of the turnover rate of soil organic carbon decomposition with temperature-sensitive (solid curve) or temperature insensitive (dashed curve) microbial growth efficiency of soil microbes in the three-pool model.

Fig. 7 .
Fig. 7. Temperature dependence of the turnover rate of soil organic carbon decomposition with temperature-sensitive (solid curve) or temperature-insensitive (dashed curve) microbial growth efficiency of soil microbes in the three-pool model.

Figure 8 .Fig. 8 .
Figure 8. Normalized equilibrium responses of the litter 809 soil warming by the three-pool model with =0.The op

Figure 9 .Fig. 9 .
Figure 9. Responses of different carbon pools to a 5 o C warming above 15 o C at time t=0 of the 818 nonlinear (dashed black curves) or linear (yellow curves) versions of the three-pool models.Panels 819 (a)-(c) for an instant acclimation to warming with a microbial growth efficiency of 0.39 as calculat 820 for a soil temperature at 15 o C and the panels (d)-(f) for no acclimation to warming with a microbial 821 growth efficiency of 0.31 as calculated for a soil temperature of 20 o C. =0.02 for all simulations 822 here.The pink curves on panels (c) and (f) are the best fitted regressions to the simulation by the 823 nonlinear model (dashed black curves).824 825 Fig. 9. Responses of different carbon pools to a 5 • C warming above 15 • C at time t = 0 of the nonlinear (dashed black curves) or linear (yellow curves) versions of the three-pool models.Panels(a-c) for an instant acclimation to warming with a microbial growth efficiency of 0.39 as calculated for a soil temperature of 15 • C, and the panels (d-f) for no acclimation to warming with a microbial growth efficiency of 0.31 as calculated for a soil temperature of 20 • C. α = 0.02 for all simulations here.The pink curves on panels (c) and (f) are the best fitted regressions to the simulation by the nonlinear model (dashed black curves).