UC Berkeley
UC Berkeley Previously Published Works
Title
Higher climatological temperature sensitivity of soil carbon in cold than warm climates
Permalink
https://escholarship.org/uc/item/7q68v57x
Journal
Nature Climate Change, 7(11)
ISSN
1758-678X
Authors
Koven, Charles D
Hugelius, Gustaf
Lawrence, David M
et al.
Publication Date
2017-11-01
DOI
10.1038/nclimate3421
Peer reviewed
Powered by the California Digital Library
University of California
Higher Climatological Temperature Sensitivity of Soil Carbon in Cold Than Warm
Climates
C. D. Koven
1*
, G. Hugelius
2,3
, D. M. Lawrence
4
, W. R. Wieder
4,5
5
Affiliations:
1
Earth and Environmental Sciences Division, Lawrence Berkeley National Laboratory, Berkeley,
CA, USA
10
2
Department of Physical Geography & Bolin Centre of Climate Research, Stockholm University,
Stockholm, Sweden
3
School of Earth, Energy, and Environmental Sciences, Stanford University, Stanford, CA, USA
15
4
Climate & Global Dynamics Laboratory, National Center for Atmospheric Research, Boulder,
CO, USA
5
Institute of Arctic & Alpine Research, University of Colorado, Boulder, CO, USA
20
*Correspondence to: cdkoven@lbl.gov
2
The projected loss of soil carbon to the atmosphere resulting from climate change is a potentially
25
large but highly uncertain feedback to warming. The magnitude of this feedback is poorly
constrained by observations and theory, and is disparately represented in Earth system models
(ESMs)
1-3
. To assess the climatological temperature sensitivity of soil carbon, we calculate
apparent soil carbon turnover times
4
that reflect long-term and broad-scale rates of
decomposition. Here, we show that the climatological temperature control on carbon turnover in
30
the top meter of global soils is more sensitive in cold climates than in warm ones and argue that
it is critical to capture this emergent ecosystem property in global-scale models. We present a
simplified model that explains the observed high cold-climate sensitivity using only the physical
scaling of soil freeze-thaw state across climate gradients. Current ESMs fail to capture this
pattern, except in an ESM that explicitly resolves vertical gradients in soil climate and C
35
turnover. An observed weak tropical temperature sensitivity emerges in a different model that
explicitly resolves mineralogical control on decomposition. These results support projections of
strong carbon-climate feedbacks from northern soils
5,6
and demonstrate a method for ESMs to
capture this emergent behavior.
40
Carbon cycle feedbacks represent a large uncertainty on the terrestrial response to climate
change
1-3
. Much of this uncertainty arises from the dynamics of decomposing soil carbon under
changing climate, in particular how the rate of carbon cycling through soils may change with
warming. Fast-timescale observations
7
and general kinetic theory
8
both suggest that
decomposition rates should increase with warming. While this temperature response has long
45
been thought to provide a positive feedback to warming
9
, its magnitude is poorly quantified due
to the many confounding factors affecting soil metabolic rates
8
. Furthermore, acclimatory
responses by soil microbiota that reduce the effect of warming on decomposition rates at longer
timescales have been proposed
10-12
to explain the reduction in temperature sensitivity observed in
experiments
13
. Given the size of global soil C stocks, especially at high latitudes
14
and the
50
potential long-term vulnerability of soil C to warming, it is critical to accurately include these
feedbacks in assessing emissions scenarios that are compatible with desired climate outcomes
15
.
The current (Coupled-Model Intercomparison Project, phase 5; CMIP5) generation of ESMs
actually show a relatively small contribution of climate-driven changes to carbon turnover times
3
on soil carbon stocks, with the majority of projected soil carbon change instead occurring due to
55
changes in plant productivity
3
. The ESMs—which couple carbon cycle and climate processes—
typically use simplified temperature sensitivities
12
, omit realistic dynamics of soil microbial
ecology
16
, show little predictive power in simulating current soil carbon stocks
17
, and
systematically overestimate the transient sensitivity of soil carbon pools to productivity
changes
18
. Equally troubling, the models omit crucial processes that may exacerbate warming-
60
related carbon losses, in particular the representation of frozen carbon in permafrost soils
6
. Thus,
although ESMs show high inter-model divergence in soil carbon predictions, they likely
underestimate the actual uncertainty
16
surrounding increased atmospheric greenhouse gas
burdens and accelerated warming under climate change scenarios.
It is difficult to directly evaluate the transient climate-response predictions made by ESM
65
soil carbon models, because dynamical observations of soil carbon responses to warming at the
relevant timescale—multi-decadal to centennial—are scarce
16
and ambiguous about both
acclimation timescales and whether changes result from productivity or turnover responses
19
. An
alternate approach is to look at model predictions across spatial gradients, since current soil
conditions reflect the long-term accumulated effects of climate, vegetation, edaphic properties,
70
and landscape changes on soil organic matter formation
20
(see Methods). Indeed, the large
spatial variation of soil carbon turnover times across climate gradients served as an early piece of
evidence supporting the idea that warming would lead to soil carbon losses
21,22
.
Because the balance of carbon inputs and decomposition determine soil carbon stocks, and
because both of these controls are mediated by climate, it is useful to separate them by defining
75
an apparent turnover time,
τ
, as the ratio of carbon losses via heterotrophic respiration to total
carbon stocks. Since, at steady state, carbon losses and inputs are equal, and because we have
more robust global estimates of productivity than of heterotrophic respiration, we assume that
soils are approximately at steady state in order to estimate
τ
as the ratio of carbon stocks to
carbon inputs
4
.
80
In figure 1 we show the global distribution of soil carbon stocks (fig. 1a), vegetation
inputs to soil (fig. 1b), and
τ
,
as a function of temperature and precipitation (fig. 1c). These
results illustrate that
τ
is
clearly sensitive to both soil temperature and moisture, but here we are
primarily interested in identifying the temperature control on soil carbon turnover. Moisture may
4
dominate turnover in places that are either highly moisture-limited
4
or saturated, so we mask all
85
points where we expect moisture to exert a dominant control (fig. S1), to isolate the temperature-
dominated soil carbon
τ
response (fig. 1d).
This relationship clearly shows a change in the
sensitivity of inferred
τ
to climatological temperature over the interval, with stronger sensitivity
in cold climates than in warm ones. We note, however, that considerable variation remains. The
residual two-fold variation in turnover times (residual variance in log(
τ
) = 0.1; table S1) is also a
90
real and important feature of the data, and this may be driven by mineralogical or other
factors
23,24
beyond the simple climate metrics used here. We recognize that further research
diagnosing the mechanisms responsible for this variation is critical, but here we focus on the
central relationship between soil
τ
and temperature that emerges from our global analysis.
Taking the derivative of the central soil
τ
to temperature relationship (fig. 1d), and placing
95
this in terms of the exponential form Q
10
, gives a “climatological Q
10
” (fig. 2), which
decreases
with temperature, from Q
10
>5 in cold climates to Q
10
=1 (i.e. no temperature sensitivity) in hot
climates. This climatological Q
10
differs from the classical short-timescale Q
10
in being
diagnosed from the
τ
, whereas short-timescale Q
10
values are diagnosed based on instantaneous
decay rates,
k
(where, at steady state,
k
= 1/
τ
). Short-timescale respiration observations show a
100
widespread Q
10
-like behavior with a value in the range of approximately 1.4
based on eddy
covariance fluxes
7
, or 1.5-2 based on soil incubations
8
. Where the climatological
Q
10
value
differs from the short-
timescale Q
10
value
, this is evidence for emergent behavior at longer
timescales that leads to the divergence between short-
and long-
timescale
temperature
sensitivities. We thus divide the world heuristically into three roughly defined regimes (fig. 2): a
105
cold-climate high-sensitivity emergent domain (climatological
Q
10
~2 to 5), a temperate non-
emergent domain (climatological
Q
10
~1.4 to 2), and possibly a warm-climate low-sensitivity
emergent domain (climatological
Q
10
<1.4).
To understand why the climatological
Q
10
exceeds short-timescale
Q
10
values in cold
climates, we explore a hierarchy of simplified decomposition models. These derive
τ
values
110
based solely on modeled soil temperature (T) dynamics and differ only in the functional form of
k
(T) and how this is scaled in space and time to calculate
τ
(fig. 3).
We consider four cases. The
first three differ only in the form of
k
(T) and all use near-surface soil temperatures (10 cm): fixed
Q
10
over the entire temperature range, (fig. 3a); temperature-sensitive Arrhenius relationship
25
over the entire temperature range, (fig. 3b); fixed Q
10
over thawed temperatures and no
115
5
respiration when soil is frozen (fig. 3c). The fourth case shown uses the same temperature
function as the third, but diagnoses
k
using depth-resolved soil temperatures over the full 0-1m
depth interval (fig. 3d). Only the fourth case is able to qualitatively capture the observed change
in slope seen in the observations.
The implication of the curve in figure 3d is that the increase in temperature sensitivity in
120
cold climates can be explained simply as a result of the combined reductions in the thawed
season length and thawed depth during the warm season. Thus, vertical variation in soil climate
must be accounted for to explain carbon stocks at high latitudes, even when considering carbon
to only 1m depth. Models that explicitly resolve this process, by replacing traditional carbon
cycle ordinary differential equations (ODEs) with a set of vertically resolved partial differential
125
equations (PDEs) that include transport, are one such approach, but suffer from high uncertainty
in the rates of cold-soil vertical mixing processes, such as cryoturbation
6
. Those long-term
mixing processes also contribute significantly to the large amounts of carbon stored even deeper,
below 1 meter depth
14,26
, which add further uncertainty and bias to global carbon cycle
projections. That a simpler ODE approach using a depth-averaged
k
approximates the observed
130
relationship suggests that, at least in the near surface, such transport processes are sufficiently
fast over long timescales for the soil to act as a well-mixed reservoir through which respiration
can occur at any depth within the 0-1m interval.
We contend that the climatological sensitivity of soil C to historical climate (fig. 1d) is an
emergent ecosystem property that models should be expected to replicate. To test whether ESMs
135
are able to match these qualitative patterns, we compare predictions of
τ
from models used in the
CMIP5
carbon cycle experiments (fig. 4a-f, table S2)
1
. Most models show a linear relationship
between log(
τ
) and MAAT, as would result from using fixed Q
10
and a single-layer model that
diagnoses
k
values from near-surface temperatures. Some models show offsets and emergent
behavior, but none are able to qualitatively capture both the increase in temperature sensitivity
140
through the entire range of cold climates as well as the reduction in temperature sensitivity in
tropical climates shown by the global data. The inability of the models to match spatial gradients
implies that the transient response to warming will likewise be biased, particularly in the cold-
climate regime where the ESMs show a systematic underestimate of the climate sensitivity on
soil carbon turnover. That the ESMs also show weak turnover-driven soil carbon responses to
145
6
warming
3
, and a net high-latitude carbon gain from warming
2
, is thus likely a shared artifact of
the weak temperature control in these ESMs.
Because none of the CMIP5 models represent permafrost carbon dynamics, we diagnose the
same relationship from a model which does represent permafrost carbon via a PDE approach,
CLM4.5
27
(fig. 4g-h). The two panels 4g and 4h differ only by a parameter, Z
τ
,
that controls
150
decomposition rates
k
at depth, beyond the resolved climatologic controls. The parameter is
defined as an e-folding depth, and a short value assumes that the base decomposition rate (
k
) of
deeper soil horizons is intrinsically slower than surface soil carbon stocks (Fig. 4g). A long e-
folding depth assumes that carbon pools in surface and deep soil horizons have similar intrinsic
decomposability, which allows climatologic controls (e.g., temperature) to more strongly
155
influence decomposition rates (Fig. 4h). The model is able to match the observed increase in
sensitivity at cold climates, but only when the base rates of deep decomposition are more similar
to those at the surface (fig. 4h). Under these conditions CLM4.5 predicts a substantial (-23 Pg
C/
˚C
) destabilizing carbon-climate feedback from the permafrost region
28
, which contrasts in
sign with the stabilizing feedback projected by the CMIP5 ESMs from this region. The patterns
160
of
τ
in fig. 4h better match observations in fig. 1 than other models analyzed here (Table S1),
suggesting that the corresponding projection of a strong permafrost carbon-climate feedback is
also more realistic.
The wide spread in
τ
(Fig 1d) and the low-sensitivity emergent domain observed in warm
climates (Fig. 2) emerges from a model that includes mineral and microbial associations (fig. 4i).
165
This mineral/microbial model predicts longer
τ
values for the clay rich tropical soils
and
therefore a reduced sensitivity to temperature. Because weathering rates increase with
temperature in sufficiently moist ecosystems, clay amounts and temperature are positively
correlated, which may explain the reduced tropical temperature sensitivity (fig. S3). However,
the reduced sensitivity of tropical soil carbon is also consistent with nonlinear models that
170
predict a temperature optimum for decomposition
29
, though only at tropical temperatures
.
Separating these potential causes is not possible with the static benchmark proposed here, but is
of great importance, as they would lead to different trajectories under global warming.
We propose that global temperature control on turnover, as expressed in spatial gradients,
is a useful benchmark on dynamic models, which must predict these static relationships if we are
175
7
to have confidence in their transient responses. We recognize that projections of soil C response
to warming are complicated by changes in plant productivity, organo-mineral stabilization, soil-
aggregate formation and shifts in belowground community structure and function
30
. Indeed,
resolving these complex interactions should be a focus for the next generation of experiments
and models
16
. The results shown here stress the importance of considering temperature effects on
180
decomposition across the full range of climates found on Earth, as well as vertically, even within
the top meter of soils. The systematic underestimation by the CMIP5 models of the sensitivity of
carbon turnover in cold climates belies their projections of weak soil turnover-driven feedbacks
to warming. Thus the relationships shown here support stronger carbon – climate feedbacks—
particularly from northern regions—than current estimates suggest.
185
8
Figures
190
Fig. 1
.
Global distributions of the inferred apparent turnover time (
τ
) of global soil organic
matter as function of climatological temperature
.
τ
is calculated as the ratio of (a) carbon
stocks to (b) net primary productivity. (c)
τ
plotted as function of Mean Annual Air Temperature
195
(MAAT). Each gridcell is colored by climatological precipitation. (d) As in (c), but after filtering
out gridcells that are likely to be dominated by either aridity (P minus PET < threshold of -1000
9
mm/yr) or saturation (peatland fraction exceeds threshold of 50%)
. Best fit regression curve in
(d) uses a quadratic regression of log(
τ
) versus MAAT, with 50% prediction intervals shown.
200
10
Fig 2.
Inferred “climatological Q
10
” as a function of temperature
. Climatological Q
10
is
calculated from the derivative of the regression relationship between
τ
and MAAT in fig. 1d. We
205
define emergent domains as those where the climatological Q
10
differs appreciably from short-
term Q
10
values (i.e Q
10
> 2 or Q
10
< 1.4).
11
210
Fig. 3
A hierarchy of simplified models to explain the cold-climate emergent regime of high
climatological temperature sensitivities
. In all cases,
τ
is calculated using daily soil
temperatures from a land surface model, applying a simplified decomposition function to derive
daily decomposition rates (
k
), and inferring
τ
as the reciprocal of the mean decomposition
215
function
k
. (a) Simple Q
10
=1.5 function evaluated. (b) Arrhenius temperature function. (c)
Thawed-only Q
10
=1.5 function evaluated at 10cm depth. (c) Thawed-only Q
10
=1.5 function
12
averaged over the surface-1m depth interval. Blue lines are the decomposition function as
evaluated on MAAT. Red lines are the best-fit curve and prediction intervals are from fig. 1d.
13
220
Fig. 4
A comparison of relationships between soil turnover times and climate as predicted
by a suite of ESMs and offline land models.
Inferred apparent turnover time,
τ
, calculated as in
figure 1 and colored by precipitation as in fig 1c, from soil models used in Earth system models.
(a-f) CMIP5 models, each of which (other than GFDL-ESM2G) use single-layer soil temperature
control on soil carbon turnover. (g-h) CLM4.5, which calculates vertically-resolved
225
14
decomposition rates. (g) and (h) differ by varying a parameter (Z
τ
)
that controls decomposition
rates with depth independently from resolved temperature, moisture, and oxygen controls (i)
MIMICS, which treats decomposition as a microbially-enabled and mineral-resolved nonlinear
model, shows the wide scatter in moist tropical climates as observed, due to its consideration of
mineralogical control on decomposition.
230
Acknowledgments:
This research was supported by the Director, Office of Science, Office of
Biological and Environmental Research of the US Department of Energy (DOE) under Contract
DE-AC02-05CH11231 as part of their Regional and Global Climate Modeling (BGC-Feedbacks
SFA), and Terrestrial Ecosystem Science Programs (NGEE-Arctic and NGEE-Tropics), and used
235
resources of the National Energy Research Scientific Computing Center, also supported by the
Office of Science of the US Department of Energy, under Contract DE-AC02- 05CH11231.
National Center for Atmospheric Research (NCAR) is sponsored by the National Science
Foundation (NSF). The CESM project is supported by the NSF and the Office of Science (BER)
of the US Department of Energy. Computing resources were provided by the Climate Simulation
240
Laboratory at NCAR’s Computational and Information Systems Laboratory, sponsored by NSF
and other agencies. G. H. acknowledges funding from the Swedish Research Council (grant
numbers E0689701 and E0641701), the EU JPI-climate COUP project and Marie Sklodowska
Curie Actions, Cofund, Project INCA 600398. D.M.L. is supported by funding from the US
Department of Energy BER, as part of its Climate Change Prediction Program, Cooperative
245
Agreement DE-FC03-97ER62402/ A010 and NSF Grants AGS-1048996 and ARC-1048987.
W.R.W. is supported by funding from the US Department of Agriculture NIFA 2015-67003-
23485 and US Department of Energy TES DE-SC0014374. We acknowledge the World Climate
Research Programme’s Working Group on Coupled Modeling, which is responsible for CMIP,
and we thank the climate modeling groups (listed in Table S2 of this paper) for producing and
250
making available their model output.
Author Contributions
15
CDK designed the study and performed analyses, based on ideas developed through discussions
with DML, GH, and WRW. WRW contributed MIMICS results, CDK and DML contributed
255
CLM4.5 results, and GH contributed NCSCD data. All authors wrote the manuscript.
16
References (Main Text)
260
1
Arora, V. K.
et al.
Carbon-concentration and carbon-climate feedbacks in CMIP5 Earth system models.
Journal of Climate
26
, 5289-5314, doi:10.1175/JCLI-D-12-00494.1 (2013).
2
Ciais, P.
et al.
in
Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to
the Fifth Assessment Report of the Intergovernmental Panel on Climate Change
(eds T.F. Stocker
et al.
)
265
(Cambridge University Press, 2013).
3
Koven, C. D.
et al.
Controls on terrestrial carbon feedbacks by productivity versus turnover in the CMIP5
Earth System Models.
Biogeosciences
12
, 5211--5228, doi:10.5194/bg-12-5211-2015 (2015).
4
Carvalhais, N.
et al.
Global covariation of carbon turnover times with climate in terrestrial ecosystems.
Nature
514
, 213-217, doi:10.1038/nature13731 (2014).
270
5
Ciais, P.
et al.
Large inert carbon pool in the terrestrial biosphere during the Last Glacial Maximum.
Nature
Geosci
5
, 74-79, doi:10.1038/ngeo1324 (2012).
6
Koven, C. D.
et al.
Permafrost carbon-climate feedbacks accelerate global warming.
Proceedings of the
National Academy of Sciences
108
, 14769--14774, doi:10.1073/pnas.1103910108 (2011).
7
Mahecha, M. D.
et al.
Global Convergence in the Temperature Sensitivity of Respiration at Ecosystem
275
Level.
Science
329
, 838--840 (2010).
8
Davidson, E. & Janssens, I. Temperature sensitivity of soil carbon decomposition and feedbacks to climate
change.
Nature
440
, 165-173, doi:DOI 10.1038/nature04514 (2006).
9
Jenkinson, D., Adams, D. & Wild, A. Model Estimates Of CO2 Emissions from soil in response to global
Warming.
Nature
351
, 304--306 (1991).
280
10
Allison, S. D., Wallenstein, M. D. & Bradford, M. A. Soil-carbon response to warming dependent on
microbial physiology.
Nature Geosci
3
, 336--340 (2010).
11
Wieder, W. R., Bonan, G. B. & Allison, S. D. Global soil carbon projections are improved by modelling
microbial processes.
Nature Clim. Change
3
, 909--912, doi:10.1038/nclimate1951 (2013).
12
Tang, J. & Riley, W. J. Weaker soil carbon-climate feedbacks resulting from microbial and abiotic
285
interactions.
Nature Clim. Change
5
, 56-60, doi:10.1038/nclimate2438 (2015).
13
Bradford, M. A.
et al.
Thermal adaptation of soil microbial respiration to elevated temperature.
Ecology
Letters
11
, 1316--1327, doi:10.1111/j.1461-0248.2008.01251.x (2008).
14
Hugelius, G.
et al.
Estimated stocks of circumpolar permafrost carbon with quantified uncertainty ranges
and identified data gaps.
Biogeosciences
11
, 6573--6593, doi:10.5194/bg-11-6573-2014 (2014).
290
15
Jones, C.
et al.
Twenty-First-Century Compatible CO2 Emissions and Airborne Fraction Simulated by
CMIP5 Earth System Models under Four Representative Concentration Pathways.
Journal of Climate
26
,
4398--4413, doi:10.1175/JCLI-D-12-00554.1 (2013).
16
Bradford, M. A.
et al.
Managing uncertainty in soil carbon feedbacks to climate change.
Nature Clim.
Change
6
, 751--758 (2016).
295
17
Todd-Brown, K. E. O.
et al.
Causes of variation in soil carbon simulations from CMIP5 Earth system
models and comparison with observations.
Biogeosciences
10
, 1717--1736, doi:10.5194/bg-10-1717-2013
(2013).
18
He, Y.
et al.
Radiocarbon constraints imply reduced carbon uptake by soils during the 21st century.
Science
353
, 1419-1424 (2016).
300
19
Crowther, T. W.
et al.
Quantifying global soil carbon losses in response to warming.
Nature
540
, 104--108
(2016).
20
Jenny, H.
Factors of soil formation: a system of quantitative pedology
. (Courier Corporation, 1941).
21
Post, W. M., Emanuel, W. R., Zinke, P. J. & Stangenberger, A. G. Soil carbon pools and world life zones.
Nature
298
, 156--159 (1982).
305
22
RAICH, J. W. & SCHLESINGER, W. H. The global carbon dioxide flux in soil respiration and its
relationship to vegetation and climate.
Tellus B
44
, 81--99, doi:10.1034/j.1600-0889.1992.t01-1-00001.x
(1992).
23
Doetterl, S.
et al.
Soil carbon storage controlled by interactions between geochemistry and climate.
Nature
Geosci
8
, 780--783 (2015).
310
24
Giardina, C. P. & Ryan, M. G. Evidence that decomposition rates of organic carbon in mineral soil do not
vary with temperature.
Nature
404
, 858--861 (2000).
17
25
LLOYD, J. & TAYLOR, J. ON THE TEMPERATURE-DEPENDENCE OF SOIL RESPIRATION.
Functional Ecology
8
, 315-323 (1994).
26
Harden, J. W.
et al.
Field Information Links Permafrost Carbon to Physical Vulnerabilities of Thawing.
315
Geophys. Res. Lett.
39
, doi:10.1029/2012GL051958 (2012).
27
Koven, C. D.
et al.
The effect of vertically-resolved soil biogeochemistry and alternate soil C and N models
on C dynamics of CLM4.
Biogeosciences
10
, 7109-7131, doi:10.5194/bg-10-7109-2013 (2013).
28
Koven, C. D., Lawrence, D. M. & Riley, W. J. Permafrost carbon
−
climate feedback is sensitive to deep soil
carbon decomposability but not deep soil nitrogen dynamics.
Proceedings of the National Academy of
320
Sciences
112
, 3752--3757 (2015).
29
Wang, Y. P.
et al.
Responses of two nonlinear microbial models to warming and increased carbon input.
Biogeosciences
13
, 887--902, doi:10.5194/bg-13-887-2016 (2016).
30
Sistla, S. A.
et al.
Long-term warming restructures Arctic tundra without changing net soil carbon storage.
Nature
497
, 615--618 (2013).
325
Methods
A central goal of this paper is to use observed spatial gradients in soil carbon turnover, in
particular the relationship between apparent turnover and temperature, to develop a benchmark
for dynamical models that are used to make projections of soil carbon storage in response to
climate change. We recognize that benchmarks derived from spatial gradients are insufficient to
constrain transient responses, which is also why the direct use of spatial gradients to extrapolate
forward in time (so-called space-for-time substitutions) is not possible. However, transient
models must also make predictions about steady-state differences in soil carbon turnover across
spatial climate gradients, which reflect a long-term climatological temperature sensitivity.
Because such gradients are observable, we seek to use this information as a test of the dynamical
models. We contend that this inverse “time-for-space” substitution can serve as a global,
observationally-derived benchmark to ask where and whether the underlying processes
represented in ESMs are consistent with the observations. We note that the transient dynamics
predicted by a given model may differ from that model’s steady-state spatial-gradient
predictions, due to slow processes such as mineral associations that are important in governing
spatial gradients, and whose transient effects may be highly timescale-dependent. We also note
that the benchmarking approach here is fundamentally a consistency test, which is a necessary
but insufficient constraint on model dynamics, that does not imply that any given model
formulation provides a unique solution. We thus advocate that future experimental work also
focus on understanding the transient soil carbon dynamics across the world’s climate regimes,
and in particular focus on monitoring both the productivity and turnover responses to climate
change.
Soil carbon stocks are estimated by combining the Harmonized World Soilds Database
(HWSD)
1
and Northern Circumpolar Soil Carbon Database (NCSCD)
2
soil carbon maps, using
NCSCD where overlap occurs, and estimating productivity via the MODIS net primary
productivity (NPP) product
3
(fig. 1a & 1b). Mean annual air temperatures (MAAT) are
estimated from the CRU dataset
4
. Both NPP and soil carbon show relationships with
temperature, however the log(NPP)-temperature relationship shows a more continuous slope (fig.
1a) than the log(soil carbon)-temperature relationship (fig. 1b), which shows a clear difference
between temperate-tropical and a cold-climate soils. Inferred
τ
is plotted on a log scale because
we expect it to have a roughly exponential relationship with temperature, following that of
19
respiration, and
spans a range
of about 2 orders of magnitude (fig. 1c & 1d). The full global
dataset (Fig. 1c and S1a) shows a main set of points that span a curve of minimum turnover
times for a given temperature, with a tail of points extending above this main population. Water
also exerts a strong control on soil carbon turnover, with reduced decomposition in soils that are
either dry or saturated. Since our main goal here is to focus on temperature controls to
decomposition, we seek to separate and mask out those soils that are most strongly affected by
having either too much or too little water.
To isolate soils whose decomposition is limited by saturation, we use the thematic
6
classification in the NCSCD and HWSD databases to exclude those gridcells that have an areal
7
coverage of more than 50% peat soils, defined as the Histosol soil order or the Histel suborder of
8
Gelisols (permafrost soils) (fig. S1b). To identify soils whose decomposition is limited by
9
aridity, we derive the rainfall regime in each gridcell using the GPCC dataset
5
(fig. 1c). We
10
note that, while long turnover times are associated with both cold climates or very dry climates,
11
only the former actually have large organic carbon stocks (fig. S1c). Furthermore, in cold
12
climates, soils may receive relatively little rainfall but still have abundant moisture because
13
evaporative demand in these climates is low. Thus, we mask out only those soils where demand
14
exceeds supply and soils are dry. We calculate the demand as a potential evapotranspiration
15
(PET) using the MODIS PET product
6
, and show the inferred
τ
-temperature relationship as
16
controlled by precipitation minus PET (fig. S1d). This supports the expectation that long-
17
turnover low latitude soils are predominantly found in dry climates with a strong moisture
18
deficit. Finally, we exclude all points where this moisture deficit falls below a threshold level (-
19
1000 mm y
-1
), to arrive at the filtered dataset that we use to define the regression curve in fig. 1d.
20
Once the points where
τ
is not primarily a function of temperature are removed, the
21
relationship between
τ
and temperature becomes clear. The relationship in log(
τ
) versus
22
temperature is negative, with positive curvature. Although a linear regression shows high
23
significance (r
2
=0.62), it also leaves a residual error with systematic structure. We thus reject a
24
linear model relating log(
τ
) to temperature, and find that a quadratic model (fig. 1d, and table S1)
25
fits the relationship well, with r
2
=0.68 and little systematic residual bias. Thus, the relationship
26
between inferred
τ
and temperature is stronger than the hypothesized exponential in temperature
27
over the range of climate conditions.
28
20
We calculate the curve in figure 2 based on the derivative of the central relationship in
29
figure 1b.
Q
10
is an exponential notation that is traditionally defined relative to the instantaneous
30
decomposition parameter
k
as:
31
𝑘
𝑇
=
𝑘
𝑇
!"#
𝑄
!"
!
!
!
!"#
!"
32
The Q
10
parameter can therefore be calculated as:
33
34
𝑄
!"
=
𝑘
𝑇
𝑘
!"#
!"
!
!
!
!"#
35
Since, by definition:
36
𝑘
=
1
𝜏
37
we can redefine Q
10
in terms of
τ
:
38
𝑄
!"
=
𝜏
!"#
𝜏
𝑇
!"
!
!
!
!"#
39
and calculate the Q
10
at any point along the curve as the derivative of log(
τ
) with respect to
40
temperature
via:
41
𝑄
!"
=
10
!
!"
!
!"#
(
!
)
!
!
42
By choosing a polynomial regression in figure 1 of the form:
43
44
log
𝜏
=
𝑎
𝑇
!
+
𝑏𝑇
+
𝑐
45
these equations combine as:
46
𝑄
!"
=
10
!
!"
(
!
!"
!
!
)
47
which is shown in figure 2.
48
21
49
We identify in figure 3 a framework to understand the implications of alternate ways of
50
representing the temperature sensitivity of soil carbon decomposition on resulting spatial
51
patterns. For results shown in figure 3, we diagnose daily decomposition rate values (
k)
as a
52
function of daily soil temperatures, which are taken from a land model, CLM4.5,
driven by bias-
53
corrected reanalysis meteorological data (CRU-NCEP, available at
54
http://dods.ipsl.jussieu.fr/igcmg/IGCM/BC/OOL/OL/CRU-NCEP/ ). We then calculate the
55
equilibrium
τ
as the reciprocal of time-averaged
k
values, and plot log(
τ
) as a function of the
56
driving MAAT. In each case, the
k
, and therefore
τ
, values are relative to an arbitrary offset to
57
align with the central estimate from figure 1d at 15 ˚C.
58
In the simplest case, we use a fixed Q
10
value (Q
10
=1.5) across all temperatures, and we
59
diagnose
k
and
τ
using near-surface (10 cm depth) soil temperatures (fig. 3a). Log(
τ
) values show
60
the expected linear relationship against air temperature with only two emergent features. First, all
61
τ
values fall slightly below the line relating MAAT to the log(
τ
) as diagnosed from that MAAT.
62
This offset is a result of the seasonal cycles in soil temperatures, which give a time-average
τ
that
63
is less than
τ
calculated from time-averaged temperatures,
because the shape of the relationship
64
k
=f(T) has a positive curvature
.
Secondly there is a step offset at the transition from temperate to
65
cold climates attributable to the insulating effect of snow, which, where seasonally present,
66
elevates mean soil temperatures above mean air temperatures
7
.
67
As a next step in complexity, we consider a temperature-dependent Arrhenius-like
68
relationship
8
:
69
𝑘
=
𝑘
!"#
𝑒
!"#
.
!"
!
!"
.
!"
!
!
!
!
!!"
.
!"
as a potentially more realistic model than a fixed
Q
10
model across all temperatures (Fig. 3b).
70
Such an approach is a better match to the climatological temperature sensitivity in the cold-
71
climate regime, but fares worse as compared to the observation-based relationship
in the warm-
72
climate regime. We note as well that, particularly in cold climates, the modeled annual-mean
τ
73
values fall well below the blue line that represents the
τ
evaluated from the annual mean
74
temperature. This suggests that such chemical-kinetic factors, on their own, are not responsible
75
for the observed climatological cold-climate sensitivity.
76
A third step in complexity is to explicitly consider the role of soil freezing. Freezing is a
powerful inhibitor of decomposition, and thus we can define the simplest freezing model as one
22
where Q
10
is fixed (Q
10
= 1.5, as in fig 3a) for thawed soils and ceases in frozen soils (
k
= 0 when
77
T
soil
< 0
˚C
). This is a simplifying assumption; in natural ecosystems, limited decomposition
78
occurs even in frozen soils, but the log(
τ
) versus MAAT
resulting
from the no-frozen-respiration
79
case
(fig. 3c)
is a useful approximation
. The result is a slight positive curvature to the
80
relationship
that compensates for the snow insulation offset. But the model still fails to match the
81
large increase in temperature sensitivity seen in the observations.
82
The fourth piece of complexity in figure 3 is to consider the vertical variation in soil
83
temperatures in addition to the seasonal variation. To do this, we diagnose
k
values at each model
84
soil level (9 levels) across the 0-1m soil depth interval, and then calculate
τ
as the reciprocal of
85
the mean
k
value averaged both in time and over the depth interval (fig. 3d). This gives a curve
86
of log(
τ
) versus MAAT that approximates the observed one. We note that this behavior arises
87
from the step-function-like behavior in
k
=f(T) due to freezing, which gives particularly high
88
sensitivity to small temperature changes around the freezing point.
89
Comparing the quadratic regression coefficients for the simplified models shown in figure 3
90
(table S1 and figure S4), the value of the quadratic parameter increases with each level of
91
complexity. While each of the simplified models shown in figure 3 underestimates the curvature
92
as seen in the observations, only the fourth model begins to approach the observation-derived
93
curvature in the temperature sensitivity. This emphasizes the importance of representing how
94
freeze/thaw state scales both temporally and vertically in governing soil carbon decomposition.
95
For the models shown in figure 4, inferred
τ
is calculated as the ratio of carbon stocks to
96
NPP for each model over the recent historical period, in order to match the observation-based
97
results discussed above, and individual gridcells are colored by a given ESM’s precipitation as in
98
figure 1c, (fig. 4a-f). The CMIP5 models used here are listed in table S2. For CMIP5 models
99
(fig. 4a-4f) and MIMICS (fig. 4i), we diagnose
τ
using total soil carbon stocks, since the models
100
do not provide depth-resolved soil carbon output. For CLM4.5 (fig. 4g-4h), we use stocks to 1m
101
as the carbon stocks in the
inferred
τ
calculation. We note that the lower boundary for soil
102
carbon both in the CMIP5 protocol and in any soil model that does not explicitly resolve depth
103
is not clearly defined, and that there may be errors in correspondence between our comparison of
104
model output and 0-1m integrated soil carbon stocks which could be reduced by models
105
explicitly defining the depth integral over which their stocks correspond; nonetheless we believe
106
this represents the best benchmark on the models.
107
23
To create a metric of how well each model captures the pattern in the observed relationship
108
(Table S1), we use a binned RMSE scoring approach following
9
.
We first filter all points in each
109
model using the same P-PET threshold using model-predicted P and MODIS PET as in the
110
observations (masked data shown fig. S5). Next, we defined bins of 1˚C MAAT over the
111
interval of -15.5C to 28.5C, and for each model, took the mean value of log(
τ
) within each bin to
112
calculate an error score,
e,
as the RMSE for each bin
i
of the model prediction
p
i
relative to the
113
observational estimate
o
i
of the central estimate from fig. 1 evaluated at the bin center:
114
𝑒
=
𝑝
!
!
−
𝑜
!
!
In addition to the RMSE score, we also report in table S1 quadratic fit parameters between
115
MAAT and log(
τ
) for each of the models (regression curves shown in fig. S5), as well as the data
116
shown in figure 1D. For the regression, we mask arid areas using the same P minus PET criteria.
117
In all cases, the quadratic coefficient
a
derived from the regression is less than that from the
118
observations, and in some cases has the wrong sign; the three models that most closely approach
119
the magnitude of the quadratic term are GFDL-ESM2G, MIMICS and CLM4.5. We also report
120
the residual variance in log(
τ
) after subtracting the regression relationship, for the models and
121
data. Each of these metrics provides different benchmarking constraints on the dynamical
122
models; we leave them as separate constraints rather than merging into a single benchmark here.
123
Each of the models shown in figure 4 has unique characteristics, and some better
approximate the observed curve than others. Four models (CCSM4
10
, MPI-ESM
11
, HadGEM2
1
,
and IPSL-CM5A-LR
12
) show temperature responses that approximate the simple Q
10
relationship. The GFDL-ESM2G
13
model captures the range in turnover in cold-climates well.
This model calculates
k
values based on a function of the mean temperature over a root-profile
weighted depth interval (E. Shevliakova, personal communication), so includes more
information about deeper soil climate than models that diagnose
k
values based only on near-
surface information; however we also note that much of the high cold-climate sensitivity in this
models appears to be due to a large offset that occurs at around
0
˚C, which may also arise from
the anomalously cold high-latitude soil temperatures in that model
14,15
. MIROC-ESM
16
also
captures some aspects of the cold-climate sensitivity, particularly over the range of temperatures
5˚C to -5˚C, which is consistent with its use of the Lloyd-Taylor equation for its temperature
sensitivity
16
; however, we note that its temperature sensitivity reverses direction in the -5˚C to -
24
10˚C range, and decreases below that. With the exception of MPI-ESM and MIMICS, all models
124
appear to underestimate the sensitivity of decomposition to dry conditions, as indicated by the
125
lack of long-turnover soils in arid conditions that is seen in fig. 1c. CLM4.5 tends to
126
underestimate turnover in the region near MAAT ~-5˚C, which indicates an underestimate on the
127
limitation of decomposition in warm permafrost conditions in that model.
128
All data and analysis scripts required to generate figures in this manuscript are available
129
online at: http://portal.nersc.gov/archive/home/c/cdkoven/www/soil_tau_temp
130
131
132
133
134
25
References (Methods)
135
136
1
FAO, IIASA, ISRIC, ISSCAS & JRC. Harmonized World Soil Database (version 1.2). (FAO and IIASA,
137
Rome, Italy and Laxenburg, Austria, 2012).
138
2
Hugelius, G.
et al.
The Northern Circumpolar Soil Carbon Database: spatially distributed datasets of soil
139
coverage and soil carbon storage in the northern permafrost regions.
Earth Syst. Sci. Data
5
, 3--13,
140
doi:10.5194/essd-5-3-2013 (2013).
141
3
Zhao, M., Heinsch, F., Nemani, R. & Running, S. Improvements of the MODIS terrestrial gross and net
142
primary production global data set.
Remote Sensing of Environment
95
, 164--176, doi:DOI
143
10.1016/j.rse.2004.12.011 (2005).
144
4
Jones, P. D. H., I. Climatic Research Unit (CRU) time-series datasets of variations in climate with
145
variations in other phenomena. (University of East Anglia Climatic Research Unit, British Atmospheric
146
Data Centre, date of citation, 2008).
147
5
Schneider, U.
et al.
FD\_M\_V6\_050
, doi:10.5676/DWD_GPCC/FD_M_V6_050 (2011).
148
6
Mu, Q., Zhao, M. & Running, S. W. Improvements to a MODIS global terrestrial evapotranspiration
149
algorithm.
Remote Sensing of Environment
115
, 1781--1800,
150
doi:
http://dx.doi.org/10.1016/j.rse.2011.02.019
(2011).
151
7
Slater, A.
et al.
The representation of snow in land surface schemes: Results from PILPS 2(d).
Journal of
152
Hydrometeorology
2
, 7-25 (2001).
153
8
LLOYD, J. & TAYLOR, J. ON THE TEMPERATURE-DEPENDENCE OF SOIL RESPIRATION.
154
Functional Ecology
8
, 315-323 (1994).
155
9
Slater, A. G., Lawrence, D. M. & Koven, C. D. Process-level model evaluation: A Snow and Heat Transfer
156
Metric.
The Cryosphere Discuss.
2016
, 1-16, doi:10.5194/tc-2016-258 (2016).
157
10
Oleson, K. W.
et al.
Technical Description of version 4.0 of the Community Land Model (CLM).
158
(NATIONAL CENTER FOR ATMOSPHERIC RESEARCH, P. O. Box 3000 BOULDER, COLORADO
159
80307-3000, 2010).
160
11
Raddatz, T. J.
et al.
Will the tropical land biosphere dominate the climate-carbon cycle feedback during the
161
twenty-first century?
Climate Dynamics
29
, 565-574, doi:DOI 10.1007/s00382-007-0247-8 (2007).
162
12
Krinner, G.
et al.
A dynamic global vegetation model for studies of the coupled atmosphere-biosphere
163
system.
Global Biogeochem. Cy.
19
, GB1015, doi:10.1029/2003GB002199 (2005).
164
13
Milly, P. C. D.
et al.
An Enhanced Model of Land Water and Energy for Global Hydrologic and Earth-
165
System Studies.
Journal of Hydrometeorology
15
, 1739--1761, doi:10.1175/JHM-D-13-0162.1 (2014).
166
14
Slater, A. G. & Lawrence, D. M. Diagnosing Present and Future Permafrost from Climate Models.
Journal
167
of Climate
26
, 5608--5623, doi:10.1175/JCLI-D-12-00341.1 (2013).
168
15
Koven, C. D., Riley, W. J. & Stern, A. Analysis of permafrost thermal dynamics and response to climate
169
change in the CMIP5 Earth System Models.
Journal of Climate
, doi:10.1175/JCLI-D-12-00228.1 (2012).
170
16
Sato, H., Itoh, A. & Kohyama, T. SEIB--DGVM: A new Dynamic Global Vegetation Model using a
171
spatially explicit individual-based approach.
Ecological Modelling
200
, 279--307,
172
doi:
http://dx.doi.org/10.1016/j.ecolmodel.2006.09.006
(2007).
173