Soil microbial respiration from observations and Earth System Models

The global soil respiration database version 2.0 (SRDB-V2) [23, 24] is based on a comprehensive literature survey of field (not laboratory) measurements with a variety of techniques for Rh partitioning from Rs (table 1), including 4387 station records over global land from 1961 to 2009. For the common period of 1961–2005 with ESMs, 2569 and 435 station records of Rs and Rh are available for this study. Besides the Rh and Rs data, we also use other variables (e.g., temperature, precipitation, net primary production, and SOC) in SRDB-V2. Each value (particularly SOC) in SRDB-V2 usually represents the average of many values measured or generated by different methods from individual study.

Despite the inherent uncertainties in measurements, data entry, and sampling, this comprehensive database is still the best available and valuable for evaluating Rh from ESMs as long as these random and bias errors remain much smaller than the errors in the models. Because the sites of SRDB-V2 are dominated by temperate biomes (1689 sites) and global forests (1823 sites) (figure S1 available at stacks.iop.org/ERL/8/034034/mmedia), here we focus on five ecosystem types (boreal forest, temperate forest, tropical forest, temperate grassland, temperate agriculture) with more than 36 records of Rh for each type. Due to the different precipitation seasonality, 'temperate forest' here does not include Mediterranean woodland.

The ESMs simulate comprehensive physical, chemical, biological, and anthropogenic processes in the Earth system, including atmosphere, ocean, land surface, sea ice, glacier, and their interactions. Based on the availability of relevant variables at the time of our original analysis, output from eight ESMs (CanESM2, GFDL-ESM2M, HadGEM2-ES, MIROC-ESM, MPI-ESM-LR, INMCM4, CCSM4, and NorESM1-M) from the Coupled Model Intercomparison Project Phase 5 (CMIP5) [25] are explored here. The last two ESMs adopt the same Community Land Model version 4 with the coupled carbon–nitrogen cycle (CLM4-CN), while the other six include the carbon cycle only. For the historical simulations (1861–2005), the CO₂ concentration and land cover/land use change are imposed based on historical observations and estimates. For the 21st Century simulations, the CO₂ concentrations are prescribed under different representative concentration pathways (RCPs) and RCPn (with n = 2.6,4.5,6.0, and 8.5 in figure 1) is used to denote the emission scenario with total radiative forcing stabilized at n W m⁻² roughly before 2100. Taylor et al [25] discussed the setup of these simulations and participating ESMs, while the information on the land component of the eight ESMs is summarized in table S1 (available at stacks.iop.org/ERL/8/034034/mmedia). As widely used in ESM evaluations against point measurements, model results are bilinearly interpolated in this study to each observational site for the same year for model-data comparison. The mismatch between point measurements and model grid box averages will be further discussed in sections 3.1 and 4.1.

Zoom In Zoom Out Reset image size

The Rh parameterization as a function of temperature, soil moisture, and other processes is similar across different ESMs; it differs usually in the exact functional forms. The temperature sensitivity in all ESMs except GFDL-ESM2M is represented by the Q₁₀ or Arrhenius equations, both of which are functionally similar. In GFDL-ESM2M, Rh increases with temperature up to an optimal temperature and then decreases. Rh in all of these ESMs either increases monotonically with increasing soil water availability or increases up to an optimal moisture level and then decreases. CCSM4 and NorESM1 include the nitrogen limitation in the Rh computation. The number of litter and soil pools and the temperature and moisture functions are summarized in the table 1 in Todd-Brown et al [26]. For instance, CLM4-CN as used in CCSM4 and NorESM1-M computes Rh through three litter pools, three soil organic matter (SOM) pools, and a coarse woody debris pool as a converging cascade [27]. At each timestep, Rh (g C m⁻² s⁻¹) from each pool is calculated as [27, 28]

$$\begin{eqnarray} &&\mathrm{Rh}\propto k C \hspace{0.167em} {f}_{\mathrm{r}}\hspace{0.167em} {f}_{\mathrm{N}}\hspace{0.167em} {f}_{\mathrm{T}}\hspace{0.167em} {f}_{\mathrm{w}} \end{eqnarray} \tag{ 1 }$$

