Stoichiometrically coupled carbon and nitrogen cycling in the MIcrobial-MIneral Carbon Stabilization model version 1.0 (MIMICS-CN v1.0)

Explicit consideration of microbial physiology in soil biogeochemical models that represent coupled carbon–nitrogen dynamics presents opportunities to deepen understanding of ecosystem responses to environmental change. The MIcrobial-MIneral Carbon Stabilization (MIMICS) model explicitly represents microbial physiology and physicochemical stabilization of soil carbon (C) on regional and global scales. Here we present a new version of MIMICS with coupled C and nitrogen (N) cycling through litter, microbial, and soil organic matter (SOM) pools. The model was parameterized and validated against C and N data from the Long-Term Inter-site Decomposition Experiment Team (LIDET; six litter types, 10 years of observations, and 13 sites across North America). The model simulates C and N losses from litterbags in the LIDET study with reasonable accuracy (C: R²=0.63; N: R²=0.29), which is comparable with simulations from the DAYCENT model that implicitly represents microbial activity (C: R²=0.67; N: R²=0.30). Subsequently, we evaluated equilibrium values of stocks (total soil C and N, microbial biomass C and N, inorganic N) and microbial process rates (soil heterotrophic respiration, N mineralization) simulated by MIMICS-CN across the 13 simulated LIDET sites against published observations from other continent-wide datasets. We found that MIMICS-CN produces equilibrium values in line with measured values, showing that the model generates plausible estimates of ecosystem soil biogeochemical dynamics across continental-scale gradients. MIMICS-CN provides a platform for coupling C and N projections in a microbially explicit model, but experiments still need to identify the physiological and stoichiometric characteristics of soil microbes, especially under environmental change scenarios.

1 Introduction

Soils contain the largest actively cycling terrestrial carbon (C) stocks on earth and also serve as the dominant source of nutrients, like nitrogen (N), that are critical for maintaining ecosystem productivity (Gruber and Galloway, 2008; Jobbágy and Jackson, 2000). Soil C cycle projections and their response to global change factors remain highly uncertain (Bradford et al., 2016; Todd-Brown et al., 2013), but recent empirical insights into microbial processing of soil C provide opportunities to update models and reduce this uncertainty (Cotrufo et al., 2013; Kallenbach et al., 2016; Lehmann and Kleber, 2015; Schmidt et al., 2011; Six et al., 2006). Several models have been developed recently with explicit representation of nonlinear microbial C processing dynamics, including the MIcrobial-MIneral Carbon Stabilization (MIMICS) model (Sulman et al., 2018; Wieder et al., 2014, 2015b) and others (Abramoff et al., 2017; Allison, 2014; Fatichi et al., 2019; Hararuk et al., 2015; Robertson et al., 2019; Sulman et al., 2014; Wang et al., 2013, 2014a, 2017). While these models serve different purposes, some can be as good as or better than models without explicit microbial pools at simulating global soil C stocks and the response of soil C to environmental perturbations (Wieder et al., 2013, 2015b), and they also predict very different long-term responses of soil C to global change (Wieder et al., 2013, 2018). Microbially explicit models have thus furthered our understanding of C cycling in the terrestrial system, but they also provide new opportunities to explore couplings between C and nutrient cycles, especially N.

Terrestrial models that couple C and N cycles reveal important ecosystem feedbacks that are absent from C-only models. For example, across ecosystems, experimental manipulations consistently indicate that N availability limits plant productivity (LeBauer and Treseder, 2008). C-only model configurations in models typically predict that CO₂ fertilization will result in a large increase in both plant productivity and the land C sink in coming decades, but nutrient limitation may constrain the magnitude of this terrestrial ecosystem C uptake (Wieder et al., 2015a; Zaehle et al., 2015; Zaehle and Dalmonech, 2011). As terrestrial models increasingly represent coupled C–N biogeochemistry, accurate model estimates of N release from soil organic matter (SOM) will become important to reducing uncertainty in the CO₂ fertilization response of the terrestrial C cycle.

Currently, most biogeochemical models that couple C and N cycles have an implicit representation of microbial activity. These conventional models represent SOM decomposition with the assumption that chemical recalcitrance of organic matter dictates the turnover of litter and SOM pools (Luo et al., 2016). Carbon and N fluxes represented in these models are directly proportional to donor pool sizes, without any explicit representation of the microbes that mediate these fluxes (Schimel, 2001, 2013). Linear decay constants and transfer coefficients determine the flow of C and N through a decomposition cascade, and rates of N immobilization and mineralization emerge from the interaction of fixed respiration fractions and the stoichiometry of donor and receiver SOM pools. The lack of plant–microbe–soil feedbacks in these models may limit their predictive capacity, especially in the face of environmental change. For example, in these models increased plant inputs to soil only build soil C and N stocks, and plants have no way to stimulate the microbial community to mine existing SOM for N without model modifications (Guenet et al., 2016; Wutzler and Reichstein, 2013). This “N mining” or “priming” effect, where increased plant inputs result in increased microbial activity and decomposition rates, has been demonstrated in experimental studies (Cheng and Kuzyakov, 2005; Dijkstra et al., 2013; Phillips et al., 2012) and may be a critical pathway for plants to obtain more N and support increased plant productivity under elevated CO₂ (Thomas et al., 2015; Zaehle et al., 2014).

Microbes are critical mediators of soil C–N couplings and the release of plant-available N. As such, models that explicitly consider microbial activity provide an opportunity to explore potential microbial control over soil C–N biogeochemical cycling and improve simulations of patterns in ecosystem C and N. Towards this end, multiple models have been introduced that explicitly consider the role of microbial activity in ecosystem C–N interactions (Averill and Waring, 2017; Fatichi et al., 2019; Huang et al., 2018; Schimel and Weintraub, 2003; Sistla et al., 2014; Sulman et al., 2014, 2017, 2018, 2019; Wang et al., 2014a, 2017, 2013). To date, the majority of these microbially explicit C–N models have been developed to explore soil biogeochemical interactions and microbial community dynamics, while only one has been validated for N dynamics across a continental-scale gradient (Fatichi et al., 2019).

Although there is great value in exploring diverse approaches to explicitly representing microbes in purely theoretical or site-specific applications, implementing these conceptual developments within larger-scale models requires convincing evidence that adding them improves model performance against large-scale data. Recent soil model comparisons report divergent responses to simulated global change experiments among microbially explicit model formulations, highlighting the large uncertainty in their underlying process-level representation and parameterization (Sulman et al., 2018; Wieder et al., 2018). The addition of explicit microbial pools may improve the predictive ability of landscape-scale models in the long run, but microbial models must be validated against landscape-scale datasets of a variety of pools and process rates before they can reasonably be expected to improve model performance and reduce uncertainty.

We developed a coupled C–N version of MIMICS (MIMICS-CN) to fill the need at the intersection of microbially explicit models, coupled C–N models, models that work well enough to be considered for use in ESMs, and models that can be validated against currently available large-scale data. The C-only iteration of MIMICS considers tradeoffs involved with microbial functional traits as well as both physicochemical (i.e., mineral associations) and chemical (i.e., recalcitrance) mechanisms of C stabilization in soil. Wieder et al. (2014, 2015b) and Sulman et al. (2018) evaluated this C-only version of MIMICS across site, continental, and global scales. Here we expand on this work, introducing MIMICS-CN, which incorporates stoichiometrically coupled C and N cycling of all microbial, litter, and SOM pools and stoichiometric constraints on microbial growth. Our core objectives were to (1) formulate a framework and parameterization for coupled C and N cycling in MIMICS; (2) validate MIMICS-CN against a continental-scale litter decomposition dataset (LIDET) and compare MIMICS-CN to a microbially-implicit, linear model (DAYCENT); and (3) evaluate equilibrium soil and microbial stocks and fluxes (and their parameter sensitivities) that are simulated by MIMICS-CN with data synthesized across published landscape-scale data. Our overarching goal was to create a microbially explicit coupled C–N model of soil that balances ecological realism with the practical considerations of large-scale simulation and to demonstrate the abilities of this model through parameterization, validation, and evaluation exercises using both dynamic and equilibrium data.

