Responses of two nonlinear microbial models to warming and increased carbon input

A number of nonlinear microbial models of soil carbon decomposition have been developed. Some of them have been applied globally but have yet to be shown to realistically represent soil carbon dynamics in the field. A thorough analysis of their key differences is needed to inform future model developments. Here we compare two nonlinear microbial models of soil carbon decomposition: one based on reverse Michaelis–Menten kinetics (model A) and the other on regular Michaelis–Menten kinetics (model B). Using analytic approximations and numerical solutions, we find that the oscillatory responses of carbon pools to a small perturbation in their initial pool sizes dampen faster in model A than in model B. Soil warming always decreases carbon storage in model A, but in model B it predominantly decreases carbon storage in cool regions and increases carbon storage in warm regions. For both models, the CO2 efflux from soil carbon decomposition reaches a maximum value some time after increased carbon input (as in priming experiments). This maximum CO2 efflux (Fmax) decreases with an increase in soil temperature in both models. However, the sensitivity of Fmax to the increased amount of carbon input increases with soil temperature in model A but decreases monotonically with an increase in soil temperature in model B. These differences in the responses to soil warming and carbon input between the two nonlinear models can be used to discern which model is more realistic when compared to results from field or laboratory experiments. These insights will contribute to an improved understanding of the significance of soil microbial processes in soil carbon responses to future climate change.


Introduction
The dynamics of soil carbon in most global biogeochemical models are modelled using first-order kinetics, which assumes that the decay rate of soil carbon is proportional to the size of soil carbon pool. This approach has been recently questioned on theoretical grounds (Schimel and Weintraub, 2003;Fontaine and Barot, 2005), and is contradicted by the observed responses of soil carbon decay to the addition of fresh organic litter (Fontaine et al., 2004;Sayer et al., 2011) or soil warming (Luo et al., 2001;Melillo et al., 2002;Bradford et al., 2008). As a result, a number of nonlinear soil microbial models have been developed (Alli-Published by Copernicus Publications on behalf of the European Geosciences Union.

888
Y.-P. Wang et al.: Responses of two nonlinear microbial models to warming and increased carbon input son et al., 2010;Manzoni and Porporato, 2007;Wutzler and Reichstein, 2008) and a few of them have been applied at global scales (Wieder et al., 2013;Sulman et al. 2014). Predictions of future soil carbon change by these nonlinear models can differ significantly from conventional linear models (Fontaine et al., 2007;Wieder et al., 2013). For example, conventional linear soil carbon models predict that soil carbon will decrease with increased temperature, all else being equal (Jenkinson et al., 1991), whereas the nonlinear models predict that the soil carbon can decrease or increase, depending on the temperature sensitivity of microbial growth efficiency and turnover rates (Frey et al., 2013;Hagerty et al., 2014;Li et al., 2014). However, the nonlinear models have yet to be validated against field measurements as extensively as the conventional linear soil carbon models (Wieder et al., 2016). They also have some undesirable features, particularly the presence of strong oscillations or bifurcations (Manzoni and Porporato, 2007;Wang et al., 2014) in their dynamics that are not observed in real-world systems. Therefore it is important to improve understanding of the behaviour of these nonlinear models before they are used in earth system models for informing climate decisions.
Nonlinear microbial models can explain why the decomposition rate of recalcitrant organic soil carbon varies after the addition of easily decomposable organic carbon to soil, which is known as the priming effect (Kuzyakov et al., 2000). This response has been observed in the field (Fontaine et al., 2004;Sayer et al., 2011) but cannot predicted by conventional linear soil carbon models without modification (Fujita et al., 2014). Theoretically, decomposition of soil organic carbon is catalysed by extracellular enzymes that are produced by soil microbes. The production rate of extracellular enzymes depends on the biomass and composition of the soil microbial population and their local environment. Therefore the decomposition rate of soil organic carbon should depend on both microbial biomass and substrate concentration (Schimel and Weintraub, 2003), rather than on substrate concentration only, as assumed in conventional linear models.
This sensitivity of soil carbon decomposition to the input of additional carbon has important implications for the storage of carbon by the biosphere in response to climate change. Soil is the largest land carbon pool and therefore the direction and magnitude of the global carbon-climate feedback strongly depends on the responses of soil carbon to future warming (Jones and Fallow, 2009;Hargety et al., 2014).
A number of nonlinear models have been developed that explicitly account for the dynamics of the soil microbial community (Parnas, 1978;Smith, 1979;Schimel and Weintraub, 2003;Wutzler and Reichstein, 2008;Allison et al., 2010;Grant, 2014;Riley et al., 2014;Tang and Riley, 2014). Parnas (1979) explored the mechanism of priming using a nonlinear soil microbial model that included both soil carbon and nitrogen dynamics. Smith (1979) developed a nonlinear model of soil carbon decomposition that included the interactions among carbon, nitrogen, phosphorus, and potas-sium. Smith's model represented multiple forms of carbon, nitrogen, and phosphorus and their transformation via abiotic (such as adsorption and desorption) and biological processes by different groups of soil microbes. The soil models developed by both Parnas (1978) and Smith (1979) were based on regular Michaelis-Menten kinetics, in which the rate of carbon decomposition depends linearly on the concentration of soil enzymes but nonlinearly on substrate concentration (Roberts, 1977). This was challenged by Schimel and Weintraub (2003), who emphasized the importance of exoenzyme limitation on soil carbon decomposition. Schimel and Weintraub (2003) used a reverse Michaelis-Menten kinetics formulation to show that the response of soil carbon decomposition to carbon substrate concentration can be nonlinear regardless of carbon supply. The reverse Michaelis-Menten kinetics for soil carbon decomposition assumes that the rate of carbon decomposition depends nonlinearly on enzyme concentration but linearly on substrate concentration.
The nonlinear soil carbon models described above have subsequently been used in a variety of studies: to explore different the fundamental mechanisms controlling soil carbon decomposition (Schimel and Weintraub, 2003, for example), to investigate the sensitivity of soil carbon and other biogeochemical processes to warming (Grant, 2014;Tang and Riley, 2014), to investigate the response of soil carbon to a small perturbation, such as priming (Wutzler and Reichstein, 2013), and to predict soil carbon responses to global change (Wieder et al., 2013;Sulman et al., 2014). Some studies have explored the mathematical properties of these nonlinear models in detail (for example, Manzoni et al., 2004;Manzoni and Porporato, 2007;Raupach, 2007;Wang et al., 2014). However, to date these have been predominantly restricted to obtaining insights for individual models and with a specific parameterization.
In this study we use mathematical analysis to improve our understanding of the key properties of nonlinear microbial models. For simplicity and analytic convenience, we choose two simple types of nonlinear microbial models: one with regular Michaelis-Menten kinetics and another with the reverse Michaelis-Menten kinetics. These models can be considered as two special cases of the more general kinetics discussed by Tang (2015). These two simple formulations are amenable to analytic approximations, whereas the formulations with more general kinetics, such as the equilibrium chemistry approximations, are not. We only represent three soil carbon pools with each model and ignore abiotic processes for simplicity, despite these being potentially important under certain conditions (see Tang and Riley, 2014 for an example). In comparing the two nonlinear microbial models, we use the standard mathematical technique to analyse their responses to a small perturbation (see Wang et al., 2014), such as a step change in soil temperature or carbon input, or whether two models exhibit oscillatory behaviour under certain conditions, and how the analytic approximations to the exact model solutions differ between the two nonlinear mod-Y.-P. Wang et al.: Responses of two nonlinear microbial models to warming and increased carbon input 889 els. We address the following questions. (1) How do the responses of these two models to soil warming differ, and why?
(2) Can both models simulate the response of soil carbon decomposition to increased carbon input as in a priming experiment and what determines the magnitude of the response in each model?