where k (h⁻¹) is the base pool-specific fixed decomposition rate, C (g C m⁻²) is the carbon content in the corresponding pool, f_r is a fixed fraction of the decomposition carbon flux released as CO₂ for specific plant functional type and tissue, f_T is based on Q₁₀ of 1.5 averaged over the top five soil layers, and f_w is proportional to soil moisture and the weighting depending on soil water potential (averaged over the top five soil layers). The nitrogen limitation f_N is computed from the competition between the plant and nitrogen immobilization demand (g N m⁻² s⁻¹) for the available soil mineral nitrogen pool. When the available nitrogen pool is insufficient to meet the combined demands, f_N is less than 1; otherwise f_N is set to 1. For T <− 10 °C or soil water potential <−10 MPa, Rh is taken as zero, consistent with the laboratory measurements [19].

The range of the global Rh estimates from eight ESMs is 41.3–71.6 Pg C a⁻¹ (figure 1), and their average is very close to the empirically inferred values centered around 60 Pg C a⁻¹ [29, 30]. However the two ESMs including coupled carbon–nitrogen dynamics (CCSM4 and NorESM1-M) underestimate the global Rh by 30% compared with these values (figure 1). Compared with the historical period from 1976 to 2005, the global annual Rh from 2071 to 2100 under the RCP2.6, RCP4.5, RCP6.0, and RCP8.5 scenarios is projected to increase by 10.6 ± 6.2%, 17.9 ± 8.6%, 21.8 ± 13.9%, and 31.0 ± 14.9%, respectively, across the eight ESMs, largely due to the increase of global mean surface temperature [31, 32] along with other confounding effects (e.g. increase in vegetation production). In absolute terms this represents widely diverging projected carbon cycle feedbacks in response to environmental change; this was greatest for the RCP8.5 condition (46.0–111.5 Pg C a⁻¹) but the variance is reduced if the two models including a nitrogen cycle are excluded (71.4–111.5 Pg C a⁻¹). The tropical latitudes represent the largest Rh source followed by the boreal latitudes in both historical and 21st century simulations in these eight models. The change of global annual Rh in response to temperature change confounded with CO₂ evolution and land use change from 1850 to 2005 varies from −0.56 to +12.6 Pg C °C⁻¹ across the ESMs. In contrast, this rate is just 2.3 Pg C °C⁻¹ in response to temperature change only (18; based on the CASA model with globally variant Q₁₀ values).

Because of the mismatch between point measurements and model grid box output, it is more appropriate to compare the observed and model-simulated ranges of Rh than the direct point-to-point comparisons. There is a substantial overlap between the magnitude of Rh estimated from observations and ESMs across all biomes (figure 2). In temperate forests, agriculture, and grasslands, the median and range of Rh predicted by the two ESMs incorporating the nitrogen cycle (CCSM4 and NorESM1-M) are closest to observed values (figures 2(b), (d) and (e)). In tropical forests, there is little difference between models with and without the nitrogen cycle (figure 2(c)), as the nitrogen supply is sufficient (f_N = 1 in equation (1)). Over the boreal forest, however, the constraint of nitrogen results in a significantly lower Rh from CCSM4 and NorESM1-M with the maximum Rh less than 200 g C m⁻² (lower than the median rate from observations), while the other six ESMs perform better on average compared with observations (figure 2(a)). Some models without the nitrogen constraint generate Rh values much higher than any observations. For instance, the minimum and maximum Rh values from the six carbon-only ESMs are 135 and 789 g C m⁻² over the boreal forest (figure 2(a)), in contrast to 25 and 591 g C m⁻² in observations. This deficiency is particularly apparent for temperate ecosystems and tropical forests; the minimum Rh from MPI-ESM-LR is as high as 379 g C m⁻², far above the observed minimum at 11 g C m⁻² (figures 2(b)–(e)).

Zoom In Zoom Out Reset image size

Rh in ESMs at each time step (usually 10–20 min) is regulated exponentially by temperature and linearly by soil moisture, and in CCSM4 and NorESM1-M limited by the availability of nitrogen (equation (1)). Averaging to longer time-steps (e.g., the annual averages) diffuses these relationships and linkages to other environmental drivers emerge. If the relationship between annual variables emerges from the mathematical structure (equation (1)) of ESMs but are absent in reality, then the mathematical structure is deficient. For example, surface temperature explains between 5% for CanESM2 and 43% for INMCM4 of the variance in annual Rh, while the corresponding value from observations is 16% (table 2). NPP explains between 63% and 97% of the variance in annual Rh, over 90% in five of the eight models, in contrast to the observational value of 23% (table S2 available at stacks.iop.org/ERL/8/034034/mmedia). Annual precipitation is a key driver of Rh in ESMs (0.19 ≤ r² ≤ 0.76) but much less important in observations (r² = 0.05, table S2). Similarly, the importance of SOC in the Rh computation is overestimated in ESMs compared with observations (table S2).