Figure 1 Overview of the pools and fluxes of C and N in MIMICS-CN. Black outlines indicate pools that contain C; green outlines indicate pools that contain N. Litter inputs (I) are determined based on site-specific net primary productivity and partitioned between metabolic and structural litter pools (LIT_m and LIT_s) using a site-specific litter quality metric (f_met) calculated using litter lignin and N content. Temperature-sensitive forward Michaelis–Menten kinetics (V_max and K_m; red lines) determine the flux of litter pool C and N and available SOM C and N (SOM_a) into microbial biomass (MIC_r and MIC_K). Fluxes of C into microbial pools result in respiration losses according to a defined carbon use efficiency (CUE). Microbes maintain biomass stoichiometry by spilling excess C as overflow respiration or excess N into the dissolved inorganic nitrogen pool (DIN) based on a prescribed biomass C : N. Microbial biomass turnover (τ, blue) varies by functional type (MIC_r and MIC_K) and is proportional to the square of microbial biomass. Microbial biomass turns over into available (SOM_a), physicochemically stabilized (SOM_p), and chemically stabilized (SOM_c) soil organic matter pools. Inorganic N (DIN) leaks from the model at a first-order rate. Numbers in parentheses indicate the equations in Appendix A that correspond to each depicted flux. Parameter values, units, and descriptions are given in Table 1.

2 Methods

2.1 Model formulation

MIMICS-CN builds upon the previous C-only version of MIMICS, described in Wieder et al. (2014, 2015b), using the same pool structure for N as for C plus an additional pool for dissolved inorganic nitrogen (DIN; Fig. 1). In-depth discussion of the reasoning behind the development of the C-only version of the model is available in these previous publications, but the general intent behind the development of MIMICS was to incorporate a simplified representation of the important aspects of microbial communities (biomass-dependent control of process rates, diversity in life history strategies, and physiological parameters) into a soil model that stabilizes organic matter through both physical (mineral-associated, protected from microbial decomposition) and chemical (recalcitrance-based, vulnerable to microbial decomposition) means. The C-only version of the model represents C flows through seven pools (Fig. 1): two litter pools, two microbial pools, and three SOM pools. Litter inputs to the model are partitioned into structural litter (LIT_s) and metabolic litter (LIT_m) pools based on estimates of litter quality for different biomes (Brovkin et al., 2012).

Temperature-sensitive forward Michaelis–Menten kinetics determine the flux of litter and SOM through microbial biomass pools that determine rates of organic matter decomposition, SOM formation, soil respiration, and nitrogen mineralization fluxes. The microbial functional groups are intended to broadly capture tradeoffs in microbial growth rates and growth efficiency, with rapidly growing microbial decomposers – low efficiency, r strategist (MIC_r) – and slower-growing microbial decomposers – higher efficiency, K strategist (MIC_K; Wieder et al., 2015b). In MIMICS-CN we extend these microbial physiological traits to include microbial stoichiometry and assume that the higher metabolic capacity of MIC_r also requires more nitrogen and, thus a lower microbial biomass C : N ratio. Fluxes of C into microbial pools result in respiration losses according to a defined carbon use efficiency (CUE) that varies by microbial functional group and substrate quality (e.g., structural or metabolic litter). Microbial pool sizes are moderated by inputs, CUE, and biomass-specific turnover rates. We implemented density-dependent microbial turnover (sensu Georgiou et al., 2017; see Appendix A) for this iteration of the model to make microbial pools behave realistically in response to small changes in C inputs (Wang et al., 2014b, 2016). The density-dependent turnover of microbial biomass dampens the oscillatory response of microbial biomass to perturbations.

Microbial biomass turns over into physicochemically stabilized (SOM_p), chemically stabilized (SOM_c), and a pool that is “available” for microbial decomposition (SOM_a). We consider the SOM_p pool to mostly consist of low C : N organic matter that is primarily composed of microbial products that are adsorbed onto mineral surfaces (e.g., mineral-associated organic matter, MAOM; Grandy and Neff, 2008). By contrast, the low-quality SOM_c pool consists of decomposed or partially decomposed litter that has more structural C compounds, such as lignin, and a higher C : N ratio (e.g., particulate organic matter, POM). Finally, SOM_a is the only SOM pool that is available for microbial decomposition; it contains a mixture of fresh microbial residues, products that are desorbed from the SOM_p pool (e.g., Jilling et al., 2018), as well as depolymerized organic matter from the SOM_c pool. We do not specifically consider soil aggregates, but we recognize that in some soils they are an important component of accruing and maintaining persistent organic matter.

The current representation of N cycling in MIMICS-CN is based on the threshold element ratio idea described in Sinsabaugh et al. (2009) and Mooshammer et al. (2014) whereby organisms maintain biomass stoichiometry by spilling excess C or N on either side of a threshold ratio. We modified the C-only iteration of MIMICS to include N by adding a parallel set of pools and fluxes for N, as well as a pool for inorganic N (Fig. 1). The C cycle drives decomposition with fluxes from litter and SOM pools to microbes based on biomass-C-based forward Michaelis–Menten kinetics. Parallel N fluxes are determined by the C : N ratio of the donor pools, which is a fixed parameter for the metabolic litter pool, varies with litter input chemistry for the structural litter pool, and depends on inputs for SOM pools. We use a fixed C : N of 15 for metabolic litter inputs, while the C : N of structural litter was allowed to vary to ensure conservation of total N inputs from litterfall (Table 1).

Table 1 Parameters used in MIMICS-CN for both LIDET and equilibrium simulations.

The coupling between C and N cycles in MIMICS-CN occurs in the microbial biomass: at each hourly time step, the total C and N in incoming fluxes available to microbes is summed and adjusted based on the C use efficiency (CUE; varies with microbial functional group and substrate) and N use efficiency (NUE; set to 0.85 for all fluxes entering microbial biomass pools in this model iteration). If the C : N of substrates being assimilated by microbial functional groups is greater or less than the C : N of the microbial biomass (defined as 6 and 10 for r and K strategists, respectively; Table 1), the microbes will spill excess C or N to maintain their biomass stoichiometry through overflow respiration or excess N mineralization. In MIMICS-CN the C : N ratio of SOM pools is flexible and determined by the inputs from microbial residues and direct inputs from litterfall fluxes (f_i; Fig. 1). All N fluxes into microbial pools leak a small quantity of N into a dissolved inorganic N pool (DIN) based on the model-defined NUE. At each time step, each microbial functional group can access a fraction of the inorganic N pool proportional to their fraction of total microbial biomass. Plant N uptake and ecosystem losses (both hydraulic and gaseous) of inorganic N are handled implicitly at this stage, with a fixed fraction (20 %) of DIN leaving the soil component model every time step.

2.2 Model parameterization and validation: cross-site litter decomposition

We parameterized and validated MIMICS-CN using C and N dynamics observed across multiple sites participating in the 10-year Long-Term Intersite Decomposition Experiment Team (LIDET) experiment (Adair et al., 2008; Harmon et al., 2009; Parton et al., 2007). The LIDET study selected standardized plant litter types with a range of litter quality (lignin and N concentration), placed litterbags containing 100 g of each litter type at sites across a continental-scale gradient of climatic conditions, and measured changes in the C and N in litterbags on an approximately annual basis for 10 years. Although the original dataset included 27 sites across North America, we utilized data from 14 sites ranging from Alaska to Puerto Rico based on the data available at those sites to drive MIMICS (see Wieder et al., 2015b, for site information). We focus our analysis on six leaf litters that were simulated across all sites that have been used previously to evaluate litter decomposition dynamics in terrestrial models (Bonan et al., 2013; Parton et al., 2007; Wieder et al., 2015b). Root litter types included in the original LIDET experiment were not included. The LIDET dataset is a robust appraisal of the impacts of climate and litter chemistry on litter decomposition and has been used as a dataset for comparing models of soil and litter decomposition in the past (Bonan et al., 2013). MIMICS has been used previously to simulate C losses in the LIDET study (Wieder et al., 2014, 2015b).

We parameterized MIMICS-CN using observations from Harvard Forest in Petersham, MA, USA. Observations included both litterbag C loss and N data from the LIDET study as well as measurements of soil C and N stocks and microbial C and N from other studies at Harvard Forest (Colman and Schimel, 2013). Multiple combinations of parameters produced equally good fits to litter decomposition data; thus ancillary data on soil and microbial C stocks were used to inform the parameter values presented here (Table 1). These ancillary data were not reported in LIDET and were not measured on identical plots to those used for the LIDET study (Harvard Forest encompasses multiple experiments and ecotypes), but these general targets were useful in distinguishing among model parameterizations. Our general targets for stocks at Harvard Forest included soil C and N (0–5 cm mineral soils, coniferous stand): 61 mg C cm⁻³ and 2.9 mg N cm⁻³; soil C : N of 21; and microbial biomass: 0.61 mg C cm⁻³ (estimated as 1 % of soil C based on Xu et al., 2013).

