Microbial dormancy and its impacts on northern temperate and boreal terrestrial ecosystem carbon budget

A large amount of soil carbon in northern temperate and boreal regions could be emitted as greenhouse gases in a warming future. However, lacking detailed microbial processes such as microbial dormancy in current biogeochemistry models might have biased the quantification of the regional carbon dynamics. Here the effect of microbial dormancy was incorporated into a biogeochemistry model to improve the quantification for the last century and this century. Compared with the previous model without considering the microbial dormancy, the new model estimated the regional soils stored 75.9 Pg more C in the terrestrial ecosystems during the last century and will store 50.4 and 125.2 Pg more C under the RCP8.5 and RCP2.6 scenarios, respectively, in this century. This study highlights the importance of the representation of microbial dormancy in earth system models to adequately quantify the carbon dynamics in the northern temperate and boreal natural terrestrial ecosystems.


Introduction 76
The land ecosystems in northern temperate and boreal regions (>45 ºN) occupy 22% of 77 the global surface and store over 40% of the global soil organic carbon (SOC) ( heterotrophic respiration data at representative sites. Third, we presented how the model was 148 applied to natural ecosystems in the region (above 45 ºN) for the 20 th and 21 st centuries and 149 discussed the dormancy effects on their regional carbon budget. with existing data of carbon stocks and fluxes, our study incorporated the microbial module into 159 an extant MIC-TEM that simulates carbon data dynamically. This coupling enables us to 160 extrapolate our model to northern temperate and boreal terrestrial ecosystems, rather than only 161 for temperate forest region in He et al. (2015). In our new model (MIC-TEM-dormancy), 162 microbial biomass pool was divided into two fractions, including the dormant and active 163 microbial biomass pools. The two microbial biomass pools and the reversible transition between 164 them have been considered explicitly in the new model (Figure 1), which was ignored in MIC-165

TEM. 166
In previous MIC-TEM, heterotrophic respiration (RH) is simply calculated as the product of 167 ASSIM and CUE, which are microbial assimilation and carbon use efficiency, respectively. For 168 detailed carbon dynamics in MIC-TEM, see Zha & Zhuang (2018). 169 Here we revised MIC-TEM by incorporating microbial dormancy dynamics according to 170 He et al. (2015). In MIC-TEM-dormancy, the soil heterotrophic respiration RH is revised to include 171 three parts: the maintenance respiration from the active and dormant microorganisms and the CO2 172 production through the process of microbial assimilation (He et al., 2015): 173 Here parameter  is maintenance weight (h -1 ), CNsoil and CNmic denotes the C:N ratios of soil and 189 that of microbial biomass. Besides, Φ is the substrate saturation level and defined as in He et al. 190 (2015) and Wang et al. (2014): 191 Where Ks is the half saturation constant for substrate uptake as indicated by the Michaelis-Menten 193 kinetic, and S is soluble C substrates that are directly accessible for microbial assimilation (Wang 194 The term 'Soluble C' denotes the state variable of soluble carbon pool. Dliq is the diffusion 198 coefficient of the substrate in the liquid phase, and is formulated as: 199 Where BD is the bulk density and PD is the soil particle density. θ is the volumetric soil moisture.  With the modification of microbial carbon dynamics by considering microbial life-history trait, 237 soil decomposition is changed since it is controlled by microbes. When microbial dormancy is 238 considered, the number of active microbes that participate in soil decomposition is much less. The 239 changes in soil decomposition directly influence the amount of soil respiration, and further 240 influence soil nitrogen (N) mineralization that determines soil N availability for plants, affecting 241 gross primary production (GPP). Since both GPP and RH can be affected by microbial dormancy, 242 net ecosystem production (NEP) will also be affected. 243 244

Model parameterization and validation 245
The detailed description of parameters that are related to microbial dormancy can be found 246 in He et al. (2015) ( Table 1). Here we calibrated the MIC-TEM-dormancy at six representative 247 sites with gap-filled monthly net ecosystem productivity (NEP, gCm -2 mon -1 ) data in northern 248 temperate and boreal regions (Table 2). Site-level climatic data and soil texture data were organized 249 for driving model. All sites information can be found on AmeriFlux network (Davidson et al., performed based on prior ranges from literature (Table 1) to minimize the difference between the 254 monthly simulated and measured NEP at the chosen sites. The cost function of the minimization 255 is: 256 Where NEP obs,i and NEP sim,i are the observed and simulated NEP, respectively. k is the number of 258 data pairs for comparison. Except for the parameters of microbial dormancy, other parameters are 259 derived directly from MIC-TEM (Zha & Zhuang, 2018). The optimized parameters were used for 260 model validation and regional simulations. 261 For model validation, we chose another six sites that containing monthly NEP data from 262 AmeriFlux network (Table 3). Four of these six sites were also used for parameterization (Table  263 2). However, we used the data of different observation periods for model validation for those 264 overlapped sites. Moreover, we also conducted site-level validations with monthly soil respiration 265 data from AmeriFlux network and Fluxnet dataset. The site information was provided in Table 4. 266 For these sites, we assumed 50% of soil respiration was heterotrophic respiration (RH) for forest 267 (Hanson et al., 2000), 60% and 70% of that was RH for grassland (Wang et al., 2009) and tundra 268 (Billings et al., 1977). Because there is a limited amount of available RH data, we could not 269 conduct a regional validation for all pixels in northern temperate and boreal regions. Instead, we from literatures. In our model, we assumed that soil texture, elevation, and potential natural 297 vegetation data only vary spatially, not vary over time (Zhuang et al., 2015). 298 In addition, regional simulations over the 21 st century were conducted under two 299 Intergovernmental Panel on Climate Change (IPCC) climate scenarios (RCP 2.6 and RCP 8.5). 300 The future climatic forcing data under these two climate change scenarios were derived from the 301 HadGEM2-ESmodel, which is a member of CMIP5project213 (https://esgf-302 node.llnl.gov/search/cmip5/). Then the regional estimations were obtained by summing up the 303 gridded outputs for our study region. The positive simulated NEP represents a CO2 sink from the 304 atmosphere to terrestrial ecosystems, while a negative value represents a source of CO2 from 305 terrestrial ecosystems to the atmosphere. 306 307

Parameter equifinality effects 308
Our previous studies using TEM has demonstrated that equifinality derived from site-level 309 parameterization will affect the uncertainty in the estimation of regional carbon dynamics (Tang 310 and Zhuang, 2008Zhuang, , 2009). Here equifinality refers to that a number of sets of parameters result in 311 model simulations that all match the data similarly well. To quantify this effect on our simulation 312 uncertainty, we conducted ensemble regional simulations with 50 sets of parameters for both 313 historical and future studies. The 50 sets of parameters were obtained according to the method in 314 Tang and Zhuang (2008). 315 316 3. Results 317

Inversed Model Parameters and model validation 318
Using SCE-UA ensemble method, 50 independent sets of parameters were converged to 319 minimize the objective function. Then the optimized parameters are calculated as the mean of these 320 50 sets of inversed parameters. The boxplot of parameter posterior distributions reflects different 321 ecosystem properties at these sites ( Figure 3). For instance, growth yield was higher in tundra types 322 than in forests, meaning microorganisms in environment with higher energy limitation tend to 323 enhance the efficiency of energy transportation. Besides, alpha, the maintenance weight, was also 324 higher in tundra types than in forests. From the plot for parameter beta, the ratio of dormant 325 maintenance rate to specific maintenance rate for active biomass in tundra types is lower than that 326 in forest types. Other microbial related parameters did not differentiate much among different 327 vegetation types. 328 After parameterization, the MIC-TEM-dormancy was validated with monthly NEP data for 329 six representative ecosystems, and the comparisons between monthly observed NEP and 330 simulated NEP were presented in Figure 4. With the optimized parameters, the dormancy-based 331 model was used to reproduce NEP to compare with the measured NEP (Table 5). The R 2 ranges 332 from 0.67 for Atqasuk to 0.93 for Bartlett Experimental Forest (Table 5). Generally, our new 333 model performs better for forest ecosystems than for tundra ecosystems. Compared with MIC-334 TEM, dormancy model performs better for alpine tundra, temperate coniferous forest, and 335 grassland. For other sites, both models show similar performance (Table 5). Besides, a set of 336 monthly soil respiration data were selected to evaluate the estimated RH. The comparisons 337 between monthly observed RH and simulated RH from two contrasting models were conducted 338 ( Figure 5). MIC-TEM-dormancy has higher R 2 and lower root mean square error (RMSE) ( Table  339 6). Sixty-one sites with average annual RH in northern temperate and boreal regions were used to 340 further evaluate the new model performance. The dormancy model has lower intercept and slope 341 with R 2 of 0.45, while R 2 of MIC-TEM is 0.3 ( Figure 6). These analyses indicate that new model 342 is more realistic in representing RH by considering microbial dormancy. 343 344 3.2 Regional carbon dynamics during the 20 th century 345 Regional extrapolation with both models estimated a regional terrestrial ecosystem carbon sink 346 but with different magnitudes (Figure 7c). With optimized parameters, MIC-TEM estimated a  Pg C yr -1 and 6.4 Pg C yr -1 , which are 5.8% and 16.3% less than the estimations from MIC-TEM, 355 respectively (Figures 7a and 7b). This pronounced difference of NEP between two models comes 356 from the disparity between the simulated NPP and RH with them since NEP is calculated as the 357 difference between NPP and RH. Without considering dormancy, MIC-TEM estimates more active 358 microbial biomass, hence overestimating both RH and NPP (due to higher simulated N 359 mineralization and uptake by plants), but resulting in lower NEP than that calculated by MIC- Pg C with an interannual standard deviation of 0.1 Pg C yr -1 (Figure 9). MIC-TEM-dormancy 389 estimates NPP and RH at 9.9 Pg C yr -1 and 8.7 Pg C yr -1 , which are 0.5 Pg C yr -1 and 1.7 Pg C yr -390 1 less than the estimations from MIC-TEM, respectively (Figure 9). Moreover, simulations under 391 the two contrasting climate scenarios (RCP 2.6 and RCP 8.5) exhibit a large difference of 81.1 392 Pg C of cumulative NEP during the 21 st century by MIC-TEM, but only 6.3 Pg C of that by 393

MIC-TEM-dormancy. 394
MIC-TEM-dormancy estimated higher RH in winter, but lower RH in summer under both 395 future scenarios in the 2090s (Figure 10). NPP is the same in winter with or without dormancy, 396 and in the late summer is higher than that without dormancy, especially in the RCP 8.5 scenario.

Regional uncertainty considering equifinality effects during 20 th and 21 st centuries 404
The ensemble simulations for the 20 th century is shown in Figure 11. Given the 405 uncertainty in parameters, MIC-TEM-dormancy predicts that the regional cumulative carbon The large bias between dormancy and non-dormancy models mainly comes from two parts.  (25°N-50°N) from the dormancy model amounted to 7.28 Pg C yr -1 436 and 8.83 Pg C yr -1 from a no-dormancy model, which is 21.3% higher than the dormancy model. 437 Although their study region and simulation period are different from our study, the results can 438 still be comparable. Both studies indicated that the magnitude of RH from no-dormancy model 439 are higher than dormancy models. Second, high soil respiration stimulates N mineralization in 440 soils (Zhuang et al., 2001(Zhuang et al., , 2002, making more nutrients for photosynthesis of plants (Raich et  Therefore, NPP will be higher due to the N enrichment from higher RH. However, how NEP  We didn't separate among functional microbial groups, but gather microbes into one "box". Our study considered the CUE sensitivity to temperature, but not nutrient availability. On the 479 other hand, some model assumptions can also cause uncertainties. For example, we assumed that 480 vegetation will not change during the transient simulation. However, over the past few decades 481 in northern temperate and boreal regions, temperature increases have led to vegetation shift from 482 one type to another (Hansen et al., 2006;White et al., 2000). The vegetation changes will affect 483 carbon cycling in these ecosystems. 484

While our analysis suggests it is important to incorporate microbial dormancy dynamics 485
into a process-based biogeochemistry model to more adequately simulate carbon dynamics in 486 northern temperate and boreal regions, we do confront modeling dilemmas. First, our process-487 based models have a relatively large number of parameters, which unavoidably creates the 488 "equifinality" problem as recognized in our previous studies for the model (e.g., 489 Zhuang, 2008, 2009). To alleviate this problem in this analysis, we have conducted parameter 490 ensemble simulations at both site and regional levels and presented our results with uncertainties, 491 which could be a standard approach for process-based complex biogeochemistry modeling 492 analyses. Second, incorporating more ecosystem processes increases the number of parameters 493 in our model, inducing even larger uncertainties for both site level and regional simulations. On 494 the one hand, the more complex model to a certain degree helps capture observations, on the 495 other hand, the model uncertainty has not been constrained or even enlarged. We highlight the 496 need to further investigate this trade-off within the modeling research community. 497 498

Conclusions 499
This study incorporated microbial dormancy into a detailed microbial-based soil 500 decomposition biogeochemistry model to examine the fate of large soil carbon storage in 501 northern temperate and boreal natural terrestrial ecosystems under changing climate conditions. 502 Regional simulations using MIC-TEM-dormancy indicated that, over the 20 th century, the region 503 is a carbon sink of 166.8 ± 97.7 Pg. This sink could decrease to 175.9 ± 105.4 Pg under the RCP 504 8.5 scenario or 125.4 ± 85.5 Pg under the RCP 2.6 scenario during the 21 st century. Whether 505 considering microbial dormancy or not can cause large differences in soil decomposition 506 estimation between two models. Meanwhile, due to available nitrogen affected by soil 507 decomposition, net primary production is consequently influenced in these two centuries. The 508 combined changes in soil decomposition and net primary production led to large differences in 509 carbon budget estimation between two models. Compared with MIC-TEM, MIC-TEM-510 dormancy projected 75.9 Pg more C stored in the terrestrial ecosystems over the last century,