As an example, these relationships are illustrated for CCSM4 and observations in figure 3. The regression slope is much greater in CCSM4 than observed, indicating the greater sensitivity of simulated Rh to Pr, NPP, and SOC. The r² value in CCSM4 is also much greater than observed. Furthermore, the r² values between Rh and other plant-related variables (e.g., leaf area index (LAI), total vegetation carbon, standing litter carbon) are much greater in CCSM4 than observed. These results are expected, since there is a close causal relationship between LAI and NPP, and vegetation is the only source of the decomposition pools in ESMs, leading to higher correlations of Rh with NPP and SOC in ESMs. In contrast, available observational data from 57 records in SRDB-V2 indicate that only 42% of total Rh is from litter and the rest originates from the SOC and mineral carbon.

Zoom In Zoom Out Reset image size

Multicollinearity in the environmental driver is investigated here using partial least squares (PLS) regression (see supplementary text S1 available at stacks.iop.org/ERL/8/034034/mmedia) of Rh against annual mean surface temperature (T, ° C) and the three variable in figure 2: annual precipitation (Pr, mm), NPP (g C m⁻²), and total SOC (g C m⁻²). Pooling observations over all sites (figure S1 available at stacks.iop.org/ERL/8/034034/mmedia), the apparent Q₁₀ of Rh is calculated as 1.48 by computing the fold change in Rh with respect to temperature variation across all sites. The corresponding values vary from 1.13 to 1.89 (median  =1.30) for the eight ESMs (table 2). When other environmental drivers of Rh are factored out, the apparent Q₁₀ is increased to 1.83 in observations, but its change is ESM-dependent (i.e. it is increased in three ESMs but decreased in the remaining five) (table 1). Conclusions are similar if the other two regression methods are used (table S3 available at stacks.iop.org/ERL/8/034034/mmedia). Even though CCSM4 and NorESM1-M adopt the same land model (CLM4-CN), their Q₁₀ values differ from each other (table 2), and interestingly, they differ from the embedded Q₁₀ value of 1.5 for CLM4-CN at each time step in equation (1).

When all four confounded variables T,ln(Pr),ln(NPP), and ln(SOC) are used for the regression of ln(Rh),r² between Rh and fitted values increases substantially in both observations and ESMs (table 2). While this regression explains 55% of the variance of Rh in observations, it accounts for 75–98% of the Rh variance in ESMs. The root mean square error (RMSE) is also significantly reduced as the covariate number is increased from one to four, but this decrease is more dramatic in ESMs than in observations (table 2). Similar results are obtained using the other two regression methods (table S3).

Because observational studies report Rs much more frequently than Rh, we have developed two minimum requirements or 'reality checks' for ESM Rh estimates based on Rs reported in the SRDB-V2 database. While all eight ESMs provide Rh and total respiration (Rh plus plant and root Ra), none of these ESMs provides Rs as a standard output. Since total respiration includes Rs and plant Ra, the first requirement is that total respiration from ESMs should not be less than observed Rs. Similarly, since Rs includes Rh and root Ra, the second requirement is that Rh from ESMs should not be greater than observed Rs. The percentage of sites where the first requirement is met varies from 40% to 90% (figure 4(a)), while it is higher and varies from 41% to 98% for the second requirement (figure 4(b)). The first requirement is met at most sites over the boreal forest in the six carbon-only ESMs, but it is met at just 41% of the sites in the two ESMs with the carbon–nitrogen cycle. In contrast, the second requirement is met at 98% of the sites over the boreal forest in these two ESMs, while the percentage is much lower (41–62%) in the other six ESMs (figure 4(b)).

Zoom In Zoom Out Reset image size

Note that more models (than the eight used in this study) and more ensemble members from each model are also available now from CMIP5 [25], but their Rh treatments are similar to those discussed here. Therefore the above conclusions are expected to be little affected, even though quantitative values (e.g., the exact values in figure 4) would be changed.