After parameterizing the model to match observations at Harvard Forest, the model was validated using data from the remaining LIDET sites. To represent litterbags in MIMICS-CN, we first spun up the underlying model to simulate steady-state soil C and N pools and fluxes across sites in the LIDET study using site-level measurements of mean annual temperature, clay content, and litter input quantity and litter chemistry (Wieder et al., 2015b). Then, we added a pulse of metabolic and structural litter based on the type of litter in the simulated litterbag. We tracked the C and N across all model pools for 10 years and calculated the C and N in litterbags as the difference between total model C and N in the simulations and total model C and N at steady state. In both the simulated and real litterbags, microbes immobilized N from the soil DIN pool, resulting in litterbag N contents for some time points in excess of the initial values. For each site, the model was sampled at time points equivalent to the real data collection dates in LIDET (approximately annually). Observed and modeled values of C and N in litterbags were compared by calculating R², root mean square error (RMSE), and bias.

To contextualize our results and better understand how our model functions compared to a widely used microbial-implicit model, we compared MIMICS-CN simulations of LIDET data against DAYCENT (Bonan et al., 2013) simulations of the same data. Bonan et al. (2013) used the full complement of 27 LIDET sites in their analysis, but here we subset those results for the 13 sites used in the MIMICS-CN validation. We calculated R², RMSE, and bias in the same way for each model and compared results across models, grouping results by biome.

2.3 Model evaluation: equilibrium C and N cycling

Building on the LIDET simulations, we independently synthesized observations to evaluate the patterns of C and N pools and fluxes across a variety of sites. Although direct, site-specific comparisons of modeled and observed values like microbial biomass would have been ideal, MIMICS-CN represents many variables that were not measured in the LIDET study and have not been synthesized across these Long-Term Ecological Research sites. Instead, we compared the range and distribution of pools (soil organic C and N, microbial biomass C and N, and total inorganic N) and fluxes (heterotrophic respiration and N mineralization) using the modeled LIDET simulations and published syntheses of observations from other sites (Cleveland and Liptzin, 2007; Colman and Schimel, 2013; Xu et al., 2013; Zak et al., 1994). To more directly compare measurements with model results, stock measurements were converted to units of percent of soil mass, and fluxes (heterotrophic respiration and net N mineralization rates) were converted to units of micrograms per cubic centimeter per hour (µg cm⁻³ h⁻¹). MIMICS reports pool values in units of grams per square centimeter (g cm⁻²) (0–30 cm); to compare MIMICS against observations we converted MIMICS values to percent by mass assuming a bulk density of 1.5 g cm⁻². Soil depth simulated by MIMICS (30 cm) is deeper than most of the observations in the compiled dataset, but the purpose of this exercise was to evaluate whether MIMICS produces realistic values for soil biogeochemical stocks and fluxes across continental-scale ecoclimatological and edaphic gradients, rather than making a direct site-specific comparison. The distribution of values produced by MIMICS across the LIDET sites was superimposed on the distributions of observed values to illustrate data–model agreement and to visualize the median and range of measurements across studies.

Finally, we documented relationships between model input variables (mean annual temperature, productivity, clay content, and litter quality) and the distribution of SOM pools that were simulated at the LIDET sites. Our aim with these analyses was to illustrate the underlying assumptions in the model and how they influence the size and distribution of C across SOM pools. Specifically, we wanted to explore how assumptions made in the model structure and parameterization of MIMICS determine the quantity and distribution of SOM pools and how they change among sites with variation in climatic, biological, and edaphic properties. To do this we looked at the absolute and relative contributions of each SOM pool simulated by MIMICS across the LIDET sites and conducted linear regressions to determine how environmental factors control their distributions. We also conducted linear regressions between soil C : N and both litter chemistry and environmental factors to assess the drivers of soil C : N in the model.

Figure 2 Litter decomposition time series simulated by MIMICS-CN (lines with shaded area) compared to observations (points and error bars) of (a) percent mass remaining and (b) percent of initial N remaining over 10 years for six different litter types at the Harvard Forest LTER. Litter decomposition data came from the LIDET study (Parton et al., 2007; Bonan et al., 2013; mean ±1 SD). Spread in the observations and model is largely generated by the effects of initial litter quality on decomposition rates and N dynamics. Model parameters were calibrated to fit MIMICS-CN to observations from Harvard Forest (Table 1).

Figure 3 MIMICS-CN simulations of percent (a, b) C remaining and (c, d) N remaining in litterbags in the LIDET study versus observed values, colored by (a, c) litter type or (b, d) biome. Dashed line shows the 1 : 1 line.

3 Results

3.1 Model parameterization and validation: cross-site litter decomposition

We parameterized MIMICS-CN to replicate litter C decay rates and N dynamics of six litter types observed in the LIDET study at the Harvard Forest Long-Term Ecological Research (LTER) site (Fig. 2). In its current parameterization, MIMICS slightly overestimates litter C loss at later stages of decay, but most time points are within uncertainty estimates of the observations (Fig. 2a). Similarly, for N, MIMICS-CN overestimates N accumulation in early stages of decay and underestimates N remaining at later stages, but most time points follow a reasonable trajectory given observations. MIMICS-CN also captures the effects of litter quality on both rates of litter decay (Fig. 2a) and litterbag N accumulation (Fig. 2b). The parameters we used to fit MIMICS-CN to Harvard Forest data also produce reasonable estimates of soil N stocks (2.0 vs. 2.9 mg N cm⁻³ for model and observations, respectively) and microbial biomass (0.65 vs. 0.61 mg C cm⁻³), although estimates of soil C (21 vs. 61 mg C cm⁻³) and soil C : N (11 vs. 21) are both lower than observations.

Parameter values used for this and subsequent simulations across all LIDET sites are shown in Table 1. Relative to the previous C-only version of the model (Wieder et al., 2014, 2015b), kinetic parameters and microbial turnover values were adjusted to account for density-dependent turnover (Georgiou et al., 2017). In addition, the fraction of structural litter that bypasses microbial biomass to enter the chemically protected pool (f_i) was increased from 5 % to 30 % as a means to produce reasonable values for total soil C : N. Finally, we adjusted the partitioning of microbial turnover to stable soil pools in order to more closely match distributions at Harvard Forest.

Applying this parameterization across all six litter types at 13 LIDET sites, MIMICS-CN simulates C losses and N dynamics from litterbags with an R² of 0.63 and 0.29, respectively (Fig. 3). MIMICS-CN captures effects of litter quality on decay rates, with faster rates of C loss and more rapid N mineralization simulated with more N-rich Drypetes glauca litter and slower rates of C loss and greater N immobilization simulated by low-quality Triticum aestivum litter (Fig. 3a, c). MIMICS-CN is best at capturing C loss rates in high- and intermediate-quality litters (Drypetes glauca, Pinus elliottii, Thuja plicata, and Acer saccharinum) but tends to underestimate litter C loss rates from the lowest-quality litter (T. aestivum). For N immobilization and loss, the model performs well especially for high-quality litters but underestimates N accumulation slightly in the lowest-quality litter. The model also captures broad climate effects on litter C loss, with slower decay rates in tundra and boreal forest sites and faster decay in tropical and deciduous forests (Fig. 3b).

Table 2 Goodness-of-fit statistics comparing MIMICS-CN and DAYCENT simulations to observations of C and N in decomposing litterbags in the LIDET study, aggregated by biome. DAYCENT results are subset from simulations in Bonan et al. (2013) to match the sites included in MIMICS-CN simulations. The values shown are the number of observations (n), Pearson's correlation coefficient squared (R²), root mean square error (RMSE), and bias calculated between observed and simulated percent C and N remaining. For more details on the sites grouped into each biome, see Wieder et al. (2015).

MIMICS-CN and DAYCENT simulations of LIDET decomposition data are compared in Table 2. Across a broad range of biomes, MIMICS-CN and DAYCENT both show good agreement with LIDET observations. Across sites MIMICS-CN has similar R² and RMSE values but lower bias compared to DAYCENT for mass loss (MIMICS-CN: R²=0.63, RMSE = 16.0, bias = −0.12; DAYCENT: R²=0.67, RMSE = 14.4, bias = 4.73) and percent N remaining (MIMICS-CN: R²=0.29, RMSE = 0.34, bias = 0.03; DAYCENT: R²=0.30, RMSE = 0.40, bias = 0.08). Broadly, MIMICS-CN outperformed DAYCENT in the warmest biomes, while DAYCENT excelled for colder sites for both C and N (Table 2), but the differences in model fit to data were slight and would be difficult to attribute to any particular differences in model structure. DAYCENT simulates decomposition based on initial litter chemistry and showed no site-specific effects on the maximum N immobilized or the relationship between C and N during decomposition for a given litter type (Figs. S1 and S2). By contrast, the amount of N that can be immobilized by a litterbag in MIMICS-CN is driven by the availability of N and the stocks and flows of N in the simulated steady-state soil, and MIMICS-CN showed site-specific variability in the shape of N immobilization and loss curves (Figs. 3 and 4).