Model description
We consider two nonlinear soil microbial models: model A, which uses reverse Michaelis-Menten kinetics, and model B, which uses regular Michaelis-Menten kinetics (specified below). Both models have three carbon pools: litter carbon, microbial biomass, and soil carbon.
Model A is based on the nonlinear microbial model of soil carbon described Wutzler and Reichstein (2013, their model A1). Their original model has four pools, modelled by and dC m dt = (µ l C l + µ s C s ) where t is time in years; C l , C s , C b , and C m represent the pool sizes of litter carbon, soil carbon, microbial biomass carbon, and assimilable soil carbon in g C m −2 , respectively; and F npp is carbon input in g C m −2 yr −1 , with the fraction a going to the soil carbon pool, and (1 − a) to the litter carbon pool. µ l , µ s , µ b , and µ m are rate constants of litter carbon, soil carbon, microbial biomass, and assimilable carbon per year, respectively (see Schimel and Weintraub, 2003); ε is microbial growth efficiency; and K b and K m are two empirical constants in g C m −2 for the dependence of the consumption of litter carbon or assimilable carbon by soil microbes.
In this study we are interested in the responses at timescales greater than 1 year. We therefore assume that C m is at steady state (dC m / dt = 0) because of its relatively fast turnover (less than a few days). Therefore the dynamics of microbial biomass, C b , can be simplified to Model A as used in this paper consists of Eqs. (1), (2), and (5) unless otherwise specified. This type of formulation was also used by Schimel and Weintraub (2003) and Drake et al. (2013).
Model B, based on the model used by Allison et al. (2010) and Wieder et al. (2013) with one additional assumption that both enzyme and dissolved organic carbon pools are at steady states, is given by and where K l and K s are Michaelis-Menten constants in g C m −2 , and V l and V s are maximum rates of substrate carbon (litter or soil) assimilation rate per unit microbial biomass per year. This type of kinetics was used by Riley et al. (2014), Wieder et al. (2014) and Wang et al. (2014). These two models make different assumptions about the rate-limiting step in carbon decomposition. Both models assume that microbes have similar access to litter and soil carbon. In model A, carbon decomposition is assumed to depend nonlinearly on the number of binding sites or the amount of substrate and linearly on enzymes or microbial biomass (Schimel and Weintraub, 2003). In model B, carbon decomposition is assumed to depend nonlinearly on enzymes or microbial biomass and linearly on the number of binding sites or the amount of substrate (Allison et al., 2010).
When carbon input, F npp , is equal to zero, the steady-state solution is zero for litter and soil carbon pools for both models (a trivial solution). When F npp > 0, the steady-state solutions to model A are and The steady-state solutions to model B are and

890
Y.-P. Wang et al.: Responses of two nonlinear microbial models to warming and increased carbon input CO 2 efflux from the decomposition of soil organic carbon (F s ) is calculated as for model A and for model B.

