The Millennial model: in search of measurable pools and transformations for modeling soil carbon in the new century

Abstract

Soil organic carbon (SOC) can be defined by measurable chemical and physical pools, such as mineral-associated carbon, carbon physically entrapped in aggregates, dissolved carbon, and fragments of plant detritus. Yet, most soil models use conceptual rather than measurable SOC pools. What would the traditional pool-based soil model look like if it were built today, reflecting the latest understanding of biological, chemical, and physical transformations in soils? We propose a conceptual model—the Millennial model—that defines pools as measurable entities. First, we discuss relevant pool definitions conceptually and in terms of the measurements that can be used to quantify pool size, formation, and destabilization. Then, we develop a numerical model following the Millennial model conceptual framework to evaluate against the Century model, a widely-used standard for estimating SOC stocks across space and through time. The Millennial model predicts qualitatively similar changes in total SOC in response to single factor perturbations when compared to Century, but different responses to multiple factor perturbations. We review important conceptual and behavioral differences between the Millennial and Century modeling approaches, and the field and lab measurements needed to constrain parameter values. We propose the Millennial model as a simple but comprehensive framework to model SOC pools and guide measurements for further model development.

Appendix: Model description

The equations we have chosen below reflect one possible mathematical expression of the Millennial conceptual model, but there are many possible numerical models for different applications. For example, decomposition of POM is here represented by a double Monod relationship, limited by both POM and microbial biomass, but for an application where competition between chemical species is particularly important, for example, ECA kinetics could be used instead (Tang 2015). Similarly, we chose temperature and moisture scalars to minimize steady-state differences between Millennial and Century for the purpose of model comparison, but for dynamic predictions one could apply Arrhenius temperature sensitivity and one of several semi-mechanistic moisture functions (Davidson et al. 2012; Manzoni et al. 2014).

The system of equations below is modeled on the conceptual figure (Fig. 1), tracking the size of and transfers between five C pools: POM, LMWC, aggregate C, MAOM, and microbial biomass. The change in POM (P) stock with time is governed by the balance between plant C input and aggregate C breakdown, aggregate C formation, and decomposition,

$$dP/dt = p_{i} F_{i} + p_{a} F_{a} - F_{pa} - F_{pl} ,$$

(1)

where F _i is aboveground plant litter, root litter and root exudates, p _i is the proportion of C input allocated to POM (1/3 of inputs to POM and 2/3 of inputs to LMWC after Oleson et al. 2013), p _a is the proportion of C in aggregate breakdown allocated to POM, F _a is aggregate C breakdown, F _pa is aggregate carbon formation from POM, and F _pl is decomposition of POM into LMWC. Decomposition of POM is governed by a double Michaelis–Menten equation,

$$F_{pl} = V_{pl} \frac{P}{{K_{pl} + P}} \frac{B}{{K_{pe} + B}} S_{t} S_{w} ,$$

(2)

where V _pl is the maximum rate of POM decomposition, K _pl is the half-saturation constant, B is the microbial biomass carbon, and K _pe is the half-saturation constant of microbial control on POM mineralization. The terms S _t and S _w refer to the temperature and moisture scalar, respectively, and are taken from DAYCENT, the daily time-step version of the Century model (Parton et al. 1998), to minimize differences in temperature and moisture effects between the Century and Millennial models due to choice of scalar,

$$S_{t} = \left( {\frac{{t_{2} + ( t_{3} /pi) {\text{atan}}(pi (T - t_{1} ))}}{{t_{2} + ( t_{3} /pi) {\text{atan}}(pit_{4} (T_{ref} - t_{1} ))}}} \right),$$

(3)

$$S_{w} = \frac{1}{{1 + w_{1} \exp ( - w_{2} RWC)}},$$

(4)

where T is the current temperature, T _ref is the reference temperature, t ₁ is the x-axis location of the inflection point (°C), t ₂ is the y-axis location of the inflection point, t ₃ is the distance from the maximum point to the minimum point, and t ₄ is the slope of the line at the inflection point. For the water scalar, RWC is the relative water content calculated as the fraction of field capacity, and w ₁ and w ₂ are empirical parameters. The temperature scalar is an arctangent function that predicts a decline in temperature sensitivity with increasing temperature and the water scalar depends on RWC, where the maximum effect on biological activity occurs at field capacity (volumetric water content = 0.35, RWC = 1.0) (Parton et al. 2010).

The formation of aggregate C (A) from POM follows Michaelis–Menten dynamics,

$$F_{pa} = \frac{{V_{pa} P}}{{K_{pa} + P}}\left( {1 - \frac{A}{{A_{max} }}} \right) S_{t} S_{w} ,$$

(5)

where V _pa is the maximum rate of aggregate formation, K _pa is the half-saturation constant of aggregate formation, and A _max is the maximum capacity of C in soil aggregates. Soil aggregate C breakdown is partitioned to POM and MAOM,