Figure 4 MIMICS-CN simulations of immobilization–mineralization thresholds across litters of different quality. Litter quality (in terms of C : N and lignin content) decreases from upper-left panel to lower-right panel. Red triangles show model simulations of C losses vs. N losses from litterbags in the LIDET study. Colored dots show observed C vs. N losses across biomes (Parton et al., 2007).

Litter quality determines the timing of N immobilization vs. mineralization in observations. This produces a functional relationship between initial litter chemistry, C loss, and N immobilization or mineralization that is fairly consistent across sites (colored dots; Fig. 4). MIMICS-CN broadly captured litter quality effects on the timing and magnitude of N immobilization and mineralization dynamics across all biomes (red triangles; Fig. 4). For example, litters with high initial chemical quality consistently mineralize N throughout all stages of litter decay, and MIMIC-CN adequately captures this functional C–N relationship (Fig. 4a, b). By contrast, litters with lower initial chemical quality immobilize N during early stages of litter decay but subsequently mineralize N as decomposition proceeds. MIMICS-CN broadly captures these patterns but without as much variation as the observations (Fig. 4c–f). The lowest-quality litter (Triticum aestivum) immobilizes N until only 40 % of C remains in litterbags. Although MIMICS-CN potentially underestimates total N immobilization in Triticum aestivum litter, it does capture the point at which net N mineralization begins (Fig. 4f).

Table 3 Ranges of MIMICS-CN estimates of steady-state values for a variety of soil pools and fluxes, compared against observed ranges from several continent-wide data synthesis studies. The ranges of values included for MIMICS-CN are derived from simulations of sites included in the LIDET study.

^* Depths simulated by MIMICS-CN are for the top 30 cm of soil, whereas published ranges represent measurements ranging from the top 5 to top 30 cm.

3.2 Model evaluation: equilibrium C and N cycling

Across all sites and litter types in the LIDET simulations, the ranges of underlying pool sizes and process rates in MIMICS-CN were compared against published ranges from similarly diverse sets of sites (Cleveland and Liptzin, 2007; Colman and Schimel, 2013; Xu et al., 2013; Zak et al., 1994). MIMICS-CN simulations produced reasonable equilibrium values for most pools and fluxes (Table 3 and Fig. 5). In general, the range of values across the 13 sites simulated by MIMICS was smaller than the ranges across the thousands of sites included in the compiled dataset of observations. For example, total soil C ranged from 7.0 to 50 mg C cm⁻³ in MIMICS simulations but ranged from 2.7 to 610 mg C cm⁻³ in observations. Despite this discrepancy, the median values of the simulations and observations were generally within reason (Fig. 5). The distributions of measured and modeled values for microbial biomass C and N as a percent of total soil C and N overlapped, providing evidence that the model reasonably represents microbial stoichiometry, microbial activity as a function of biomass, and microbial biomass as a function of SOM. For soil C : N, the model tended to produce low values with a relatively narrow range, relative to observed values.

Finally, we explored the environmental controls on the distribution of SOM across physicochemically protected, chemically protected, and available pools in MIMICS-CN by examining the correlations between pool sizes and salient input variables (mean annual temperature, productivity, clay content, and litter lignin content). The results are shown in Fig. 6. The absolute concentration of SOM simulated across the LIDET sites was most strongly correlated with annual net primary productivity (ANPP; R²=0.52), but it also tended to increase with MAT, albeit inconsistently (Fig. 6a; R²=0.15). The distribution of SOM across stabilized pools strongly favored chemically protected SOM at sites with lower temperatures, while the relative proportion of physicochemically protected SOM increased with increasing temperature (Fig. 6b). The relative proportion of SOM in the available pool remained fairly consistent across simulated sites. Physicochemically protected SOM was tightly positively correlated with the product of ANPP and clay content (R²=0.96; Fig. 6c), while chemically protected and available SOM were negatively correlated with MAT (Fig. 6d; R²=0.40 and 0.47, respectively) and positively correlated with litter lignin content (Fig. 6e; R²=0.68 and 0.32, respectively). The C : N of individual pools was fairly consistent across sites and tended to be higher for chemically protected SOM (∼ 15) than available (∼ 8) or physicochemically protected SOM (∼ 10). As a result, soil C : N was largely driven across sites by the distribution of SOM across pools, especially the absolute size of the SOM_p pool (Fig. 6f; R²=0.79). Given that clay content was an important driver of physicochemically protected SOM in the model, clay content was tightly correlated with soil C : N (R²=0.88). Other litter characteristics and environmental factors were not strong drivers of soil C : N (R² for MAT: 0.42; litter lignin: 0.03; litter C : N: 0.005).

4 Discussion

Terrestrial models are increasingly representing coupled C–N biogeochemistry, and MIMICS-CN is among the first attempts to do so with a microbially explicit soil biogeochemical model that can be used to project C and N dynamics across continental-scale gradients. Our formulation and parameterization of MIMICS-CN captures site level observations of litter C loss and N immobilization at the Harvard Forest LTER site (Fig. 2). Cross-site validation of the model demonstrates that it broadly captures climate and litter quality effects on rates of C and N transformations from the LIDET observations (Figs. 3–4). Notably, the results simulated by MIMICS-CN represent N dynamics during litter decomposition about as well as a first-order model that implicitly represents microbial activity (Table 2). It also generates steady-state pools and fluxes of C and N that seem reasonable compared to published syntheses (Table 3; Fig. 5). Below we discuss these dynamic and equilibrium model simulations in greater detail, as well as some of the limitations of MIMICS-CN that will be addressed in future work.

4.1 Model parameterization and validation: cross-site litter decomposition

We first parameterized and validated MIMICS-CN using the cross-site litter decomposition study, LIDET. Previous LIDET simulations using MIMICS have successfully replicated observed C loss patterns, and adding coupled N cycling to MIMICS neither improved nor degraded simulations of LIDET litter C losses relative to the C-only model (Figs. 2–3; Wieder et al., 2015b report global RMSE for the C-only model is 14.6 vs. 16.0 in this study). Our results show higher-than-observed rates of litter C mass loss in deciduous and coniferous forest (Figs. 2a, 3b; Table 2). This suggests that the partitioning of plant detrital inputs into litter pools that are chemically defined works well for initial stages of litter decay but may not consider the changes in substrate chemistry or microbial community succession that occur in later stages of decomposition that slow rates of mass loss (Berg, 2000; Melillo et al., 1989). Models that implicitly represent microbial activity capture this phenomena by using a three pool structure (Adair et al., 2008), and future studies can consider how to more mechanistically understand interactions between initial litter quality, decomposer communities, climate, nutrient availability, and late-stage litter decay rates (e.g., Craine et al., 2007; Hobbie et al., 2012; Wickings et al., 2012) in models like MIMICS-CN. In MIMICS-CN, carbon and nitrogen move together through model pools, but model dynamics are primarily driven by C, with N dynamics following suit based on pool stoichiometry. The N dynamics do, however, constrain C cycling in the model if microbes are N limited, in which case microbes lose excess C through overflow respiration. At equilibrium, microbes in our MIMICS-CN simulations primarily obtained N through recycling of SOM pools with favorably low C : N ratios, with the result that modeled microbes were almost always C limited at equilibrium and rarely exhibited overflow respiration. Large pulses of low-quality litter can perturb this equilibrium and induce N limitation, but in the absence of losses of or plant competition for inorganic and dissolved organic N, C cycling in MIMICS proceeds in essentially the same way with or without accounting for N.

MIMICS-CN accurately captured the stoichiometric relationships between C and N during litter decomposition (Fig. 4). This stoichiometric relationship has been well defined in the past using theoretical microbial stoichiometry and CUE (Parton et al., 2007), but comparable soil models without explicit microbial physiology have tended to overpredict N accumulation in litterbags (Bonan et al., 2013). Moreover, models without microbially explicit physiology also show N immobilization mineralization dynamics that are completely determined by initial litter quality, whereas MIMICS simulations show greater site-level variation (Figs. 4, S2). In MIMICS-CN, stoichiometric relationships drive litterbags to accumulate soil N until they reach a threshold C : N, after which litterbags become net sources of N. This threshold, representing the balance between microbial N requirements and availability, is a function of changes in litter stoichiometry during decomposition, as well as of the stoichiometry of microbes and their nutrient use efficiencies. By explicitly considering these dynamics, MIMICS-CN has a similar or lower RMSE for N remaining in litter bags than a model that implicitly represents microbes, DAYCENT (Table 2).