Parameter values
We allow all model parameters to vary with soil temperature (T s ) with the exception of parameter a. Based on the work of Allison et al. (2010) and Hagerty et al. (2014), we model the temperature dependence of parameters as and for both models, where T R is reference soil temperature in • C (i.e. 15 • C), ε R and µ bR are the values of ε and µ b at T s = T R , respectively, and x and b are two empirical constants (see Table 1 for their default values). Previously there has been debate about the temperature sensitivities of ε and µ b (see Frey et al., 2013;Hargety et al., 2014). The microbial models as developed by Allison et al. (2010) and used by Wieder et al. (2013), and  assumed that ε was temperature-sensitive and µ b was temperature-insensitive (or b = 0). This assumption was recently challenged by Hargety et al. (2014), who found that µ b was temperature-sensitive and ε was not, based on a laboratory soil-warming experiment. Here we will explore the consequence of different assumptions about the temperature sensitivities of ε and µ b on the simulated response of soil carbon to warming by the two models (see Sect. 3.2).
We also assume that three additional model parameters in model A, K b , µ l , and µ s depend on soil temperature exponentially, with and where K bR , µ lR , and µ sR are the values of K b , µ l and µ s when soil temperature (T s ) is equal to the reference temperature, T R (15 • C in this study), and α k , α l , and α s are three empirical constants with their default values listed in Table 1.
For model B, we assume that K l , K s , V l , and V s increase with soil temperature exponentially: and where K lR , K sR , V lR , and V sR are the values of K l , K s , V l , and V s at the reference soil temperature (T R ), respectively, and β kl , β ks , β vl and β vs are four empirical constants for model B (see Table 1). As found by Wang et al. (2014), the microbial biomass as simulated by model B using the parameter values of Wieder et al. (2013) was low (< 1 % of total soil carbon). We therefore reduced the turnover rate of microbial biomass to 1.1 yr −1 by assuming that 2 % of total soil organic carbon is microbial biomass carbon at a soil temperature of 15 • C. Some parameter values in model A at the reference temperature were obtained by calibrating the equilibrium litter and soil carbon pool sizes against those from model B for a soil temperature of 15 • C and carbon input of 400 g C m −2 yr −1 , as used in Wang et al. (2014).

Analytic solutions and numerical simulations
We derived and used analytic solutions whenever possible for comparing the two models. Specifically, we mathematically analysed the temperature dependence of steady-state soil carbon pool size, and derived an analytic approximation of soil temperature at which equilibrium soil carbon is at a minimum (e.g. Eq. B4 for model B). We also derived an approximate solution for the maximum CO 2 loss from soil carbon decomposition after the increased carbon input for each model (e.g. Eq. C12 for model A and Eq. C15 for model B). When an analytic solution was not possible or too cumbersome, we used numerical simulations to show the differences between the two models in their responses of carbon pools to a small perturbation in litter or microbial carbon pool sizes, and the response of CO 2 efflux from soil carbon decomposition to litter addition at a tropical forest site (Sayer et al., 2011).

Results
Before comparing the responses of our models to soil warming and increased carbon input, we first analyse some key properties of their responses to a small perturbation, i.e. whether both models oscillate in response to a small change in their initial pool sizes and what determines the period and amplitude of the oscillation. As a step change in soil temperature or carbon input can be considered to be a perturbation, Table 1. Default values of model parameters and their temperature sensitivities ( • C −1 ). Four parameters were tuned: 1 tuned using the microbial biomass data measured from a tropical forest site (see Sayer et al., 2011), and 2 tuned against the soil carbon pool size simulated by model B by Wang et al. (2014). identifying differences in those key properties will help us understand the differences in the responses of the two models to soil warming and increased carbon input. The response of model B to perturbation has already been analysed by Wang et al. (2014), and will not be elaborated here, but the results from that analysis will be used to compare the period and amplitude of the response to perturbation to that of model A.

Comparison of the perturbation responses of both models
Perturbation analysis is a standard mathematical technique for analysing the behaviour of a dynamic system near its equilibrium state (see Drazin, 1992, for further details). There are two kinds of perturbation responses: stable or unstable. The system states, or carbon pool sizes in this study, will always approach their equilibrium states for a stable response, or otherwise for an unstable response. For both stable and unstable responses, the transient change in a carbon pool size over time can be oscillatory or monotonic. As shown in Appendix A, the response of a carbon pool to a small perturbation is always stable, and oscillatory only if F npp < or monotonic otherwise for model A. This region of oscillation in the two-dimensional space of carbon input and soil temperature is shown in black in Fig. 1. The response of model A to a small perturbation is oscillatory under most conditions; the conditions with low soil temperature and high carbon input are uncommon in terrestrial ecosystems.
The results of a singular perturbation analysis are strictly applicable only when the perturbation is small. However, our simulations show that the predictions from the perturbation analysis approximate well the responses of our two models to any realistic perturbation (see Appendix A of this paper and Appendix B in Wang et al., 2014). Therefore we can predict how soil carbon or other carbon pools change over time in response to a change in carbon inputs or soil warming (i.e. a perturbation of the external environment) and explain why the responses of the carbon pools are different between the two models.
To illustrate how the responses of carbon pools to a small perturbation differ between the two models, we numerically simulated the recovery of all three carbon pools in each model after a 10 % reduction at time t = 0 in both litter and microbial carbon from their respective steady-state values, while no perturbation was applied to soil carbon at t = 0 (see Fig. 2). The amplitude of the initial oscillation is about 70 g C m −2 for the litter pool (see Fig. 2b) and 7 g C m −2 for the microbial carbon pool (see Fig. 2d) in model B, compared to about 25 g C m −2 (see Fig. 2a) for the litter pool and 4 g C m −2 for the microbial pool (see Fig. 2c) in model A. After 20 years, both the litter and microbial carbon pools are very close to their respective steady-state values in model A, but continue to oscillate in model B.
The oscillatory response can be mathematically characterized by its half-life (t 0.5 ) and period (p). For a stable oscillatory response, the amplitude of the oscillation decays ex- ponentially. The time for the amplitude to reach 50 % of its initial value is defined as the half-life (t 0.5 ). The smaller t 0.5 is, the faster the oscillation dampens. As explained in Appendix A, values of t 0.5 and p for model A are much smaller than model B for any given soil temperature and perturbation. This explains why the oscillatory response of model A dampens much faster than model B.
There are significant differences in the response of soil carbon between the two models. While there is no response of soil carbon to a small perturbation in litter carbon and microbial biomass in model B, soil carbon in model A decreases initially to a minimum value at 5 years after the perturbation, then gradually increases to its steady-state value. These differences in the response of soil carbon between the two models can be explained by the differences in the structure of eigenvectors for litter carbon and microbial biomass between the two models (see Appendix A for further details).