There is a horizontal scale mismatch between point measurements and model grid box output. There is also a vertical scale mismatch between point measurement of environmental variables (e.g. temperature) and model output based on multiple soil layers. The latter mismatch is less serious for annual (mean or accumulated) variables. For instance, annual mean temperature is relatively close to each other in different soil layers. The Rh data in SRDB-V2 also contain uncertainty due to the use of different methods to partition observed Rs into heterotrophic and autotrophic respirations. Table 1 summarizes these partitioning methods and indicates that the mean Rh:Rs ratio in the SRDB-V2 database varies from 0.45 to 0.65 (and the median ratio varies from 0.44 to 0.70) depending on the method used. These factors are partially responsible for the above model-data inconsistencies, particularly for the direct Rh comparison between observations and model output (figure 2). Their impact on the comparison of the relationship between Rh and other variables (figure 3, tables 2, S2, and S3), however, is expected to be much less. The correlation between annual Rh and Pr, NPP, and SOC and between the logarithm of annual Rh and T are much higher in ESMs than in observations (tables 2 and S2). The annual Rh is also more sensitive to these climatic, vegetation, and edaphic variables in ESMs than in observations, as indicated by the regression slopes in figure 3.

While soil respiration can be observed directly, 'observations' of Rh are derived from a range of methods (table 1). For this reason, there are many more Rs observations than Rh data. Furthermore, both annual and seasonal Rs data are available from SRDB-V2, while seasonal Rh is not available. To make use of Rs observations, we have created two 'reality checks' for ESMs. All models (with or without the nitrogen cycle) fail to pass these tests in all biomes a significant percentage of the time (figure 4). To facilitate more explicit model validation and improvements, modeling groups should provide Rs as a standard output in the future.

The Rh computation in ESMs is dependent on various factors (equation (1)), including the availability of nitrogen [27, 33]. Without the nitrogen constraint, some ESMs tend to produce Rh values much higher than observations (figure 2). The two models that include an explicit nitrogen cycle perform very well in temperate biomes and might be expected to show more realistic mechanistic Rh responses. However, the nitrogen constraint appears to be too severe in the boreal forest (figure 2), which may explain the particularly low global Rh estimates for CCSM4 and NorESM1-M relative to the other models (figure 1). To test the robustness of this conclusion, we have evaluated another ESM incorporating the nitrogen cycle from China (BCC-CSM1.1 [34, 35]) and found that this model compares well with Rh derived from observations. For example, the 25th, 50th (median) and 75th percentiles of annual Rh over the boreal forest are 221, 323 and 341 g C m⁻², in good agreement with observations in figure 2(a). Furthermore, the lower Rh than expected over the boreal forest in CCSM4 and NorESM1-M is unlikely to be related to the turnover times for global soil carbon in these two models: the turnover times across the CMIP5 ESMs range from 11 to 39 years, in which CCSM4 (about 12 years) is in the lower end but NorESM1-M is in the middle of the range [26]. These analyses and the high correlation between Rh and NPP (figure 3(c), table S2) suggest that low Rh estimates in both models are most probably caused by low NPP which is associated with high Ra. Therefore, the treatment of Ra in both models needs to be improved, including the maintenance respiration (as a function of temperature and tissue nitrogen concentration) and growth respiration (proportional to the difference between leaf photosynthesis and maintenance respiration).