MIMICS-CN and DAYCENT capture N dynamics during decomposition with similar overall degrees of fit but for different reasons. In DAYCENT, N immobilization and loss dynamics are driven by initial litter chemistry, and good model fit to data is achieved by capturing the average N immobilized for a given litter type regardless of biome and climate conditions (see Figs. S1 and S2). By contrast, litterbag N immobilization in MIMICS-CN is driven by the availability of N in the underlying modeled soil and by site-specific effects (e.g., climate, clay content) on the simulated stocks and fluxes of N. As a result, MIMICS-CN generates greater variation in the amount N immobilized for a given litter type across sites (Figs. 3 and 4). Site-specific variability in N immobilization patterns is also clearly visible in LIDET observations (colored dots; Fig. 4), but the introduction of site-specific variability in MIMICS-CN does not substantially improve model fit to data relative to DAYCENT. Spatial variability in ecosystem processes, like N mineralization rates, may be linked to factors like local-scale microbial community composition, soil moisture, or mineralogy (Graham et al., 2016; Smithwick et al., 2005; Soranno et al., 2019; Doetterl et al., 2015). While more work needs to be done to understand the factors controlling within and among site variation in soil C–N dynamics (Bradford et al., 2017), these results highlight that the explicit representation of microbial activity in MIMICS-CN may present opportunities to explore factors responsible for biogeochemical heterogeneity across scales.

Although MIMICS-CN broadly captures appropriate climate and litter quality effects on leaf litter decomposition patterns, the model underestimates N accumulation in the highest C : N ratio litter (Triticum aestivum; Fig. 4f). Microbes in MIMICS-CN recycle nitrogen from necromass and necromass-derived SOM, which might allow microbes to scavenge the N required to decompose high C : N litter without having to accumulate it from the inorganic soil pool. In a real litterbag, necromass might be lost through leaching, and microbial access to recycled biomass might be limited, and some microbial-derived compounds may require extensive depolymerization and proteolysis before the N is available for recycling (Schulten and Schnitzer, 1997), thus favoring N uptake from the soil pool. Alternatively, N inputs to real litterbags in the LIDET study may have come from atmospheric deposition or other unintended sources that MIMICS-CN does not address. Nonetheless, the high C : N ratio of Triticum aestivum is not typical of the majority of litter inputs across diverse biomes (Brovkin et al., 2012), which are well within the range that MIMICS-CN can simulate.

Figure 5 Distributions of MIMICS-CN estimates of steady-state values for a variety of soil pools and fluxes, compared against observed ranges from several continent-wide data synthesis studies. Black lines show the median value across all observations; red lines show median value of MIMICS-CN simulations.

4.2 Model evaluation: equilibrium C and N cycling

We conducted additional model evaluation by comparing model pools and fluxes at equilibrium to published observations. The parameter values used in the LIDET simulations produced reasonable estimates of equilibrium pools (soil organic C and N, microbial biomass C and N, and total inorganic N) and fluxes (heterotrophic respiration and N mineralization) (Table 3; Fig. 5). In combination with the LIDET results, these results indicate that MIMICS-CN can produce realistic simulations of both the short-term dynamic processes involved in litter decomposition and the soil-forming processes that produce equilibrium pools and fluxes over much longer timescales. In addition, MIMICS-CN simulates microbial stoichiometry, microbial growth and turnover, and microbially mediated decomposition, rather than using prescribed values as in models that lack explicit representation of microbes. This increases the power of MIMICS-CN to explore the microbial and biogeochemical processes underpinning model predictions.

Figure 6 Variation in steady-state SOM pools and environmental factors controlling their distribution in MIMICS-CN simulations across LIDET sites. Top panels show the (a) total C stocks in physicochemically protected, chemically protected, and available SOM pools (SOM_p, SOM_c, SOM_a pools, respectively) arranged by the site mean annual temperature (MAT) or the (b) relative fraction of each SOM pool arranged in the same way. Upper-right and bottom panels show the correlations between C in each SOM pool and environmental drivers including: (c) SOM_p vs. the product of annual net primary productivity (ANPP) and clay content, (d) SOM_c and SOM_a vs. MAT, and (e) SOM_c and SOM_a vs. lignin content of litter inputs at each site. Finally, (f) soil stoichiometry is largely determined by the fraction of total SOM pools that are considered physicochemically protected.

Continent-wide observations of soil pools and fluxes range over several orders of magnitude (Table 3), but MIMICS simulations agreed well with the median of those ranges. Observations tended to be spread over a much larger range of values than the MIMICS-CN simulations, but these simulations only included information from 13 sites, while the observations included thousands of locations. The median values of observed and simulated values were within a factor of 2.5 for all pools (Fig. 5). Differences in measurement depth or error in estimated bulk density values could account for some of the differences between measurements and simulations and for the spread across observed values. This is less of a concern for three of the variables used here (soil C : N, microbial biomass C as a percent of total soil C, and microbial biomass N as a percent of total soil N), which are ratios that are comparable across sites. Microbial biomass C as a percent of total soil C and microbial biomass N as a percent of total soil N were highly conserved across sites, relative to soil stocks or microbial C or N, and may be particularly useful metrics for evaluating microbially explicit soil biogeochemical models since the size of the microbial biomass pool directly controls rates of SOM turnover and formation in models like MIMICS-CN. For these ratios, MIMICS-CN reproduced distributions and median values that overlapped well with observations. In future work, direct comparisons of modeled and measured values for these ratios at specific sites may shed light on the limitations of the model and the origins of data–model disagreement. However, even the simple range comparisons included here provide evidence that the mechanistic representation of soil biogeochemistry in MIMICS-CN is ecologically realistic. Examinations of model realism like this are a crucial step in transitioning from theory and small-scale model tests to applications in ESMs or at larger scales where evaluation data are more sparse.

Besides representing appropriate soil biogeochemical stocks, fluxes simulated by the models also agree well with observations. Specifically, MIMICS-CN simulations of heterotrophic respiration and net N mineralization rates fell within observed bounds, although the variation in observations was much greater than the variation in simulated values. Our simulations calculated rates at equilibrium assuming constant temperature and other factors, while real rates of these processes are driven by seasonally and diurnally variable temperature, soil moisture, and other factors, so predictably, our simulations produced smaller-than-observed variability in rates. MIMICS-CN produced total soil C : N values that fall within observed ranges, although observations again show greater variation in soil C : N ratios and have maximum values that are much higher than the maximum C : N ratios simulated by MIMICS-CN. SOM pools in MIMICS-CN are mostly comprised of microbial necromass, in addition to a small proportion of litter that enters SOM pools directly without first passing through microbial biomass. Increasing this proportion in the model is one way to increase the C : N of SOM pools and the overall system at equilibrium. At some sites, litter may contribute more directly to SOM pools than microbial necromass (Jilling et al., 2018). For example, forests often have a higher proportion of total soil C in the light fraction, which is almost entirely made up of plant residues, compared to agroecosystems and many grasslands (Grandy and Robertson, 2007). For those sites with large, direct contributions of plant matter to SOM, increasing the fraction of litter that passes directly into SOM in MIMICS may be appropriate.

4.3 Exploring emergent SOM dynamics

The distribution of SOM across simulated pools in MIMICS-CN (Fig. 6) illustrates how model-defined assumptions about pool stabilization mechanisms drive potential responses to environmental variables. The wide variation in SOM pool distributions among contrasting environments in our simulations provides support for experimental efforts aimed at distinguishing between SOM pools to understand SOM responses to environmental changes and potential ecosystem feedbacks. For example, global change factors like warming can cause a range of different responses among SOM pools (Conant et al., 2008; Li et al., 2013; von Lützow and Kögel-Knabner, 2009; Plante et al., 2010). Experimental studies also show that increases in SOM resulting from increased inputs are not typically evenly distributed across different SOM pools (Lajtha et al., 2017; Stewart et al., 2009), which can influence feedbacks to productivity as well as the persistence of soil C gains in response to shifts in climate. Thus, while our broad-scale projections of how and why SOM differs among pools need to be evaluated with experiments and data synthesis across environments, they can provide a starting point for understanding SOM responses to global change factors across environments.

In MIMICS, the turnover of chemically protected and available SOM pools is based on temperature-sensitive Michaelis–Menten kinetics and litter chemistry (the latter controlling allocation of litter pools to the different microbial functional groups). This results in SOM_C pools (analogous to light fraction or POM pools) that are negatively correlated with MAT and positively correlated with litter lignin content (Fig. 6d, e). Turnover of the physicochemically protected SOM pool, on the other hand, occurs via first-order kinetics with a rate constant modified by clay content, and the equilibrium values of this pool are determined by inputs that largely come from microbial biomass and biomass turnover rates (Fig. 1). Therefore, the equilibrium values of SOM_p (analogous to heavy fraction or MAOM pools) were strongly positively correlated with the product of ANPP and clay content (Fig. 6c). This relationship broadly reflects the expected importance of total soil C inputs and their potential to be preserved after microbial processing by association with clays (Kleber et al., 2015). However, these two variables are also likely to covary with others, especially MAT, highlighting the difficulty of isolating individual mechanisms that regulate SOM.