Response of soil carbon to warming
Here we explore how soil carbon responds to a step increase in soil temperature, as in many soil-warming experiments (Luo et al., 2001;Mellilo et al., 2002), and ignore the response of carbon input to warming.
As explained in Appendix A, the response of soil carbon to warming is always stable in both models and is likely to be weakly oscillatory in model A and monotonic in model B. The transient change in soil carbon after warming can be pre-dicted using the generalized solution for soil carbon for each model (see Eq. B1 of Wang et al., 2014). Therefore the directional change in soil carbon in response to warming, i.e. increasing or decreasing only, depends on the sensitivity of the equilibrium soil carbon pool to soil temperature in both models.
As shown in Appendix B, the equilibrium pool size of soil carbon of model A always decreases with soil warming if carbon input does not increase with warming. For model B, the equilibrium pool size of soil carbon can increase or decrease in response to warming, depending on soil temperature and model parameter values. In Appendix B, we show that a soil temperature (T x ) may exist at which the equilibrium soil carbon is at a minimum for model B. Identifying T x is important for predicting the directional change in soil carbon by model B in a warmer world, because soil carbon will decrease if the warmed soil temperature is below T x , and will increase otherwise.
The value of T x for model B depends on three parameters: the fraction of carbon input directly into the soil pool (a), microbial biomass turnover rate (µ b or its temperature sensitivity b), and microbial growth efficiency (ε or its temperature sensitivity x). Figure 3a shows that T x for model B decreases with an increase in a or x. Over the ranges of values of x and a, T x can vary across the range of air temperature experienced by most terrestrial ecosystems. For example, T x is > 40 • C when x < 0.005 • C −1 and a < 0.5; therefore the equilibrium soil carbon predicted by model B decreases with warming when the warmed soil temperature is below 40 • C. When a > 0.4 and x>0.02 • C −1 , T x is < 0 • C (the black region on the top left corner of Fig. 3a); therefore the simulated equilibrium soil carbon by model B increases with warming if the warmed soil temperature is above 0 • C. Figure 3b shows that T x for model B decreases with an increase in b or x. When the turnover rate of microbial biomass is not sensitive to soil temperature (b = 0) and x = 0.016 • C −1 as the default value for model B, T x is about 35 • C. For b = 0.063, as estimated by Hagerty et al. (2014), T x < 0 • C; therefore the equilibrium soil carbon pool size as simulated by model B always increases with soil warming for most terrestrial ecosystems, irrespective of the value of x.
Therefore the simulated responses of the soil carbon pool to warming by the two models can be quite different: the equilibrium soil carbon pool size always decreases with soil warming in model A, but can increase or decrease in model B, depending on the temperature sensitivities of microbial growth efficiency and microbial turnover rate and the fraction of carbon input entering soil carbon pool directly.