Change in CO₂ release from terrestrial ecosystems is dependent on the temperature sensitivity (Q₁₀) of autotrophs and heterotrophs, which has received much attention in recent years with some controversy [18, 22, 36–38]. However, much less attention has been paid to the use of observed intrinsic and apparent Q₁₀ values in ESMs. The intrinsic (or unconfounded) temperature sensitivity (quantifying the inherent kinetic properties of substrate decomposition) is found to converge to 1.4 ± 0.1 based on flux tower measurements across different biomes [37], and the embedded Q₁₀ value of 1.5 used in the Rh computation at each time step in CLM4-CN (equation (1)) falls into this range. The apparent Q₁₀ is typically calculated from the short-term temporal variability of measured Rh and temperature, and it represents the temperature sensitivity confounded by other factors that covary with temperature. The apparent Q₁₀ values in CCSM4 and NorESM1-M (adopting the same CLM4-CN) differ from each other and from the embedded Q₁₀ value of 1.5. The apparent Q₁₀ values vary broadly from 1.43 to 3.43 in the SRDB-V2 database (figure S1 available at stacks.iop.org/ERL/8/034034/mmedia) and differ from the intrinsic values of 1.4 ± 0.1. In contrast, the apparent Q₁₀ using annual Rh and T across all sites in figure S1 is 1.48, coincidentally close to the embedded value in CLM4-CN and the intrinsic Q₁₀ values from flux tower measurements [37]. The Q₁₀ values are also sensitive to the number of covariates considered in observations and ESMs (tables 2 and S2). For the same climate models, the Q₁₀ values in table 2 are also different from those computed based on the linear regression of ln(SOC) with ln(NPP) and T in the table 4 of Todd-Brown et al [26]. Because Q₁₀ can be, and has been, calculated in a variety of ways, it is crucial to clearly explain methods used to calculate Q₁₀ to avoid confusion between the modeling and measurement communities.

In ESMs, the Rh computation is proportional to the base pool-specific fixed decomposition rate, the carbon content in the corresponding pool, and the fraction of the decomposition carbon flux as CO₂ (equation (1)). With such simple treatments, it is not surprising that the annual Rh is much more sensitive to SOC in ESMs than in observations, as indicated by the regression slopes in figures 3(e) and (f). In contrast, there is essentially no statistically significant relationship between SOC and Rh in SRDB-V2 (figure 3(f)) and in an empirical study [10] for several possible reasons. Not all organic compounds are equally accessible to microorganisms [12]. Labile carbon from root exudates and nutrient content in SOC are important in controlling Rh [38, 39]. Rh is differentially controlled within and outside of the rhizosphere and inhibited in low oxygen conditions [15] and is strongly regulated by extracellular microbial enzymes dependent on the dynamics of microbial function, population size and nutrient availability [31, 33, 36, 40–43]. ESMs omit these mechanisms (figure 3(e)) despite including enzyme-driven process in photosynthesis [44]. In some models the problem is more fundamental. For instance, the soil organic matter carbon to nitrogen ratio is nearly constant at 10 in CLM4-CN, whereas the ratio varies in a broad range of 7–26 in different ecosystems in SRDB-V2. To improve the Rh computation, one promising approach is to incorporate a microbial biomass pool and constraints of extracellular enzymes (e.g. [31]) to represent direct coupling between microbes and soil carbon turnover, mimicking the thermal adaptation and acclimation of Rh by predicting the microbial growth efficiency and respiratory capacity [45].

Note that the lack of statistically significant spatial correlation between annual SOC and Rh over various sites in SRDB-V2 is not necessarily counter to the conventional wisdom on the positive temporal correlation between SOC and Rh over a given site; rather a strong linear relationship is not expected [11]. While more observational data are still needed to establish complex and likely non-linear relationships between annual SOC and Rh in contrasting conditions, the simple treatment of the above three factors used in the Rh computation in ESMs needs to be revised. For instance, the constant decomposition rate in equation (1) may be relaxed to be dependent on soil depth, and the simple carbon allocation rule may need to consider ecosystem acclimation and adaptation.

Because Rh parameterizations in current ESMs lack key regulatory processes like the role of different carbon substrates, microbial dynamics and nutrient constraints, the global Rh range is as large as 17.4 Pg C a⁻¹ (=71.6–54.2) from 1976 to 2005 using the six carbon-only ESMs, and it increases to 40.1 Pg C a⁻¹ from 2071 to 2100 under the RCP8.5 scenario. When the two models with the carbon–nitrogen cycle are also included, these ranges become even larger, reaching 30.3 and 65.5 Pg C a⁻¹, respectively. To reduce such large uncertainties in ESM Rh projections, the observational and ESM community need to work together to design field experiments for testing and improving the Rh parameterization as a function of soil temperature, moisture, and other processes (e.g. equation (1)), rather than just the Q₁₀ factor. Seasonal Rh data, which are not available from SRDB-V2 but are very useful for model Rh evaluations, need to be developed and observations of Rh and likely drivers of microbial metabolism in contrasting conditions would be instructive. In addition to comparing measured flux rates to those predicted by models, model development would be aided by benchmarking mean residence times of C pools used by models against estimates from the field or laboratory based on isotopic or incubation experiments.