Across the sites included in these simulations, chemically protected SOM formed a higher proportion of total SOM at lower MAT, while physicochemically protected SOM was favored at warmer sites (Fig. 6b). In global simulations with the carbon-only version of MIMICS, these assumptions result in MIMICS projecting longer soil C turnover times in soil C pools and larger soil C pools in the tropics than other models (Koven et al., 2017; Wieder et al., 2018) and a higher vulnerability of high-latitude soil C stocks (Wieder et al., 2015b, 2019). Evaluating the accuracy of our model assumptions and the resulting patterns in soil C and N cycling requires coupling process-level studies of the fate of decomposing litter (e.g., using isotope tracers) to broad-scale evaluation of SOM pool distributions across environmental gradients.

Soil C : N ratios simulated by MIMICS-CN across sites were highly correlated with soil clay content (R²=0.88), suggesting that, in the model, soil stoichiometry emerges from the relative contributions of SOM across physicochemically and chemically protected pools (Fig. 6). Although the spread of C : N values across the sites simulated by MIMICS-CN was small (Fig. 6f), C : N tended to decrease with increasing temperature, and simulated soil C : N was more correlated with site temperature (R²=0.42) than any of the litter characteristics used to drive the model, such as litter lignin (R²=0.03) or litter C : N (R²=0.005). This result directly contradicts a recent study using a first-order linear model which presumed that litter quality and soil quality at equilibrium were directly proportional (Menichetti et al., 2019). Although many soil biogeochemical models prescribe soil C : N ratios for individual pools, the stoichiometry of SOM in MIMICS-CN is an emergent property of the model.

The lack of correlation between simulated soil C : N and litter C : N in MIMICS-CN simulations suggests an intriguing followup question: in the field, is SOM stoichiometry correlated with litter quality, or is it better explained by climate, edaphic, and mineralogical gradients that impact soil microbial community composition, microbial activity, and mineral-mediated mechanisms of SOM persistence? Various regional studies provide limited support for the relationships generated by MIMICS-CN between soil C : N and MAT (Miller et al., 2004) or clay content (Hassink et al., 1993; Homann et al., 2007; Jenny, 1941), though a large-scale synthesis of measurements across all of these variables is still needed. Presently, MIMICS-CN assumes that microbial biomass stoichiometry largely controls the C : N ratios of stable SOM, with relatively minor contributions from litter quality. However, a small proportion of litter inputs become stabilized in MIMICS-CN without first passing through the stoichiometric filter of microbial biomass, and increasing this fraction in the model is a means to increase the C : N of simulated stable SOM. The strength of the mineral sink for microbial necromass in the model also impacts the relative balance of microbe- or plant-derived stable SOM, which in turn impacts modeled soil C : N. This result implies that in the field, C : N stoichiometry might be used as a means to differentiate the degree to which a given soil fraction is derived from direct plant inputs or microbial biomass, and mineralogical variables might be useful for explaining differences in fraction distributions across soils that impact C : N. Studies like Mikutta et al. (2019) illustrate the way that C : N can be used to assess the relative contributions of plant matter or microbial residues to stable SOM. Future work will use measured C : N of soils and soil fractions and isotopic insights into the plant or microbial origins of stable SOM to improve the parameterization of this aspect of the model and better understand the relationship between mechanisms of SOM stabilization and soil stoichiometry.

4.4 Limitations and future work

MIMICS-CN combines reasonable biogeochemical simulations with the option to explore underlying microbial processes, but limitations remain. For example, MIMICS only represents two microbial groups with different stoichiometric and physiological parameters, but real soils contain a much more diverse array of microbial functional groups with different responses to environmental conditions and different couplings between C and N cycles. CUE and NUE are critical microbial parameters in MIMICS-CN, but the relationships between CUE and microbial community composition (Maynard et al., 2017), microbial growth rate (Molenaar et al., 2009; Pfeiffer et al., 2001), temperature (Allison, 2014; Dijkstra et al., 2011; Frey et al., 2013; Steinweg et al., 2008), substrate quality (Blagodatskaya et al., 2014; Frey et al., 2013; Sinsabaugh et al., 2013), or any number of other aspects of microbial metabolism are complex, difficult to quantify, and challenging to represent at the scale of a whole soil community (Geyer et al., 2016). In its current configuration, MIMICS-CN also simplifies a number of ecosystem biogeochemical processes, and there are several important pathways of N cycling currently absent from the model. For example, MIMICS-CN does not currently represent free living biological N fixation, direct mycorrhizal exchanges for plant C for microbial N, dissolved organic C or N losses, denitrification/nitrification/other inorganic N transformation and loss pathways, plant uptake of N, or inorganic N leaching beyond a simple linear decay rate. Some of these shortcomings may be remedied by integrating MIMICS with a full ecosystem biogeochemical model that represents the greater complexity of the plant–soil continuum. MIMICS-CN provides a pathway to reconcile mechanistic explanations for phenomena like priming and plant–soil feedbacks with emergent patterns in terrestrial biogeochemistry across landscapes. MIMICS-CN and microbial models like it are a good first step towards representing the complex ecological factors that drive the coupling of soil C and N biogeochemistry, including the distribution of SOM among functionally relevant pools and SOM C : N ratios. Future work could compare model formulations that take different approaches to microbial community and stoichiometric parameters (e.g., flexible microbial parameters like C : N or CUE, additional microbial groups, partitioning microbial metabolism into a greater number of pathways) and refinement of mechanisms that confer SOM persistence. These efforts should also assess the ramifications of different choices for simulating existing data and predicting the long-term response of soil C and N cycles to global change. Our work demonstrates that MIMICS-CN can reproduce site and litter quality effects on litter decomposition C and N dynamics at a landscape scale, while also pointing to the importance of underlying, interacting microbial and biogeochemical factors in regulating SOM dynamics. Future work coupling MIMICS-CN to experiments and syntheses relating the distribution of SOM across pools to their underlying controls across gradients will improve our confidence in our ability to understand and project SOM dynamics.

Appendix A: Model equations

The structure and assumptions in the C-only version of MIMICS have been described previously (Wieder et al., 2014, 2015b), and the structure and assumptions in MIMIC-CN are described in Sect. 2.1 (Model formulation) of the methods section of this paper. The C fluxes (mg C cm⁻³ h⁻¹) from donor to receiver pools in MIMICS-CN, numbered in Fig. 1, are defined by the following:

$$\begin{array}{}\text{(A1)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}={\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}\times V_{max\left[r\mathrm{1}\right]}\times {\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}}{(K_{\mathrm{m}\left[\mathrm{r}\mathrm{1}\right]}+{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}})}},\text{(A2)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}={\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}\times V_{max\left[\text{r2}\right]}\times {\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}}{(K_{\mathrm{m}\left[\mathrm{r}\mathrm{2}\right]}+{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}})}},\text{(A3)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}={\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}\times V_{max\left[\text{r3}\right]}\times {\mathrm{SOM}}_{\mathrm{a},\mathrm{C}}}{(K_{\mathrm{m}\left[\mathrm{r}\mathrm{3}\right]}+{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}})}},\text{(A4)}& {\displaystyle }& {\displaystyle }{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}\mathit{\_}{\mathrm{SOM}}_{\text{C}}={\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}^{\mathit{\beta }}\times {\mathit{\tau }}_{\left[r\right],}\text{(A5)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}={\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}\times V_{max\left[K\mathrm{1}\right]}\times {\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}}{(K_{\mathrm{m}\left[\mathrm{K}\mathrm{1}\right]}+{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}})}},\text{(A6)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}={\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}\times V_{max\left[K\mathrm{2}\right]}\times {\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}}{(K_{\mathrm{m}\left[\mathrm{K}\mathrm{2}\right]}+{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}})}},\text{(A7)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}={\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}\times V_{max\left[K\mathrm{3}\right]}\times {\mathrm{SOM}}_{\mathrm{a},\mathrm{C}}}{(K_{\mathrm{m}\left[\mathrm{K}\mathrm{3}\right]}+{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}})}},\text{(A8)}& {\displaystyle }& {\displaystyle }{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}\mathit{\_}{\mathrm{SOM}}_{\text{C}}={\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}^{\mathit{\beta }}\times {\mathit{\tau }}_{\left[K\right]},\text{(A9)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{p},\mathrm{C}}\mathit{\_}{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}}={\mathrm{SOM}}_{\mathrm{p},\mathrm{C}}\times D,\end{array}$$