The response of soil carbon to an increased litter input
We compare the simulated responses of soil carbon to litter addition by the two models with field measurements from an experiment described by Sayer et al. (2011). used three treatments: litter removal (L − ), with aboveground litter being removed regularly; increased litter input (L + ) with the added litter from the litter removal treatment; and a control (C). Measurements of CO 2 efflux from soil were made and the contribution of root-rhizosphere respiration to soil respiration was estimated using a δ 13 C technique. Sayer et al. (2011) found that the CO 2 efflux from the decomposition of soil organic carbon in the L + treatment was 46 % higher than in the control. Therefore, increased litter addition accelerated the decomposition of soil organic carbon. Here we assess whether the observed response of soil carbon decomposition to increased litter input can be reproduced by running both models for L + and C treatments. Inputs to each model, including the monthly data of soil temperature and litter input from 2002 to 2008 for two treatments (C and L + ) at the site, were compiled from Sayer and Tanner (2010a, b; see Fig. 4 for monthly litter input as an example). We also assumed that the contribution of fine-root respiration to total soil respiration (root respiration plus heterotrophic respiration) was 35 % for the control treatment and 21 % for the litter addition treatment, based on the estimates by Sayer et al. (2011).
The initial sizes of all pools were obtained by running each model with the monthly inputs for the first 2 years repeated until all pools reached steady state (i.e. the change in pool size between two successive cycles is less than 0.01 %).
Using the initial pool sizes for each model and the monthly input from 2002-2008, we numerically integrated both models and calculated the average contributions to total soil CO 2 efflux from the decomposition of litter and soil organic carbon for the last 2 years (2007)(2008) and compared the simulated results with the estimates from field measurements by Sayer et al. (2011).
By tuning values of two model parameters (µ bR and K bR ) (see Table 1), we obtained an initial microbial biomass carbon 240 g C m −2 for both models, very close to the measured microbial biomass carbon of 219 g C m −2 by Sayer et al. (2007). The simulated initial soil carbon is 6715 g C m −2 for model A and 6945 g C m −2 for model B, which is higher than the estimated soil carbon of 5110 g C m −2 in the top 25 cm (Cavelier et al., 1992) and lower than the estimated soil carbon of 9272 g C m −2 in the top 50 cm soil (Grimm, 2007).
The estimated total soil CO 2 efflux from the control treatment by Sayer et al. (2011) was 1008 g C m −2 yr −1 from 2007 to 2008, which was closely simulated by both models (1004 g C m −2 yr −1 by model A and 1008 g C m −2 yr −1 by model B). However, both models overestimated the total soil CO 2 efflux from the litter addition treatment.  B (b). The dark-grey bar and black bars represent CO 2 effluxes from litter and soil organic carbon decomposition, respectively. The light-grey bar for the litter addition treatment represents the additional CO 2 efflux from soil organic carbon decomposition due to additional litter input. mated efflux by Sayer et al. (2011) was 1380 g C m −2 yr −1 , as compared with the simulated flux of 1425 g C m −2 yr −1 by model A and 1502 g C m −2 yr −1 by model B (see Fig. 5).
The additional CO 2 efflux from the decomposition of soil carbon in the litter addition treatment was estimated to be 180 ± 50 g C m −2 year −1 by Sayer et al. (2011), which was quite well simulated by model B (105 g C m −2 yr −1 ) (see Fig. 5b) but was underestimated by model A (29 g C m −2 yr −1 ) (see Fig. 5a).
The difference in the simulated response of soil organic carbon decomposition to increased litter input by the two models can be explained by differences in their substrate kinetics. The rate of carbon loss from the decomposition of soil carbon depends on both soil carbon and microbial biomass in both models. Because soil carbon is unlikely to change significantly within a few years, the rate of CO 2 emission from soil carbon decomposition will largely depend on microbial biomass, and that dependence is nonlinear following the reverse Michaelis-Menten equation in model A (see Eq. 2), but is linear in model B (see Eq. 7). Therefore the simulated response of soil organic carbon decomposition to increased litter input by model B is more sensitive to microbial biomass than model A.

Response to priming: maximum CO 2 efflux from soil carbon decomposition
Results from the above comparison of the responses of two models to the increased litter input are likely dependent on soil temperature, carbon input, and model parameter values.
To understand the differences in the responses of our two models to litter addition at different rates and soil temperatures for any parameter value, we use the analytic approximations to maximum CO 2 efflux from the priming treatment for each model to identify key differences in their response to priming.
Priming is defined as the change in organic carbon decomposition rate after the addition of an easily decomposable organic substance to soil (Kuzyakov et al., 2000). In lab priming experiments, a given amount of isotopically labelled C substrate is added to the primed treatment only at the beginning of the experiment (t = 0) and no substrate is added to the control. CO 2 effluxes from soil carbon decomposition are estimated from measurements for the following weeks or longer (Cheng et al., 2014). The effect of priming, p, is calculated as (R p −R c )/R c , where R c and R p are the CO 2 efflux from the decomposition of soil organic carbon in the control and primed treatments, respectively. Maximum values of p are usually reported in most priming studies (see Cheng et al., 2014).
However, analytic approximations to p for both models are quite cumbersome for analysing their differences in the responses to priming. Another way to quantify the priming effect is by measuring the maximum CO 2 efflux from soil organic carbon decomposition after carbon addition at time t = 0 (Jenkinson et al., 1985;Kuzyakov et al., 2000). This quantity can be easily measured in the laboratory or field.
In both models, the equilibrium soil microbial biomass is proportional to carbon input (see Eqs. 11 and 13). In the primed treatment, the amount of carbon added at t = 0 usually is well above the rate of the carbon input under natural conditions, and no further carbon is added. Therefore the microbial biomass will increase until reaching a maximum value, then decreases with time after t = 0.
As shown in Appendix C, the maximum CO 2 efflux from soil carbon decomposition in the primed treatment, F max , depends on the maximum microbial biomass and microbial growth efficiency for both models, as well as on soil carbon turnover rate for model A (see Eq. C12 for F A ) and the microbial turnover rate for model B (see Eq. C15 for F B ). Figure 6 shows that F max (or F A for model A, F B for model B) increases with carbon input, and decreases with an increase in soil temperature for both models. However, the sensitivity of F max to carbon input at different soil temperatures is different between the two models. For model A, the sensitivity of F max to carbon input is greatest around 25 • C, and is quite small at < 5 • C. For model B, the sensitivity of F max to carbon input decreases with an increase in soil temperature (see Fig. 6).
The sensitivity of F max to soil temperatures in both models can be explained by the analytic approximations (Eq. C12 for model A and C15 for model B). Maximum CO 2 efflux is proportional to soil carbon in model A, and to the maximum microbial biomass in model B. Both soil carbon and maximum microbial biomass in both models decrease with an increase in soil temperature for the parameter values we used (see Fig. 6c); therefore F max also decreases with an increase in soil temperature.
Differences in the sensitivity of F max to carbon input at different soil temperatures in the two models can also be explained by their respective analytic approximations, par- ticularly the dependence of maximum microbial biomass on both carbon input and initial microbial biomass in model A (see Eq. C11) and on equilibrium litter carbon pool size in model B (see Eq. C14), because F max depends on the maximum microbial biomass in both models. In model A, F A nonlinearly varies with maximum microbial biomass (see Eq. C12), which increases linearly with carbon addition at t = 0 ( C l ) and varies nonlinearly with the initial pool size of microbial biomass (C * b ) (see Eq. C11). Because C * b increases with a decrease in soil temperature or an increase in C l (see Fig. 6c), F A increases with an increase in C l (either directly (Eq. C11) or via the effect on C * b ), and with a decrease in soil temperature (via the temperature dependence of C * b ). In model B, the sensitivity of F B to carbon input is determined by the maximum microbial biomass (C b max, B ), which varies with equilibrium litter pool size (C * l ) following the regular Michaelis-Menten equation (C b max, B ∝ M l in Eq. C14) for a given amount of carbon input ( C l ). The equilibrium litter carbon pool size increases with soil temperature, and is independent of carbon input based on Eq. (12) (see Fig. 6d). When soil temperature is low, C * l is low, and therefore sensitivity of F B to carbon input is high. When soil temperature is high, C * l is high and the sensitivity of F B in model B to carbon input is low because of saturating response in the regular Michaelis-Menten equation.