$$F_{a} = k_{b} S_{t} S_{w} A,$$

(6)

where k _b is the rate of breakdown.

The change in LMWC (L) depends on LMWC input, the leaching rate, decomposition of POM, adsorption to minerals, and microbial uptake. In a multilayer version of the Millennial model LMWC would also depend on leaching input, but in this single layer version we assume that the leaching input is included in the LMWC input,

$$dL/dt = F_{i} \left( {1 - p_{i} } \right) - F_{l} + F_{pl } - F_{lm} - F_{lb} ,$$

(7)

where F _l is the LMWC leaching loss,

$$F_{l} = k_{l} S_{t} S_{w} L$$

(8)

and where k _l is the leaching rate, F _lm is the adsorption of LMWC to MAOM, and F _lb is the uptake of LMWC by microbial biomass. Adsorption of LMWC to minerals is controlled by a Langmuir saturation function,

$$F_{lm} = S_{t} S_{w} L\left( {\frac{{K_{lm} Q_{max} L}}{{1 + (K_{lm} L)}} - M} \right)/Q_{max} ,$$

(9)

$$K_{lm} = 10^{( - 0.186pH - 0.216)} ,$$

(10)

$$Q_{max} = BD10^{{(c_{1} \log (\% clay) + c_{2} )}} ,$$

(11)

where K _lm is the binding affinity that is adjustable based on the pH. Q _max is the maximum sorption capacity (mg C kg⁻¹ dry soil) that is converted to C density (g C m⁻²) by multiplying soil bulk density (BD = 1350 kg m⁻³), assuming a 1 m soil profile. The parameters c ₁ and c ₂ are the coefficients for computing Q _max from the clay content in percent, derived from Mayes et al. (2012). The Langmuir function parameters were derived from measurements of DOC sorption on over 200 soils in the eastern US. The measurements demonstrate a nonlinear saturation with respect to DOC concentrations in soils, and several recent models have used approaches that also impose a mechanism for DOC saturation on mineral surfaces (Wang et al. 2013; Riley et al. 2014; Ahrens et al. 2015; Dwivedi et al. 2017).

Microbial uptake of LMWC is a function of microbial biomass and LMWC concentration, temperature, water, and temperature-dependent CUE,

$$F_{lb} = V_{lm} S_{t} S_{w} L \frac{B}{{B + K_{lb} }}\left( {CUE_{ref} - CUE_{T} \left( {T - T_{{ae{ - }ref}} } \right)} \right),$$

(12)

$$F_{gr} = V_{lm} S_{t} S_{w} L \frac{B}{{B + K_{lb} }} \left( {1 - \left( {CUE_{ref} - CUE_{T} \left( {T - T_{{ae{ - }ref}} } \right)} \right)} \right),$$

(13)

where V_lm is the potential uptake rate of LMWC. F _gr is microbial growth-related respiration, K_lb is the half-saturation constant for microbial activity, CUE_ref is the reference CUE, and CUE_T is the CUE dependence on temperature. T_ae-ref and T are the reference and current temperature, respectively. Both MAOM and POM can be incorporated into the aggregate C pool,

$$dA/dt = F_{ma} + F_{pa} - F_{a} ,$$

(14)

$$F_{ma} = \frac{{V_{ma} M}}{{K_{ma} + M}}\left( {1 - \frac{A}{{A_{max} }}} \right)S_{t} S_{w} ,$$

(15)

where F _ma is the carbon flow from MAOM to aggregate C, V _ma is the maximum rate of aggregate formation, and K _ma is the half-saturation constant of aggregate formation. MAOM is formed by adsorption of LMWC and microbial necromass, and is affected by transfer into and out of the aggregate C pool,

$$dM/dt = F_{lm} + F_{bm} - F_{ma} + F_{a} \left( {1 - p_{a} } \right),$$

(16)

$$F_{bm} = k_{mm} S_{t} S_{w} B,$$

(17)

where F _bm is the carbon flow from microbial biomass to MAOM, namely adsorption of necromass, and k _mm is the adsorption rate of microbial biomass. In this particular iteration of the Millennial model, we assume that adsorbed microbial biomass is no longer alive, but by allowing adsorbed microbial biomass to take up LMWC and perform growth and maintenance, one could modify the model to accommodate the assumption that live microbial biomass can sorb to minerals, or even to other microbes (i.e., biofilms). Microbial biomass changes as a result of uptake, adsorption to minerals, and loss via maintenance,

$$dB/dt = F_{lb} - F_{bm} - F_{mr} ,$$

(18)

$$F_{mr} = k_{m} S_{t} S_{w} B,$$

(19)

where F _mr is the maintenance respiration of microbial biomass, and k _m is the microbial turnover rate (Table 3).

Table 3 Parameters, pools, fluxes and other variables used in the Millennial model