$$\begin{array}{}\text{(A10)}& \begin{array}{rl}{\mathrm{SOM}}_{\mathrm{c},\mathrm{C}}\mathit{\_}{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}}& =\left({\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}\times V_{max\left[\text{r2}\right]}\times {\mathrm{SOM}}_{\mathrm{c},\mathrm{C}}}{({\mathrm{KO}}_{\left[r\right]}\times K_{\mathrm{m}\left[\mathrm{r}\mathrm{2}\right]}+{\mathrm{SOM}}_{\mathrm{c},\mathrm{C}})}}\right)\\ & +\left({\displaystyle \frac{{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}\times V_{max\left[K\mathrm{2}\right]}\times {\mathrm{SOM}}_{\mathrm{c},\mathrm{C}}}{({\mathrm{KO}}_{\left[K\right]}\times K_{\mathrm{m}\left[\mathrm{K}\mathrm{2}\right]}+{\mathrm{SOM}}_{\mathrm{c},\mathrm{C}})}}\right),\end{array}\end{array}$$

where pools and parameters are described in Sect. 2.1 and Table 1, respectively. The N fluxes (mg N cm⁻³ h⁻¹) from donor to receiver pools in MIMICS-CN are calculated based on the C fluxes between pools and the C : N ratio of donor pools. These fluxes are numbered in Fig. 1 and defined by the following:

$$\begin{array}{}\text{(A11)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{m},\mathrm{N}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{r},\mathrm{N}}=\mathrm{A}\mathrm{1}\times {\mathrm{LIT}}_{\mathrm{m},\mathrm{N}}/{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}},\text{(A12)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{s},\mathrm{N}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{r},\mathrm{N}}=\mathrm{A}\mathrm{2}\times {\mathrm{LIT}}_{\mathrm{s},\mathrm{N}}/{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}},\text{(A13)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{a},\mathrm{N}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{r},\mathrm{N}}=\mathrm{A}\mathrm{3}\times {\mathrm{SOM}}_{\mathrm{a},\mathrm{N}}/{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}},\text{(A14)}& {\displaystyle }& {\displaystyle }{\mathrm{MIC}}_{\mathrm{r},\mathrm{N}}\mathit{\_}\mathrm{SOM},\mathrm{N}=\mathrm{A}\mathrm{4}\times {\mathrm{MIC}}_{\mathrm{r},\mathrm{N}}/{\mathrm{MIC}}_{\mathrm{r},\mathrm{C}},\text{(A15)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{m},\mathrm{N}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{K},\mathrm{N}}=\mathrm{A}\mathrm{5}\times {\mathrm{LIT}}_{\mathrm{m},\mathrm{N}}/{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}},\text{(A16)}& {\displaystyle }& {\displaystyle }{\mathrm{LIT}}_{\mathrm{s},\mathrm{N}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{K},\mathrm{N}}=\mathrm{A}\mathrm{6}\times {\mathrm{LIT}}_{\mathrm{s},\mathrm{N}}/{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}},\text{(A17)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{a},\mathrm{N}}\mathit{\_}{\mathrm{MIC}}_{\mathrm{K},\mathrm{N}}=\mathrm{A}\mathrm{7}\times {\mathrm{SOM}}_{\mathrm{a},\mathrm{N}}/{\mathrm{SOM}}_{\mathrm{a},\mathrm{C}},\text{(A18)}& {\displaystyle }& {\displaystyle }{\mathrm{MIC}}_{\mathrm{K},\mathrm{N}}\mathit{\_}\mathrm{SOM},\mathrm{N}=\mathrm{A}\mathrm{8}\times {\mathrm{MIC}}_{\mathrm{K},\mathrm{N}}/{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}},\text{(A19)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{p},\mathrm{N}}\mathit{\_}{\mathrm{SOM}}_{\mathrm{a},\mathrm{N}}=\mathrm{A}\mathrm{9}\times {\mathrm{SOM}}_{\mathrm{p},\mathrm{N}}/{\mathrm{SOM}}_{\mathrm{p},\mathrm{C}},\text{(A20)}& {\displaystyle }& {\displaystyle }{\mathrm{SOM}}_{\mathrm{c},\mathrm{N}}\mathit{\_}{\mathrm{SOM}}_{\mathrm{a},\mathrm{N}}=\mathrm{A}\mathrm{10}\times {\mathrm{SOM}}_{\mathrm{c},\mathrm{N}}/{\mathrm{SOM}}_{\mathrm{c},\mathrm{C}}.\end{array}$$

Each time step, the microbial pools in MIMICS-CN take up inorganic N from the DIN pool proportional to the biomass in each pool. Subsequently, the C : N ratio of all the inputs to each microbial pool is calculated, and the microbial pools spill either excess C or excess N to maintain a model-defined C : N ratio of microbial biomass. The algorithm that determines the release of excess C or N is determined using the following equations:

$$\begin{array}{}\text{(A21)}& {\displaystyle }& {\displaystyle }{\mathrm{DINup}}_{\mathrm{r}}={\displaystyle \frac{(\mathrm{1}-N_{\mathrm{leak}})\times \mathrm{DIN}\times {\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}}{({\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}+{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}})}},\text{(A22)}& {\displaystyle }& {\displaystyle }{\mathrm{DINup}}_{\mathrm{K}}={\displaystyle \frac{(\mathrm{1}-N_{\mathrm{leak}})\times \mathrm{DIN}\times {\mathrm{MIC}}_{\mathrm{K},\mathrm{C}}}{({\mathrm{MIC}}_{\mathrm{r},\mathrm{C}}+{\mathrm{MIC}}_{\mathrm{K},\mathrm{C}})}},\text{(A23)}& {\displaystyle }& {\displaystyle }{\mathrm{upMIC}}_{\mathrm{r},\mathrm{C}}={\mathrm{CUE}}_{\left[\mathrm{1}\right]}\times (\mathrm{A}\mathrm{1}+\mathrm{A}\mathrm{3})+{\mathrm{CUE}}_{\left[\mathrm{2}\right]}\times \left(\mathrm{A}\mathrm{2}\right),\text{(A24)}& {\displaystyle }& {\displaystyle }{\mathrm{upMIC}}_{\mathrm{r},\mathrm{N}}=\mathrm{NUE}\times (\mathrm{A}\mathrm{11}+\mathrm{A}\mathrm{13}+\mathrm{A}\mathrm{12})+\mathrm{A}\mathrm{21},\text{(A25)}& {\displaystyle }& {\displaystyle }{\mathrm{CNup}}_{\mathrm{r}}={\displaystyle \frac{\mathrm{A}\mathrm{23}}{\mathrm{A}\mathrm{24}}},\text{(A26)}& {\displaystyle }& {\displaystyle }{\mathrm{Overflow}}_{\mathrm{r}}=\mathrm{A}\mathrm{23}-(\mathrm{A}\mathrm{24}\times min({\mathrm{CN}}_{\mathrm{r}},\mathrm{A}\mathrm{25}\left)\right),\text{(A27)}& {\displaystyle }& {\displaystyle }{\mathrm{Nspill}}_{\mathrm{r}}=\mathrm{A}\mathrm{24}-\left({\displaystyle \frac{\mathrm{A}\mathrm{23}}{max({\mathrm{CN}}_{\mathrm{r}},\mathrm{A}\mathrm{25})}}\right),\text{(A28)}& {\displaystyle }& {\displaystyle }{\mathrm{upMIC}}_{\mathrm{K},\mathrm{C}}={\mathrm{CUE}}_{\left[\mathrm{3}\right]}\times (\mathrm{A}\mathrm{5}+\mathrm{A}\mathrm{7})+{\mathrm{CUE}}_{\left[\mathrm{4}\right]}\times \left(\mathrm{A}\mathrm{6}\right),\text{(A29)}& {\displaystyle }& {\displaystyle }{\mathrm{upMIC}}_{\mathrm{K},\mathrm{N}}=\mathrm{NUE}\times (\mathrm{A}\mathrm{15}+\mathrm{A}\mathrm{17}+\mathrm{A}\mathrm{16})+\mathrm{A}\mathrm{22},\text{(A30)}& {\displaystyle }& {\displaystyle }{\mathrm{CNup}}_{\mathrm{K}}={\displaystyle \frac{\mathrm{A}\mathrm{28}}{\mathrm{A}\mathrm{29}}},\text{(A31)}& {\displaystyle }& {\displaystyle }{\mathrm{Overflow}}_{\mathrm{K}}=\mathrm{A}\mathrm{28}-(\mathrm{A}\mathrm{29}\times min({\mathrm{CN}}_{\mathrm{K}},\mathrm{A}\mathrm{30}\left)\right),\text{(A32)}& {\displaystyle }& {\displaystyle }{\mathrm{Nspill}}_{\mathrm{K}}=\mathrm{A}\mathrm{29}-\left({\displaystyle \frac{\mathrm{A}\mathrm{28}}{max({\mathrm{CN}}_{\mathrm{K}},\mathrm{A}\mathrm{30})}}\right).\end{array}$$