Discussion
Here we analysed the responses of different carbon pools to perturbation, soil warming, and increased carbon input in two nonlinear microbial soil carbon models. Table 2 lists the key differences in those responses.
Some of the differences between the two models also depend on the chosen parameter values for each model. For example, there has been debate about the temperature sensitivities of microbial biomass turnover rate and microbial growth efficiency (Frey et al., 2013;Hargety et al., 2014), and the simulated sensitivity of soil carbon to warming (Hagerty et al., 2014). Regardless of the temperature sensitivity of microbial growth efficiency, model A always simulates a decrease in the equilibrium soil carbon under warming, whereas model B can simulate an increase or a decrease in the equilibrium soil carbon under warming, depending on the temperature sensitivities of microbial growth efficiency and turnover rate. If microbial growth efficiency is sensitive to soil temperature and microbial turnover rate is not, as found by Frey et al. (2013), the simulated responses of equilibrium soil carbon to warming by the two nonlinear models are quite similar in the direction of response over temperate and boreal regions, but different in the tropical regions. This is because the minimum soil carbon temperature, T x , for model B is about 25 • C for x = 0.015 K −1 and a = 0.05, the values used by Allison et al. (2010) and German et al. (2012) (see Fig. 3a). In that case the equilibrium soil carbon, as simulated by model B, will decrease over most temperate and boreal regions, for which the mean soil temperature within the rooting zone is below 25 • C for most of the growing season, and will increase in tropical regions, for which the mean soil temperature in the top 100 cm of soil is close to 25 • C for most of the year. However, if microbial turnover rate is sensitive to soil temperature and microbial growth efficiency is not, as found by Hargety et al. (2014), then T x is < 0 • C at α s > 0.055 ( • C) −1 for model B, causing equilibrium soil carbon to increase in model B with warming, but decrease in model A with warming. Therefore, the predicted responses of soil carbon to warming by the two nonlinear models differ significantly across all major global biomes where mean rooting zone soil temperature over the growing season is above 0 • C.
Some of the key differences in the responses of the two nonlinear models can be used to discern which model is more applicable to the real world. For example, the oscillatory response of model A generally is quite small (<1%), which is quite consistent with the results from litter removal experiments (Sayer et al., 2007, for example). The relatively large and more persistent oscillation in model B has not been observed in the field, and the insensitivity of soil carbon to 896 Y.-P. Wang et al.: Responses of two nonlinear microbial models to warming and increased carbon input a perturbation in the litter or soil microbial carbon pool in model B also needs to be assessed against long-term field experiments such as the DIRT experiment (Nadelhoffer et al., 2004). Model B in its present form may not be applicable under field conditions. It has been argued that the influences of microbial community structure and their activities on mineral soil carbon decomposition at field scale may be much smaller than at the rhizosphere scale (Schimel and Schaeffer, 2012), because substrate concentration rather microbial activity is the rate-limiting step for the decomposition of soil organic matter in mineral soils. A recent study by Sulman et al. (2014) clearly showed the importance of physical protection of microbial by-products in forming stable soil organic matter, and its implications for the response of global soil carbon to carbon inputs. This mechanism has been recently incorporated into a nonlinear soil microbial carbon model . Whether the large oscillatory responses of model B will be significantly dampened by the addition of such physical protection mechanism is yet to be studied.
The two models also have quite different sensitivities to soil warming (see Table 2), particularly in warm regions. Results from a decade-long soil-warming experiment showed that warming did not reduce soil carbon, because plant carbon production increased as a result of the increased availability of soil mineral nitrogen in a nitrogen-limited forest (Melillo et al., 2002). However, this is quite a different mechanism because model B in our study includes neither a nitrogen cycle nor the response of carbon input to warming.
Overall both models can simulate the priming response to a change in carbon inputs, although model A simulates a weaker response than model B and the sensitivities to carbon input at different soil temperature are different between the two models, particularly under cool climate conditions (see Table 2). So far, results from litter manipulation experiments in the field have not been analysed for their sensitivity to soil temperature. The differences in the responses of soil carbon decomposition to an increased carbon input we identified between the two models can also be used to assess which model is more applicable in the field using experiments with different carbon input under cool (mean annual air temperature < 10 • C) and warm (mean annual temperature > 20 • C) conditions. If the sensitivity of soil carbon decomposition to an increased carbon input under cool conditions is greater than that under warm conditions, then model B is more appropriate than model A.
This has yet to be tested.
Our analysis here does not include some other key processes, such as the transformations of different forms of organic carbon substrates by different microbial communities as included in some models (see Grant, 2014;Riley et al., 2014, for example). Therefore the conclusions from this study about the two nonlinear models should be interpreted with some caution. As shown by Tang and Riley (2014), interactions among soil mineral sorption, carbon substrate, and microbial processes can generate transient changes in the apparent sensitivity of soil carbon decomposition to soil temperature; therefore the static dependence of microbial processes on soil temperature as used in our study may not be applicable. Our simplification of the soil microbial community and soil carbon fractions is necessary for analytic tractability, but may also limit the applicability of our results to field experiments. For example, Allison (2012) showed that the apparent kinetics of soil carbon decomposition can vary with the spatial scale: the regular Michaelis-Menten kinetics at microsites coupled with an explicit representation of different strategies for facilitation and competition among different microbial taxa generated litter carbon decomposition kinetics similar to the reverse Michaelis-Menten equation. Therefore, the identified differences between the two models should vary with spatial scale.
The regular and reverse Michaelis-Menten kinetics can be considered as two special cases of a more general kinetics, as discussed by Tang (2015). Both models use different mass balance constraints (see Tang, 2015), which are unlikely to hold across a wide range of conditions. In the real world, the kinetics and parameter values of carbon decomposition likely depend on a number of other factors, such as soil physical properties, substrate quality, and soil nutrient availability (Manzoni and Porporato, 2009). Future studies of soil carbon decomposition kinetics need to include those factors and the role of root growth dynamics and photosynthetic activities in rhizosphere priming (see Kuzyakov, 2002).  , 2010). The values of those parameters may be quite different under field conditions. Evaluation of their applicability under a wide range of field conditions will require an integrated approach, such as applications of model-data fusion using a range of field experiments (Wieder et al., 2016). This will eventually lead to a better understanding of the significance of microbial activity on soil carbon decomposition and more accurate predictions of carbon-climate interactions.