Inorganic N leaches slowly from the model according to a model-defined rate:

$$\begin{array}{}\text{(A33)}& {\displaystyle }{\displaystyle }\mathrm{LeachingLoss}=N_{\mathrm{leak}}\times \mathrm{DIN}.\end{array}$$

Given the fluxes defined above, the changes in C and N pools in each hourly time step (mg C or N cm⁻³) are described by the following:

$$\begin{array}{}\text{(A34)}& {\displaystyle \frac{{\mathrm{dLIT}}_{\mathrm{m},\mathrm{C}}}{\mathrm{d}t}}={\mathrm{I}}_{{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}}\times \left(\mathrm{1}-{\mathrm{f}}_{\mathrm{i},\mathrm{met}}\right)-\mathrm{A}\mathrm{1}-\mathrm{A}\mathrm{5},\text{(A35)}& {\displaystyle \frac{{\mathrm{dLIT}}_{\mathrm{s},\mathrm{C}}}{\mathrm{d}t}}={\mathrm{I}}_{{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}}\times \left(\mathrm{1}-{\mathrm{f}}_{\mathrm{i},\mathrm{struc}}\right)-\mathrm{A}\mathrm{2}-\mathrm{A}\mathrm{6},\text{(A36)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dMIC}}_{\mathrm{r},\mathrm{C}}}{\mathrm{d}t}}& =\mathrm{CUE}\left[\mathrm{1}\right]\times (\mathrm{A}\mathrm{1}+\mathrm{A}\mathrm{3})+\mathrm{CUE}\left[\mathrm{2}\right]\\ & \times \left(\mathrm{A}\mathrm{2}\right)-\mathrm{A}\mathrm{4}-{\mathrm{Overflow}}_{\mathrm{r}},\end{array}\text{(A37)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dMIC}}_{\mathrm{K},\mathrm{C}}}{\mathrm{d}t}}=& \mathrm{CUE}\left[\mathrm{3}\right]\times (\mathrm{A}\mathrm{5}+\mathrm{A}\mathrm{7})+\mathrm{CUE}\left[\mathrm{4}\right]\\ & \times \left(\mathrm{A}\mathrm{6}\right)-\mathrm{A}\mathrm{8}-{\mathrm{Overflow}}_{\mathrm{K}},\end{array}\text{(A38)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dSOM}}_{\mathrm{p},\mathrm{C}}}{\mathrm{dt}}}& ={\mathrm{I}}_{{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}}\times {\mathrm{f}}_{\mathrm{i},\mathrm{met}}+({\mathrm{f}}_{\mathrm{p},\mathrm{r}}\times \mathrm{A}\mathrm{4})\\ & +({\mathrm{f}}_{\mathrm{p},\mathrm{K}}\times \mathrm{A}\mathrm{8})-\mathrm{A}\mathrm{9},\end{array}\text{(A39)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dSOM}}_{\mathrm{c},\mathrm{C}}}{\mathrm{d}t}}& =I_{{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}}\times {\mathrm{f}}_{\mathrm{i},\mathrm{struc}}+({\mathrm{f}}_{\mathrm{c},\mathrm{r}}\times \mathrm{A}\mathrm{4})\\ & +(f_{\mathrm{c},\mathrm{K}}\times \mathrm{A}\mathrm{8})-\mathrm{A}\mathrm{10},\end{array}\text{(A40)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dSOM}}_{\mathrm{a},\mathrm{C}}}{\mathrm{d}t}}& =({\mathrm{f}}_{\mathrm{a},\mathrm{r}}\times \mathrm{A}\mathrm{4})+({\mathrm{f}}_{\mathrm{a},\mathrm{K}}\times \mathrm{A}\mathrm{8})\\ & +\mathrm{A}\mathrm{9}+\mathrm{A}\mathrm{10}-\mathrm{A}\mathrm{3}-\mathrm{A}\mathrm{7},\end{array}\text{(A41)}& {\displaystyle \frac{\mathrm{dLI}{T}_{\mathrm{m},\mathrm{N}}}{\mathrm{d}t}}={\displaystyle \frac{I_{{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}}\times \left(\mathrm{1}-f_{\mathrm{i},\mathrm{met}}\right)}{{\mathrm{CN}}_{\mathrm{m}}}}-\mathrm{A}\mathrm{11}-\mathrm{A}\mathrm{15},\text{(A42)}& {\displaystyle \frac{\mathrm{dLI}{T}_{\mathrm{s},\mathrm{N}}}{\mathrm{d}t}}={\displaystyle \frac{I_{{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}}\times \left(\mathrm{1}-f_{\mathrm{i},\mathrm{struc}}\right)}{{\mathrm{CN}}_{\mathrm{s}}}}-\mathrm{A}\mathrm{12}-\mathrm{A}\mathrm{16},\text{(A43)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dMIC}}_{\mathrm{r},\mathrm{N}}}{\mathrm{d}t}}& =\mathrm{NUE}\times (\mathrm{A}\mathrm{11}+\mathrm{A}\mathrm{13}+\mathrm{A}\mathrm{12})-\mathrm{A}\mathrm{14}\\ & +{\mathrm{DINup}}_{\mathrm{r}}-{\mathrm{Nspill}}_{\mathrm{r}},\end{array}\text{(A44)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dMIC}}_{\mathrm{K},\mathrm{N}}}{\mathrm{d}t}}& =\mathrm{NUE}\times (\mathrm{A}\mathrm{15}+\mathrm{A}\mathrm{17}+\mathrm{A}\mathrm{16})\\ & -\mathrm{A}\mathrm{18}+{\mathrm{DINup}}_{\mathrm{K}}-{\mathrm{Nspill}}_{\mathrm{K}},\end{array}\text{(A45)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dSOM}}_{\mathrm{p},\mathrm{N}}}{\mathrm{d}t}}& ={\displaystyle \frac{I_{{\mathrm{LIT}}_{\mathrm{m},\mathrm{C}}}\times \left(f_{\mathrm{i},\mathrm{met}}\right)}{{\mathrm{CN}}_{\mathrm{m}}}}\\ & +(f_{\mathrm{p},\mathrm{r}}\times \mathrm{A}\mathrm{14})+(f_{\mathrm{p},\mathrm{K}}\times \mathrm{A}\mathrm{18})-\mathrm{A}\mathrm{19},\end{array}\text{(A46)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dSOM}}_{\mathrm{c},\mathrm{N}}}{\mathrm{d}t}}& ={\displaystyle \frac{I_{{\mathrm{LIT}}_{\mathrm{s},\mathrm{C}}}\times \left(f_{\mathrm{i},\mathrm{struc}}\right)}{{\mathrm{CN}}_{\mathrm{s}}}}\\ & +(f_{\mathrm{c},\mathrm{r}}\times \mathrm{A}\mathrm{14})+(f_{\mathrm{c},\mathrm{K}}\times \mathrm{A}\mathrm{18})-\mathrm{A}\mathrm{20},\end{array}\text{(A47)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{{\mathrm{dSOM}}_{\mathrm{a},\mathrm{N}}}{\mathrm{d}t}}& =(f_{\mathrm{a},\mathrm{r}}\times \mathrm{A}\mathrm{14})+(f_{\mathrm{a},\mathrm{K}}\times \mathrm{A}\mathrm{18})\\ & +\mathrm{A}\mathrm{19}+\mathrm{A}\mathrm{20}-\mathrm{A}\mathrm{13}-\mathrm{A}\mathrm{17},\end{array}\text{(A48)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{\mathrm{dDIN}}{\mathrm{d}t}}& =\left(\mathrm{1}-\mathrm{NUE}\right)\times \left(\mathrm{A}\mathrm{11}+\mathrm{A}\mathrm{12}+\mathrm{A}\mathrm{13}+\mathrm{A}\mathrm{15}\right.\\ & \left.+\mathrm{A}\mathrm{16}+\mathrm{A}\mathrm{17}\right)+{\mathrm{Nspill}}_{\mathrm{r}}+{\mathrm{Nspill}}_{\mathrm{K}}-{\mathrm{DINup}}_{\mathrm{r}}\\ & -{\mathrm{DINup}}_{\mathrm{K}}-\mathrm{LeachingLoss}.\end{array}\end{array}$$