Conclusions
This study analysed the mathematical properties of two nonlinear microbial soil carbon models and their responses to soil warming and carbon input. We found that the model using the reverse Michaelis-Menten kinetics (model A) has shorter and more frequent oscillations than the model using regular Michaelis-Menten kinetics (model B) in response to a small perturbation.
The responses of soil carbon to warming can be quite different between the two models. Under global warming, model A always simulates a decrease in soil carbon, but model B will likely simulate a decrease in soil carbon in temperate and boreal regions, and an increase in soil carbon in tropical regions, depending on the sensitivities of microbial growth efficiency and microbial biomass turnover rate.
The response to carbon input varies with soil temperature in both models. The simulated maximum response to priming by model A generally is smaller than that by model B. The maximum rate of CO 2 efflux from SOC decomposition (F max ) to carbon input in the primed treatment decreases with an increase in soil temperature in both models, and the sensitivity of F max to the amount of carbon input increases with soil temperature in model A but decreases monotonically with an increase in soil temperature in model B.
Based on those differences between the two models, we can design laboratory or field experiments to assess which model is more applicable in the real world and, therefore, advance our understanding of the importance of microbial processes at regional to global scales. The Jacobian at the equilibrium pool sizes, J , is given by where a 1 = µ l g, a 2 = µ s g, a 3 = µ l C * l ∂ g and C * s are the equilibrium pool sizes of litter carbon, microbial biomass, and soil carbon in g C m −2 , respectively.
The three eigenvalues of J are given by where . These correspond to three carbon pools (λ 1 for litter carbon, λ 2 for microbial biomass, and λ 3 for soil carbon). If the eigenvalue of a carbon pool is complex, then the response of that pool to a small perturbation is oscillatory, or monotonic otherwise. If the real part of the eigenvalue is negative, then the response is stable. Therefore, the responses of all three carbon pools to a small perturbation are monotonic if F > 0 or F npp > 4 (1−ε) 2 ε µ l µ 2 b (µ b −µ l ) 2 K b , or oscillatory otherwise (or F < 0). The responses of all carbon pools are always stable because The corresponding eigenvectors of J are given by When the responses of carbon pools to a small perturbation are oscillatory and stable, the amplitude of oscillation decreases exponentially after t = 0. The oscillatory response can be characterized by its half-life (t 0.5 ) and period (p) (both in years) calculated from their eigenvalues. The amplitude of a stable oscillation decreases exponentially over time, and time when the amplitude is half as much as the amplitude at t = 0 is defined as t 0.5 . t 0.5 and p are calculated for model A as Wang et al. (2014) gave the formulae for t 0.5 and p for model B (their Eqs. 24 and 25). As shown in Fig. A1, the half-life is longest for both models when soil temperature is high and carbon input is low, conditions often experienced in arid ecosystems, implying a strong oscillation at these conditions. At a given soil temperature and carbon input, the half-life for model A is about half as much as that for model B (see Appendix Fig. A1a and  b). When carbon input is > 1000 g C m −2 yr −1 , as in tropical rainforests, the half-life is less than 1 year for model A at a soil temperature between 20 and 30 • C, and for model B at a soil temperature between 0 and 20 • C only.
Over the range of realistic carbon inputs and soil temperatures, the values of both t 0.5 and p of model A are less than half as much as those of model B (see Fig. A1). Therefore the responses of carbon pool sizes to a small perturbation in model A oscillate faster and those oscillations also dampen faster than model B.
As shown by Wang et al. (2014) (their Appendix B, Eq. B1; there a i is eigenvalue and v i is eigenvector), the evolution of each carbon pool after a small perturbation can be mathematically represented using the eigenvalues, eigenvectors and initial pool sizes (Eq. B1 in Appendix B of . The third elements of the eigenvectors corresponding to litter carbon (v 1 in Eq. A3) and microbial biomass (v 2 in Eq. A3) represent the influences of those two carbon pools at any time on soil carbon. Because those elements are equal to 1 (see the matrix in Eq. A3), the oscillation of litter carbon and microbial biomass will also cause the response of soil carbon to be oscillatory, although the oscillation is small and dampens very quickly. In model B, the third elements of the eigenvectors corresponding to litter carbon and microbial biomass are zero (see the bottom row of the matrix in Eq. A4 of Wang et al., 2014) and therefore oscillatory responses of litter carbon and microbial biomass have no effect on the response of soil carbon. The eigenvalue of the soil carbon in model B is negative real; therefore the response of soil carbon to a small perturbation always is monotonic and stable in model B (see Appendix A in Wang et al., 2014).

Appendix B: Soil temperature at which equilibrium soil carbon pool is minimum (T x )
The steady-state soil carbon pool size of model A is The first term on the right-hand side of Eq. (B1) always decreases with an increase in T s , and the second term has two parts: 1 + a ε −1 − 1 and µ b K b µ s . Because both K b and µ s increase with T s exponentially, and the sensitivity µ s to T s is much greater than K b , K b µ s always decreases with an increase in T s , and that decrease is much greater than the increase in 1 + a ε −1 − 1 with T s . As a result, the second term also decreases with an increase in soil temperature, independent of temperature sensitivity of µ b . In summary for model A, dC * s dT s < 0.
The steady-state pool of soil carbon in model B is Assuming that V s µ b ε ε+a(1−ε) >>1, we can therefore approximate C * s as It can be easily shown that T x can only exist only when β k + b − β v ≤ 0 and 0 < a < 1 and when a = 0, T x does not exist and for model B.

Appendix C: Derivation of an analytic approximation for the timing and magnitude of the maximum microbial biomass after priming
Both models can be used to simulate the response of soil carbon to priming by specifying different initial pool sizes for the primed and control treatments. The initial values are C l (t = 0) = C * l + C l ; C b (t = 0) = C * b and C s (t = 0) = C * s for the priming treatment; C l (t = 0) = C * l ; C b (t = 0) = C * b and C s (t = 0) = C * s for the control. Here we assume that all pools are at equilibrium just before the priming treatment at t = 0. C * l , C * b , and C * s are equilibrium pool sizes, and C l is the amount of litter carbon added at time t = 0. No carbon is added to both treatments after t = 0.
The CO 2 efflux from soil carbon decomposition is calculated using Eq. (15) for model A and Eq. (16) for model B. Therefore we need to solve the three equations for C b and C s for t>0. Observations show that maximum priming response occurs soon after priming treatment (Kuzyakov et al., 2000); therefore, maximum priming response can be considered as a short-timescale phenomenon. At a short timescale, C s can be considered as being constant, and the maximum CO 2 efflux from the priming treatment will occur when the microbial biomass reaches a maximum after t = 0. Therefore we will use a second-order Taylor expansion to obtain the approximate solutions to the timing and magnitude of maximum CO 2 efflux from the soil carbon decomposition in the priming treatment for each model.
For model A, Eqs. (1) and (2) for both treatments after t>0 becomes dC l dt = −µ l C l As the litter pool size at time t = 0 is above its equilibrium value, the microbial biomass will likely increase after t = 0 and then reach its maximum value.
Equations (C1), (C2), and (C3) can be simplified using variable substitution. Let Then those three equations can be written as with the initial pool sizes of C b (0) = a 3 a 1 ε 1−ε , C s (0) = a 3 a 1 ε 1−ε + a + 1 + a 1−ε ε for both treatments, and C l (0) = (1 − a)( a 3 a 1 + 1−ε ε ) + C l for the primed treatment, and C l (0) = (1 − a)( a 3 a 1 + 1−ε ε ) for the control treatment. At relatively short timescales, a 2 1, C s (t) → C s (t = 0) Microbial biomass carbon after t = 0 can be approximated using the second-order Taylor expansion (Abramowitz and Stegun, 1972) Differentiating both sides of Eq. (C5) with respect to t, we have Assuming that C b is maximum at t = t max,A , then C b t max,A = 0. Equation (C6) becomes Both C b (0) and C b (0) can be obtained differentiating Eq. (C4) at t = 0, giving Substituting Eqs. (C8) and (C9) into (C7), and solving for t max,A , we have Substituting Eq. (C10) into (C5), we have the maximum microbial biomass at t max,A , or C bmax,A for the primed treatment as follows: The maximum rate of CO 2 release from decomposition of soil organic carbon, F CO 2 at t = t max,A is given by Similarly we derived the approximations for the timing (t max,B ) and magnitude of maximum microbial biomass (C b max, B ) in the primed treatment at t>0 as where M l = V l C * l + C l C * l + C l + K l .
The rate of CO 2 release from decomposition of soil carbon, F B , for model B at time t = t max,B is given by Comparison with numerical simulations shows that the relative error of Eq. (C12) is < 3 % across soil temperature and carbon input within their realistic ranges. However, errors in Eq. (C15) for model B can be quite large, particularly at high carbon input. Equation (C15) is only reasonably accurate (relatively error < 10 %) at low carbon input < 700 g C m −2 .