Georgiou et al. (2024). Emergent temperature sensitivity of soil organic carbon driven by mineral associations

Nature Geoscience

| Volume 17 | March

2024 | 205–212

205

nature geoscience

https://doi.org/10.1038/s41561-024-01384-7

Article

Emergent temperature sensitivity of

soil organic carbon driven by mineral

associations

Katerina Georgiou

, Charles D. Koven

, William R. Wieder

3,4

Melannie D. Hartman

, William J. Riley

, Jennifer Pett-Ridge

1,6

Nicholas J. Bouskill

, Rose Z. Abramoff

2,7

, Eric W. Slessarev

1,8

Anders Ahlström

, William J. Parton

, Adam F. A. Pellegrini

Derek Pierson

, Benjamin N. Sulman

, Qing Zhu

Robert B. Jackson

13,14

Soil organic matter decomposition and its interactions with climate depend
on whether the organic matter is associated with soil minerals. However, data
limitations have hindered global-scale analyses of mineral-associated and
particulate soil organic carbon pools and their benchmarking in Earth system
models used to estimate carbon cycle–climate feedbacks. Here we analyse
observationally derived global estimates of soil carbon pools to quantify
their relative proportions and compute their climatological temperature
sensitivities as the decline in carbon with increasing temperature. We find
that the climatological temperature sensitivity of particulate carbon is
on average 28% higher than that of mineral-associated carbon, and up to
53% higher in cool climates. Moreover, the distribution of carbon between
these underlying soil carbon pools drives the emergent climatological
temperature sensitivity of bulk soil carbon stocks. However, global models
vary widely in their predictions of soil carbon pool distributions. We show
that the global proportion of model pools that are conceptually similar
to mineral-protected carbon ranges from 16 to 85% across Earth system
models from the Coupled Model Intercomparison Project Phase 6 and offline
land models, with implications for bulk soil carbon ages and ecosystem
responsiveness. To improve projections of carbon cycle–climate feedbacks,
it is imperative to assess underlying soil carbon pools to accurately predict
the distribution and vulnerability of soil carbon.

Soil carbon–climate feedbacks represent a major uncertainty in the
response of the terrestrial biosphere to climate change

. This uncer-

tainty stems, in part, from poorly constrained relationships between
the temperature sensitivities of the decomposition and stabilization
processes of soil organic carbon (C), which are key parameters in soil
biogeochemical models

. Whereas many empirical and modelling

studies have explored the effects of temperature on bulk soil C stocks
and turnover rates at ecosystem to global scales

–

, they largely ignore

heterogeneity in soil organic matter. This is a key oversight given that
underlying soil organic matter fractions or ‘pools’ (for example, par-
ticulate or mineral-associated organic matter) can exhibit different
responses to warming in laboratory- and field-based manipulations

–

Received: 23 October 2022

Accepted: 23 January 2024

Published online: 20 February 2024

Check for updates

A full list of affiliations appears at the end of the paper. e-mail:

georgiou1@llnl.gov

Page 2

Nature Geoscience

| Volume 17 | March

2024 | 205–212

206

Article

https://doi.org/10.1038/s41561-024-01384-7

and particulate soil C (derived from an observational synthesis of
soil fractions

) to quantify the distribution of soil organic C between

these two pools and their respective climatological temperature sen-
sitivities at the global scale. We then evaluate the distributions and
temperature sensitivities of protected (mineral-associated-like) and
unprotected (particulate-like) soil C stocks across an ensemble of
12 global models—namely, nine Earth system models (ESMs) from
the Coupled Model Intercomparison Project Phase 6 (CMIP6)

and

three offline land models

, including two microbially explicit soil C

models (Supplementary Tables 1–3). Although there are important
considerations regarding the interpretability of modelled pools and
measured fractions (see ‘Model pool interpretability’ in the Methods),
pool-specific benchmarks offer a critical opportunity to refine model
formulations and reduce uncertainties. We argue that resolving the
distribution and temperature sensitivities of underlying soil C pools is
imperative for developing a predictive understanding of the emergent
temperature dependence of bulk soil C decomposition and constrain-
ing global biogeochemical models.

Climatological temperature sensitivity of soil

carbon pools

Soil organic C stocks are known to broadly decrease with increasing
temperature across climate gradients

. However, the magnitude

of this climatological temperature sensitivity can vary across soils

Indeed, this relationship appears to be modulated by clay and silt con-
tent, where fine-textured soils (that is, soils containing higher amounts
of clay and silt minerals) have a lower climatological temperature
sensitivity compared with coarse-textured sandy soils (inferred from
the relative strengths of relationships in Fig.

). These contrasting

temperature sensitivities can be attributed in part to the mineralogi-
cal capacity of each soil

, but more specifically to its effect on the

distribution of soil C among mineral-associated and particulate C pools.
Since mineral-associated C often has a lower temperature sensitivity
than particulate C

, the proportion of C that is mineral-associated

can play a dominant role in driving the emergent bulk soil temperature
sensitivity (see conceptual schematic in Supplementary Fig. 1).

Globally, our results indicate that the proportion of soil C that is

mineral-associated (hereafter termed ‘protected’ for consistency with
model pools

) increases with increasing mean temperature (Fig.

; a

trend that is not captured in most global models (Supplementary Figs. 2
and 3)). Furthermore, fine-textured soils have greater proportions
of protected C than do coarse-textured soils across all temperatures.
We therefore expect fine-textured soils to exhibit a lower climato-
logical temperature sensitivity compared with coarse-textured sandy

Understanding and quantifying the temperature sensitivity of these
underlying pools at the global scale is thus critical for accurately pre-
dicting emergent feedbacks and soil organic C vulnerabilities to a
changing climate.

Warmer soil temperatures accelerate the rate of organic matter

decomposition, often leading to losses in bulk soil C stocks

. How-

ever, the potential decomposability of organic matter in a given soil
can depend on its association with minerals

. Mineral-associated

organic matter consists of small plant- or microbial-derived monomers
and biopolymers that are chemically or physically bound to clay- and
silt-sized minerals, whereas particulate organic matter is largely com-
prised of coarser, partially decomposed fragments of plant material

Mineral associations can limit access of microbial constituents and
their enzymes to organic matter, and as a result, mineral-associated C is
often older than particulate C

–

. In soil biogeochemical models, these

pools can be represented explicitly

–

or implicitly

–

(see ‘Model pool

interpretability’ in the Methods), with implications for their transient
dynamics and long-term response to novel conditions. For consist-
ency across the data and models, we refer to mineral-associated C
(and slowest-cycling pools) as ‘protected’ and particulate C (and all
other model pools) as ‘unprotected’, while acknowledging that the bio-
availability and persistence of each pool is more complex in practice

Importantly, much uncertainty remains in the responses of protected
and unprotected soil C to warming across broad spatial scales

, hinder-

ing model benchmarking efforts.

Evaluating the transient dynamics of soil biogeochemical models

to warming remains a challenge, in part because observations of soil C
at multi-decadal to centennial timescales are limited

. Laboratory

incubations and field studies demonstrate that decomposition rates
increase, and C stocks often decrease, with warming

. However, at

multi-decadal timescales, acclimatory responses in microbial activity
may limit the effect of warming on decomposition rates

. Further-

more, many experiments focus on the response of bulk soil C stocks,
and neither measure nor report the relative responses of underlying soil
C fractions. Given these data limitations, an alternative approach is to
benchmark global models using spatial gradients

—although this

approach is insufficient for constraining transient dynamics (that is,
space-for-time assumptions)

, climate gradients do reflect long-term

climatological temperature sensitivities that should, in principle, be
captured by soil C models.

Here we assess the climatological temperature sensitivity of bulk

soil C stocks, and the role of mineral-associated and particulate soil
C pools in driving this emergent property in data and models. We
first leverage a globally gridded data product of mineral-associated

Clay + silt (%)

100

MAT (°C)

Tot

al soil C stock

s (k

gC m

−2

)

Fine-textur

Coar

se-te

xtured

100

MAT (°C)

opor

tion of pr

otected C (%)

Fine-textur

Coarse-te

xtured

Fig. 1 | Global distribution of observationally derived soil C stocks as a
function of climatological temperature and clay and silt minerals. a,b, Total
soil C stocks (

a) and the proportion of total soil C that is protected (mineral-

associated) (

b) as a function of the mean annual temperature (MAT) globally.

Each gridcell is coloured by the percentage of clay and silt minerals, and best-fit
trends are depicted for fine- and coarse-textured soils; here, fine-textured soils
were classified as those with >70% clay + silt content and coarse-textured soils
with <20% clay + silt content.

Page 3

Nature Geoscience

| Volume 17 | March

2024 | 205–212

207

Article

https://doi.org/10.1038/s41561-024-01384-7

soils where higher unprotected C stocks drive a higher temperature
sensitivity. The climatological temperature sensitivity of bulk soil C
emerges from the underlying pool distributions and their respective
temperature dependencies.

To assess the climatological temperature sensitivity of each soil

C pool, we calculated the proportional decline in C stocks for every
10 °C increase in mean temperature while controlling for potential
confounding factors that include primary productivity, precipitation,
and clay and silt content (Methods). Our approach builds upon previous
studies

and provides a benchmarking metric analogous to a climato-

logical Q

(that is, the increase in process rates for each 10 °C increase in

temperature) but for C stocks instead of respiration or turnover rates

We find that, at the global scale, unprotected C stocks have a 28% (95%
confidence interval (CI): 26, 30%) greater climatological temperature
sensitivity than protected C stocks (Fig.

). This pattern is especially

pronounced in cool regions (<15 °C) where unprotected C is 53% (46,
60%) more temperature sensitive than protected C, compared with
warm regions (≥15 °C) where the relative difference is only 15% (11, 20%).
Both soil C pools show weak climatological temperature sensitivities in
warm regions (Fig.

). This latter result suggests that decomposition

in cool regions is more temperature-limited, and that warm regions do
not experience as stark a decline in soil C for a given increase in mean
temperature. This difference is supported by gradient studies in tropi-
cal soils

, and is consistent with incubation studies showing that

warm-adapted microbial communities typically show lower temperature
sensitivities (namely, Q

) and higher microbial growth efficiencies

Furthermore, warm regions contain higher proportions of protected C
(Fig.

), such that the emergent bulk soil C temperature sensitivity more

closely mirrors that of the protected C pool. Therefore, the proportion
of protected C—in addition to its climatological temperature sensitiv-
ity—may serve as an important benchmark for evaluating global models.

Global patterns in soil carbon pools across data

and models

Soil C models vary in their underlying assumptions, mathematical
representations and parameterizations. However, several studies have
proposed that the slowest-cycling (‘passive’) pool of soil C in many
ESMs corresponds broadly to mineral-associated (protected) C

–

This reflects a common conceptualization of donor-pool sources,

mineralogical controls, turnover times and stoichiometry of protected
C that are used to approximate the formation and characteristics of
mineral-associated organic matter (see ‘Model pool interpretability’
in the Methods). We thus adopt this convention here, while recogniz-
ing that there are still inherent uncertainties in the interpretability
of modelled pools and operationally defined fractions. We posit that
benchmarks based on relative pool distributions, in tandem with radio-
carbon ages

, will be critical for refining future model parameteriza-

tions and constraining transient dynamics.

Global soil C models are typically benchmarked on the basis of bulk

soil C metrics

. However, even models that predict similar bulk soil C

stocks can vary widely in their underlying pool distributions (Supple-
mentary Fig. 4). Whereas bulk soil C stocks in non-permafrost mineral
soils vary 2.5-fold across the CMIP6 ESMs and offline land models in this
study, protected soil C stocks vary nearly sevenfold. As a result, the pro-
portion of protected C ranges from 16 to 85% across the global models
(Fig.

and Supplementary Table 4). This large spread in model predic-

tions presents a critical opportunity to apply data-driven constraints
and reduce uncertainty. We find that about half of the analysed ESMs
underestimate the proportion of protected C in soils globally compared
with the observationally derived data product, and consequently,
contain too much C in faster-cycling, unprotected C pools. This result
is especially true in non-permafrost mineral soils, where several of the
global models underestimate the proportion of protected C in both
cool and warm regions (Supplementary Fig. 5 and Supplementary Table
5). In particular, six models (CESM2, CNRM-ESM2-1, E3SM-1-1-ECA,
IPSL-CM6A-LR, NorESM2 and the microbially explicit model CORPSE)
roughly capture the proportion of protected C, whereas the other six
models (ACCESS-ESM1-5, BCC-CSM2-MR, MIROC-ES2L, MRI-ESM2-0,
CASA-CNP and the microbially explicit model MIMICS) considerably
underestimate this proportion (Fig.

). However, we note that all

of the models in this study are amenable to matching pool-specific
constraints in future work—indeed, both first-order and microbially
explicit formulations with refined parameterizations have been able
to better match observations in past studies

–

. Nevertheless, our

results generally suggest that soil C is skewed towards faster-cycling
C pools in several global models, consistent with radiocarbon-based
benchmarking studies that have shown that current ESMs underesti-
mate the mean age of soil C across latitudes

0.1

1.0

10.0

100.0

MAT (°C)

otected soil C stock

s (k

gC m

−2

)

0.1

1.0

10.0

100.0

MAT (°C)

npr

otected soil C stock

s (k

gC m

−2

)

1.0

1.5

2.0

2.5

3.0

Cool

All

arm

opor

tional decline in C (10 °C)

Unprotected C
Protected C

Fig. 2 | Global distribution of observationally derived protected and
unprotected soil C stocks as a function of climatological temperature.
a,b, Protected (mineral-associated) (a) and unprotected (particulate) (b) soil
C stocks as a function of MAT in regions above 0 °C globally. Linear regressions
are shown for cool (<15 °C) and warm (≥15 °C) regions (solid lines) and across
all temperatures above 0 °C (dashed lines).

c, Climatological temperature

sensitivity calculated as the proportional decline of each soil C pool for a 10 °C
increase in MAT (analogous to a climatological Q

, but for C stocks) while

accounting for potential confounding variables (see ‘Data analysis’ in the
Methods). Points depict the proportional decline of each soil C pool (derived
from the regression slopes; n = 32,528 gridcells for non-permafrost mineral
soils above 0 °C globally, with n = 13,410 in cool regions and n = 19,118 in warm
regions), and error bars denote the 95% CIs (Supplementary Table 6). Higher
values (>1) indicate greater decreases in C with increasing climatological
temperature.

Page 4

Nature Geoscience

| Volume 17 | March

2024 | 205–212

208

Article

https://doi.org/10.1038/s41561-024-01384-7

At higher latitudes, the two microbially explicit models depict a

lower proportion of protected C (Supplementary Fig. 4), consistent
with recent studies

. The observationally derived data product

also supports this latitudinal trend (Fig.

). By contrast, several CMIP6

models show a larger proportion of protected C at higher latitudes
(particularly above 0 °C; Supplementary Fig. 3), which suggests that
the interpretability of the protected pool as mineral-associated may
break down in these regions within certain CMIP6 models. For instance,
in models with vertically resolved soil biogeochemistry (CESM2, E3SM-
1-1-ECA and NorESM2) these larger proportions of protected C may be
caused by the relative timescales of decomposition and vertical mixing
processes. This is because high-latitude active-layer C is transported
to permafrost layers via a slow constant diffusion rate representing
cryoturbation

, which consequently transports a greater fraction

of slower-cycling C relative to fast-cycling C pools. Whereas soil C is
indeed old in these high latitudes

, its decomposition is limited more

by low temperatures that reduce the availability of liquid water than by
association with minerals

. Given that temperature and hydrology can

change at decadal timescales, representing long-lived organic matter at
high latitudes as an inherently slow-cycling pool (as opposed to unpro-
tected, particulate forms that can be more sensitive to disturbance

)

has important implications for the responsiveness of vast permafrost
soil C pools to projected climate change (Supplementary Figs. 6–8).

Temperature sensitivity of soil carbon in global

models

The climatological temperature sensitivity of bulk soil C is an emergent
property that ultimately depends on the distribution of C within under-
lying pools and their respective temperature sensitivities. Whereas

some models may appear to capture the observationally derived clima-
tological temperature sensitivity of bulk soil C stocks, this agreement
does not always occur for the right reasons and may break down when
considering the contribution of underlying protected and unprotected
soil C pools. For example, six of the CMIP6 ESMs (ACCESS-ESM1-5,
CESM2, CNRM-ESM2-1, E3SM-1-1-ECA, IPSL-CM6A-LR and NorESM2)
and two microbially explicit models (MIMICS and CORPSE) roughly
capture the temperature sensitivity of bulk soil C stocks (Fig.

and

Supplementary Table 6). However, these same models show divergent
distributions and climatological temperature sensitivities of underly-
ing protected and unprotected soil C pools (Fig.

). Only a few models

roughly capture the climatological temperature sensitivity of protected
C (CESM2, CNRM-ESM2-1, IPSL-CM6A-LR, NorESM2 and CORPSE) and
even fewer of unprotected C (ACCESS-ESM1-5, MIMICS and CORPSE).

The observationally derived data product broadly suggests that

protected C has a lower climatological temperature sensitivity than
unprotected C (Fig.

). The models generally agree with this trend,

although the difference in climatological temperature sensitivities
of pools is less pronounced in the CMIP6 models than in the microbi-
ally explicit models (Fig.

). This difference may arise because many

first-order soil models impose the same temperature-dependence
parameters across all pools

, and consequently the emergent cli-

matological temperature sensitivities are often similar between pools.
There are a few notable exceptions among the CMIP6 models, in which
unprotected C appears less temperature sensitive than protected C
(namely, MRI-ESM2-0 and to a lesser degree E3SM-1-1-ECA (Fig.

), as

well as CESM2, CNRM-ESM2-1, MIROC-ES2L and NorESM in cool regions
(Supplementary Fig. 9)). This could arise due to temperature-sensitivity
parameterizations based on a presumed greater molecular complexity

Dat

a pr

oduct

Latitude (deg)

opor

tion of pr

otected C (%)

60° N

–20

–40

100

30° N

30° S

60° S

120° W

Proportion of protected C (%)

100

Mineral-

associated
Passive
Physicochemically

protected

60° E 120° E

60° W 0°

120° W

60° E 120° E

60° W 0°

Dat

a pr

oduct

ACCES

S-ESM1-5

BCC-CSM2

-MR

CESM2

CNRM-ESM2

-1

E3SM-

1-E

-CM6A

-LR

MIR

OC-ES2L

MRI-ESM2

-0

NorESM2

-CNP

MIMICS

CORP

100

Proportion of protected C (%)

0°

60° N

30° N

30° S

60° S

0°

Multi-model me

Multi-model

mean

Data

product

Fig. 3 | Proportion of soil C in protected pools globally across data
and models. a–c, The percentage of total soil C that is protected in the
observationally derived data product and in the multi-model mean of the
ensemble of CMIP6 ESMs and offline land models, shown globally (

a,c) and as

latitudinal means (

b; solid lines). Shaded uncertainty ranges in b correspond to

90% prediction intervals for the data product

, and 5th and 95th quantiles of the

multi-model ensemble across each latitudinal band.

d, The global percentage of

total soil C that is protected across the individual ESMs and offline land models.

Protected C was estimated as the mineral-associated fraction for the data
product, physicochemically protected pools in MIMICS and CORPSE, and the
slowest-cycling passive pool in the CMIP6 ESMs and CASA-CNP (see ‘Model pool
interpretability’ in the Methods). Global means are given for non-permafrost
mineral soils, excluding tundra, deserts and peatlands (n = 39,620 gridcells for
0.5° resolution), and the data product uncertainty range corresponds to the 90%
prediction interval

(summarized in Supplementary Table 5).

Page 5

Nature Geoscience

| Volume 17 | March

2024 | 205–212

209

Article

https://doi.org/10.1038/s41561-024-01384-7

of protected C in select models

. Ultimately, the distinct climatologi-

cal temperature sensitivities of underlying soil C pools observed in the
data product and process-based models can have important implica-
tions for transient projections, especially with regard to C ages and
ecosystem responsiveness, and constitute an important benchmark
for evaluating model performance.

Since protected C makes up the majority of bulk soil C globally in

the data product (Fig.

), it also dominates the climatological tem-

perature sensitivity of bulk soil C stocks (Fig.

). However, many models

do not capture this pattern. In fact, unprotected C dominates the bulk
soil C stocks and climatological temperature sensitivity in a few models
(BCC-CSM2-MR, MRI-ESM2-0, CASA-CNP and MIMICS), whereas most
other ESMs have similar proportions and temperature sensitivities of
underlying C pools. Only a couple of global models have high enough
proportions of protected C to drive bulk climatological temperature
sensitivities that are dominated by the protected C pool (for example,
E3SM-1-1-ECA and CORPSE), as observed in the data product (Fig.

and

Supplementary Fig. 9).

Furthermore, climatological temperature sensitivities of bulk soil

C stocks and turnover can vary across climates, where temperature
dependencies are typically higher in cooler climates than in warmer
climates

. However, the temperature sensitivities of underlying

soil C pools within these climate regimes have not been explored in
global-scale models or observations. We find that both protected and
unprotected C pools show greater climatological temperature sensitiv-
ities in cool than warm climates, but interestingly, this contrast across
climates is more pronounced for unprotected C than for protected C,
and especially in the observationally derived data product compared
with the global models (Figs.

and

and Supplementary Figs. 9–11).

Overall, our analysis shows that global models tend to overesti-

mate the climatological temperature sensitivity of protected C stocks in

both cool and warm climates, and of both C pools in warm climates (Fig.

). This is especially true of MIROC-ES2L, CASA-CNP and MRI-ESM2-0,

which overestimate the climatological sensitivity of protected C by
more than a factor of two. However, the global models tend to under-
estimate the climatological temperature sensitivity of unprotected
C in cool climates (Fig.

). The inability of most ESMs to match these

regional patterns suggests that the transient response to warming may
also be biased, particularly in cool climates. Protected (slower-cycling)
C in ESMs may be overly sensitive to warming in cool regions, whereas
unprotected (faster-cycling) C may not be sensitive enough to warm-
ing, with implications on soil C dynamics and age distributions under
future scenarios.

Implications and future perspectives

We propose that the distributions and climatological temperature sen-
sitivities of underlying soil C pools are ecosystem properties that global
models should be expected to reproduce. Our results highlight that
even ESMs that predict similar bulk soil C stocks and climatological tem-
perature sensitivities can vary widely in their predictions of underlying
soil C pools, with implications for C ages and ecosystem responsiveness
under changing conditions. In particular, owing to an overestimation of
C in unprotected (faster-cycling) pools, several CMIP6 ESMs may be too
responsive to C inputs and may overestimate the effects of increased
primary production on potential soil C sequestration (for example, with
large increases in unprotected C; Supplementary Fig. 6); a finding that
agrees with independent radiocarbon-based studies

. This implies

that model estimates of the amount of C sequestration that may result
from climate change and CO

fertilization of terrestrial ecosystems

may be too strong. Furthermore, most ESMs tend to underestimate
the climatological temperature sensitivity of unprotected C in cool
climates, potentially leading to an underestimation of C losses from

0.5

1.0

1.5

2.0

2.5

3.0

3.5

Dat

a pr

oduct

ACCES

S−ESM1−5

BCC−

CSM2−MR

CESM2

CNRM−ESM2−1

E3SM−1−1−E

−CM6A

−LR

MIR

OC−ES2L

MRI−ESM2−0

NorESM2

CAS

−C

MIMICS

CORP

opor

tional decline in C (10 °C)

Unprotected C

Protected C

Total soil C

Proportion of C in each pool

14%

86%

C decrease

C increase

Thermal

independence

Fig. 4 | Climatological temperature sensitivity of total soil C and the
contribution of protected and unprotected pools. The climatological
temperature sensitivity of total soil C (denoted by a diamond) is an emergent
property of the climatological temperature sensitivity of protected (brown
circles) and unprotected (green circles) soil C pools, as well as the proportion
of C within each of these pools (denoted by the marker area). The proportional
decline in C for a 10 °C increase in MAT was calculated for each C pool, controlling

for potential confounding variables (see ‘Data analysis’ in the Methods). Results
are shown for the data product and an ensemble of ESMs and offline land models
across gridcells with MAT > 0 °C globally; 95% CIs are given in Supplementary
Table 6. Higher values (>1) indicate greater decreases in C with increasing
climatological temperature, lower values (<1) indicate increases in C and values
equal to 1 (grey dashed line) indicate thermal independence.

Page 6

Nature Geoscience

| Volume 17 | March

2024 | 205–212

210

Article

https://doi.org/10.1038/s41561-024-01384-7

these unprotected pools with warming, an underestimation of the
global carbon cycle–climate feedback and compounding projections of
higher productivity-driven C accumulation in these regions. Whereas
this bias may be partially offset in ESMs that overestimate the clima-
tological temperature sensitivity of protected C across climates, the
effect of this latter overestimate on projections of carbon cycle–cli-
mate feedbacks is probably weaker because these pools are intrinsi-
cally slow-cycling in first-order soil models and consequently change
very little on timescales relevant to climate policy (Supplementary
Figs. 6–8). Divergent pool-specific responses can ultimately lead to
differences in soil C composition and functioning. It is therefore impera-
tive that models benchmark not only bulk climatological sensitivities
but also those of underlying soil C pools, to confidently project soil C
dynamics into the future.

The CMIP protocol requests but does not require the reporting

of soil C pool distributions

, and, thus, this information is currently

available only for a subset of CMIP6 ESMs. Our results show that, with
recent advances in global estimates of protected and unprotected soil
C stocks

, redefining and benchmarking these underlying soil C pools

can serve as an important constraint for ESMs. We therefore urge that
all global models report the distribution of soil C among underlying
pools in CMIP Phase 7 and beyond, so that future studies may lever-
age this benchmark for model evaluation across climate gradients.
Although further study of transient responses is also imperative for
constraining dynamic models

, the static relationships presented

here serve as important benchmarks of long-term temperature effects
on the formation and decomposition of soil C pools.

Online content

Any methods, additional references, Nature Portfolio reporting sum-
maries, source data, extended data, supplementary information,
acknowledgements, peer review information; details of author con-
tributions and competing interests; and statements of data and code
availability are available at

https://doi.org/10.1038/s41561-024-01384-7

References

1. Friedlingstein, P. et al. Uncertainties in CMIP5 climate

projections due to carbon cycle feedbacks. J. Clim.

, 511–526

(2014).

2. Davidson, E. A. & Janssens, I. A. Temperature sensitivity of soil

carbon decomposition and feedbacks to climate change. Nature

440

, 165–173 (2006).

3. García-Palacios, P. et al. Evidence for large microbial-mediated

losses of soil carbon under anthropogenic warming. Nat. Rev.
Earth Environ.

, 507–517 (2021).

4. Koven, C. D., Hugelius, G., Lawrence, D. M. & Wieder, W. R. Higher

climatological temperature sensitivity of soil carbon in cold than
warm climates. Nat. Clim. Chang.

, 817–822 (2017).

5. Crowther, T. et al. Quantifying global soil C losses in response to

warming. Nature

104

, 104–108 (2016).

6. Todd-Brown, K., Zheng, B. & Crowther, T. W. Field-warmed soil

carbon changes imply high 21st-century modeling uncertainty.
Biogeosciences

, 3659–3671 (2018).

7. Van Gestel, N. et al. Predicting soil carbon loss with warming.

Nature

554

, E4–E5 (2018).

8. Hicks Pries, C. E., Castanha, C., Porras, R. & Torn, M. S. The

whole-soil carbon flux in response to warming. Science

355

1420–1423 (2017).

9. Rocci, K. S., Lavallee, J. M., Stewart, C. E. & Cotrufo, M. F. Soil

organic carbon response to global environmental change
depends on its distribution between mineral-associated and
particulate organic matter: a meta-analysis. Sci. Total Environ.

793

, 148569 (2021).

10. Conant, R. T. et al. Temperature and soil organic matter

decomposition rates – synthesis of current knowledge and a way
forward. Glob. Chang. Biol.

, 3392–3404 (2011).

11. Soong, J. L. et al. Five years of whole-soil warming led to loss

of subsoil carbon stocks and increased CO

efflux. Sci. Adv.

eabd1243 (2021).

ACC

BCC

CES

CNR

E3S

IPS

MIR

MRI

Nor

CAS

MIM

COR

ACC

BCC

CES

CNR

E3S

IPS

MIR

MRI

Nor

CAS

MIM

COR

Protected C

0.5

1.0

1.5

2.0

2.5

3.0

3.5

4.0

4.5

5.0

Cool climates

Warm climates

Pro

portional decline in C (10 °C)

Pro

portional decline in C (10 °C)

ACC

BCC

CES

CNR

E3S

IPS

MIR

MRI

Nor

CAS

MIM

COR

ACC

BCC

CES

CNR

E3S

IPS

MIR

MRI

Nor

CAS

MIM

COR

Unprotected C

0.5

1.0

1.5

2.0

2.5

3.0

3.5

4.0

4.5

5.0

Cool climates

Warm climates

Data product
Multi-model mean
Model output

ACC = ACCESS-ESM1-5

BCC = BCC-CSM2-MR

CES = CESM2

CNR = CNRM-ESM2-1

E3S = E3SM-1-1-ECA

IPS = IPSL-CM6A-LR

MIR = MIROC-ES2L

MRI = MRI-ESM2-0

Nor = NorESM2

CAS = CASA-CNP

MIM = MIMICS

COR = CORPSE

Fig. 5 | Climatological temperature sensitivity of protected and unprotected
soil C across temperature regimes. a,b, The proportional decline in protected
(

a) and unprotected (b) soil C stocks for a 10 °C increase in MAT, controlling

for potential confounding variables. Proportional declines are shown for cool
climates (0 °C < MAT < 15 °C; blue circles) and warm climates (MAT ≥ 15 °C; red
circles). Results are shown for the data product (black triangles), the multi-
model ensemble mean (asterisks) and each ESM and offline land model (solid

circles; labelled with abbreviations); 95% CIs are given in Supplementary Table 6.
Kernel density estimates are shown for each panel, illustrating the distribution
of climatological temperature sensitivities across the model ensemble for each
pool and temperature regime. Higher values (>1) indicate greater decreases in C
with increasing climatological temperature, lower values (<1) indicate increases
in C and values equal to 1 (grey dashed line) indicate thermal independence.

Page 7

Nature Geoscience

| Volume 17 | March

2024 | 205–212

211

Article

https://doi.org/10.1038/s41561-024-01384-7

12. Nottingham, A. T., Meir, P., Velasquez, E. & Turner, B. L. Soil carbon

loss by experimental warming in a tropical forest. Nature

584

234–237 (2020).

13. Pellegrini, A. F. A. et al. Low

‐

intensity frequent fires in coniferous

forests transform soil organic matter in ways that may offset
ecosystem carbon losses. Glob. Chang. Biol.

, 3810–3823 (2021).

14. Lugato, E., Lavallee, J. M., Haddix, M. L., Panagos, P. & Cotrufo, M. F.

Different climate sensitivity of particulate and mineral-associated
soil organic matter. Nat. Geosci.

, 295–300 (2021).

15. Sokol, N. W. et al. Global distribution, formation and fate of

mineral-associated soil organic matter under a changing climate:
a trait-based perspective. Funct. Ecol.

, 1411–1429 (2022).

16. Lavallee, J. M., Soong, J. L. & Cotrufo, M. F. Conceptualizing soil

organic matter into particulate and mineral-associated forms to
address global change in the 21st century. Glob. Chang. Biol.

261–273 (2020).

17. Kögel-Knabner, I. et al. Organo-mineral associations in temperate

soils: integrating biology, mineralogy, and organic matter
chemistry. J. Plant Nutr. Soil Sci.

171

, 61–82 (2008).

18. Heckman, K. et al. Beyond bulk: density fractions explain

heterogeneity in global soil carbon abundance and persistence.
Glob. Chang. Biol.

, 1178–1196 (2022).

19. Wieder, W. R., Allison, S. D., Davidson, E. A., Georgiou, K. &

Hararuk, O. Explicitly representing soil microbial processes in
Earth system models. Glob. Biogeochem. Cycles

, 1782–1800

(2015).

20. Sulman, B. N., Phillips, R. P., Oishi, A. C., Shevliakova, E. & Pacala,

S. W. Microbe-driven turnover offsets mineral-mediated storage
of soil carbon under elevated CO

. Nat. Clim. Chang.

, 1099–1102

(2014).

21. Ahrens, B., Braakhekke, M. C., Guggenberger, G., Schrumpf, M.

& Reichstein, M. Contribution of sorption, DOC transport and
microbial interactions to the

C age of a soil organic carbon

profile: Insights from a calibrated process model. Soil Biol.
Biochem.

, 390–402 (2015).

22. Sulman, B. N. et al. Multiple models and experiments underscore

large uncertainty in soil carbon dynamics. Biogeochemistry

141

109–123 (2018).

23. Wieder, W. R., Sulman, B. N., Hartman, M. D., Koven, C. D. &

Bradford, M. A. Arctic soil governs whether climate change drives
global losses or gains in soil carbon. Geophys. Res. Lett.

14486–14495 (2019).

24. Berardi, D. et al. 21st-century biogeochemical modeling:

challenges for Century-based models and where do we go from
here? Glob. Chang. Biol. Bioenergy

, 774–788 (2020).

25. Kleber, M. et al. Dynamic interactions at the mineral–organic

matter interface. Nat. Rev. Earth Environ.

, 402–421 (2021).

26. Melillo, J. M. et al. Long-term pattern and magnitude of soil

carbon feedback to the climate system in a warming world.
Science

358

, 101–105 (2017).

27. Walker, T. W. N. et al. A systemic overreaction to years versus

decades of warming in a subarctic grassland ecosystem. Nat.
Ecol. Evol.

, 101–108 (2020).

28. Bradford, M. Thermal adaptation of decomposer communities in

warming soils. Front. Microbiol.

, 333 (2013).

29. Hartley, I. P., Hill, T. C., Chadburn, S. E. & Hugelius, G. Temperature

effects on carbon storage are controlled by soil stabilisation
capacities. Nat. Commun.

, 6713 (2021).

30. Doetterl, S. et al. Soil carbon storage controlled by interactions

between geochemistry and climate. Nat. Geosci.

, 780–783

(2015).

31. Abramoff, R. Z., Torn, M. S., Georgiou, K., Tang, J. & Riley, W. J. Soil

organic matter temperature sensitivity cannot be directly inferred
from spatial gradients. Glob. Biogeochem. Cycles

, 761–776

(2019).

32. Luo, Z. et al. Convergent modelling of past soil organic carbon

stocks but divergent projections. Biogeosciences

, 4373–4383

(2015).

33. Georgiou, K. et al. Global stocks and capacity of

mineral-associated soil organic carbon. Nat. Commun.

, 3797

(2022).

34. Eyring, V. et al. Overview of the Coupled Model Intercomparison

Project Phase 6 (CMIP6) experimental design and organization.
Geosci. Model Dev.

, 1937–1958 (2016).

35. Jobbágy, E. G. & Jackson, R. B. The vertical distribution of soil

organic carbon and its relation to climate and vegetation. Ecol.
Appl.

, 423–436 (2000).

36. Ahrens, B. et al. Combination of energy limitation and sorption

capacity explains

C depth gradients. Soil Biol. Biochem.

148

107912 (2020).

37. Gentsch, N. et al. Temperature response of permafrost soil

carbon is attenuated by mineral protection. Glob. Chang. Biol.

3401–3415 (2018).

38. Giardina, C. P. & Ryan, M. G. Evidence that decomposition rates

of organic carbon inmineral soil do not vary with temperature.
Nature

404

, 858–861 (2000).

39. Giardina, C. P., Litton, C. M., Crow, S. E. & Asner, G. P.

Warming-related increases in soil CO

efflux are explained

by increased below-ground carbon flux. Nat. Clim. Chang.

822–827 (2014).

40. Wang, C. et al. The temperature sensitivity of soil: microbial

biodiversity, growth, and carbon mineralization. ISME J.

2738–2747 (2021).

41. Shi, Z. et al. The age distribution of global soil carbon inferred

from radiocarbon measurements. Nat. Geosci.

, 555–559

(2020).

42. He, Y. et al. Radiocarbon constraints imply reduced carbon uptake

by soils during the 21st century. Science

353

, 1419–1424

(2016).

43. Todd-Brown, K. E. O. et al. Causes of variation in soil carbon

simulations from CMIP5 Earth system models and comparison
with observations. Biogeosciences

, 1717–1736 (2013).

44. Ito, A. et al. Soil carbon sequestration simulated in CMIP6-LUMIP

models: implications for climatic mitigation. Environ. Res. Lett.

124061 (2019).

45. Abramoff, R. Z. et al. Improved global-scale predictions of soil

carbon stocks with Millennial Version 2. Soil Biol. Biochem.

164

108466 (2022).

46. Zimmermann, M., Leifeld, J., Schmidt, M. W. I., Smith, P. & Fuhrer,

J. Measured soil organic matter fractions can be related to pools
in the RothC model. Eur. J. Soil Sci.

, 658–667 (2007).

47. Pierson, D. et al. Optimizing process-based models to predict

current and future soil organic carbon stocks at high-resolution.
Sci. Rep.

, 10824 (2022).

48. Koven, C. et al. On the formation of high-latitude soil carbon

stocks: effects of cryoturbation and insulation by organic
matter in a land surface model. Geophys. Res. Lett.

https://doi.

org/10.1029/2009GL040150

(2009).

49. 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

, 7109–7131 (2013).

50. Plaza, C. et al. Direct observation of permafrost degradation and

rapid soil carbon loss in tundra. Nat. Geosci.

, 627–631

(2019).

51. Parton, W. J., Schimel, D. S., Cole, C. V. & Ojima, D. S. Analysis

of factors controlling soil organic matter levels in Great Plains
grasslands. Soil Sci. Soc. Am. J.

, 1173–1179 (1987).

52. Haaf, D., Six, J. & Doetterl, S. Global patterns of geo-ecological

controls on the response of soil respiration to warming. Nat. Clim.
Chang.

, 623–627 (2021).

Page 8

Nature Geoscience

| Volume 17 | March

2024 | 205–212

212

Article

https://doi.org/10.1038/s41561-024-01384-7

53. Jones, C. D. et al. C4MIP – the Coupled Climate–Carbon Cycle

Model Intercomparison Project: experimental protocol for CMIP6.
Geosci. Model Dev.

, 2853–2880 (2016).

54. Bouskill, N. J., Riley, W. J., Zhu, Q., Mekonnen, Z. A. & Grant, R. F.

Alaskan carbon-climate feedbacks will be weaker than
inferred from short-term experiments. Nat. Commun.

5798 (2020).

Publisher’s note

Springer Nature remains neutral with regard to

jurisdictional claims in published maps and institutional affiliations.

Open Access

This article is licensed under a Creative Commons

Attribution 4.0 International License, which permits use, sharing,

adaptation, distribution and reproduction in any medium or format,
as long as you give appropriate credit to the original author(s) and the
source, provide a link to the Creative Commons licence, and indicate
if changes were made. The images or other third party material in this
article are included in the article’s Creative Commons licence, unless
indicated otherwise in a credit line to the material. If material is not
included in the article’s Creative Commons licence and your intended
use is not permitted by statutory regulation or exceeds the permitted
use, you will need to obtain permission directly from the copyright
holder. To view a copy of this licence, visit

http://creativecommons.

org/licenses/by/4.0/

Physical and Life Sciences Directorate, Lawrence Livermore National Laboratory, Livermore, CA, USA.

Climate and Ecosystem Sciences Division,

Lawrence Berkeley National Laboratory, Berkeley, CA, USA.

Climate and Global Dynamics Laboratory, National Center for Atmospheric Research,

Boulder, CO, USA.

Institute of Arctic and Alpine Research, University of Colorado, Boulder, CO, USA.

Natural Resource Ecology Laboratory, Colorado

State University, Fort Collins, CO, USA.

Department of Life and Environmental Sciences, University of California Merced, Merced, CA, USA.

Ronin

Institute, Montclair, NJ, USA.

Department of Ecology and Evolutionary Biology, Yale University, New Haven, CT, USA.

Department of Physical Geography

and Ecosystem Science, Lund University, Lund, Sweden.

Department of Plant Sciences, University of Cambridge, Cambridge, UK.

Rocky Mountain

Research Station, United States Forest Service, Boise, ID, USA.

Climate Change Science Institute and Environmental Sciences Division, Oak Ridge

National Laboratory, Oak Ridge, TN, USA.

Department of Earth System Science, Stanford University, Stanford, CA, USA.

Woods Institute for the

Environment and Precourt Institute for Energy, Stanford University, Stanford, CA, USA. e-mail:

georgiou1@llnl.gov

Page 9

Nature Geoscience

Article

https://doi.org/10.1038/s41561-024-01384-7

Methods

Data sources and processing

The overarching goal of this study was to use observed spatial gradients
in soil C stocks—particularly in underlying unprotected (particulate)
and protected (mineral-associated) C pools—and their relation-
ships with temperature, to develop benchmarks for ESMs. Although
such benchmarks derived from spatial gradients are insufficient
for constraining transient responses or making projections (that is,
space-for-time assumptions)

, they do reflect steady-state differ-

ences and long-term climatological temperature sensitivities of soil
C pools. These observed differences across spatial climate gradients
should, in principle, be predicted by dynamical soil C models, and,
thus, these benchmarks serve as important tests

. However, meet-

ing this constraint is not sufficient for ensuring the ability of models
to produce accurate projections under global change, and so we urge
experimental work that quantifies the transient responses of soil C
pools across diverse climates and biomes

Observationally derived climate and edaphic datasets for all global

analyses were used at 0.5 × 0.5° resolution. MAT was estimated from
the CRU dataset (version 3.10)

and mean annual precipitation from

the GPCC dataset

, both as 30 year annual averages. Land cover was

obtained from the MODIS MCD12C1 product

, and productivity was

estimated from the MODIS net primary productivity product

. Fol-

lowing Shi et al.

, soil organic C stocks to 1 m depth were estimated as

the mean of Harmonized World Soil Database

and SoilGrids

maps

to reflect a best estimate and uncertainty across data products. Cor-
responding mineral-associated C stocks and 90% prediction intervals
were estimated in Georgiou et al.

using a machine-learning (random

forest) algorithm, and multiple cross-validation approaches were used
to provide rigorous assessments of predictability

. We also note that

our conclusions were robust to the selected fractionation method
(namely, size and density; Supplementary Fig. 5 and Supplementary
Table 5) and reflect a global best estimate and uncertainty given exist-
ing observations. However, we encourage further measurements in
cooler climates and sandy soils, where the relative contribution of fine
particulate organic matter may be more substantial

. In the present

study, we focused our analyses on non-permafrost mineral soils with
MAT > 0 °C given data limitations at higher latitudes

. To exclude

soils whose decomposition is limited by water saturation or aridity,
we also removed gridcells that contain a coverage of more than 50%
peat, defined as the Histosol soil order or Histel suborder, and those
that are hyperarid (with an aridity index of <0.05) or receive less than
100 mm yr

−1

of precipitation. These gridcells were masked in white in

all maps (Fig.

) and were excluded from subsequent analyses.

Global land model output

Model output was sourced from CMIP6 ESMs

and three offline bio-

geochemical testbed models

for present-day (‘historical’) and future

projections (‘ssp585’, or RCP8.5). Only CMIP6 models that reported soil
C pool distributions were used in this study. Model details are listed
in Supplementary Tables 1–3. Carbon stocks and climate covariates
from each model were taken over the historical period as ten year
annual averages for C stocks and 30 year annual averages for cli-
mate (spanning 2005–2015 and 1985–2015, respectively) to match
the observation-based quantities. Specifically, MAT (‘tas’ for CMIP6
models), mean annual precipitation (‘pr’) and net primary productivity
(‘npp’) were used from each model. Since the models do not provide
depth-resolved C stocks, we used total soil C stocks (‘cSoil’ for CMIP6
models) in our analyses. Although the depth of soil C in CMIP6 models
is not clearly defined, previous studies have compared with 0–1 m
integrated soil C stocks

. We believe that this represents the best

benchmark on current models and follow this depth correspondence
in comparison with the data, but note that any errors in the corre-
spondence with depth could be reduced with future models that are
depth-resolved and/or explicitly define a depth interval

For the offline models of the biogeochemical testbed (that is,

CASA-CNP, MIMICS and CORPSE), mineral-associated (protected) C
was estimated by the explicit physicochemically protected pools for
MIMICS and CORPSE, and the slowest-cycling (passive) pool for the
first-order CASA-CNP model (see Wieder et al.

for details on physico-

chemical protection within the three biogeochemical testbed models).
Analogously, for the first-order CMIP6 models, the slowest-cycling
C pool (‘cSoilSlow’) was used to represent mineral-protected C (see
‘Model pool interpretability’ in the Methods). We note that the CMIP6
model output repository uses the terminology cSoil for total soil
C stocks and ‘cSoilFast’, ‘cSoilMedium’ and cSoilSlow for the three
soil C pools that correspond to the ‘active’, ‘slow’, and passive pools,
respectively.

Of the ten CMIP6 ESMs that reported soil C pool distributions,

only one model constituted a conceptualization that did not align
with mineral protection; namely, in TaiESM1, the soil C dynamics
were parameterized as a decomposition cascade based on a synthesis
of microcosm studies using radio-labelled substrates

–

. We thus

excluded TaiESM1 from our analyses, and note that the model can
be redefined and parameterized for comparison with pool-specific
benchmarks in future work. We also emphasize that not all first-order
soil C models were formulated as a decomposition cascade based on
recalcitrance. Rather, several first-order soil C models implicitly repre-
sent key processes by rerouting flows and imposing particular controls
(for example, lignin, soil texture and climate) on the formation of select
pools

(see ‘Model pool interpretability’ in the Methods; Supplemen-

tary Fig. 12). To the best of our knowledge, the nine remaining CMIP6
ESMs (namely, ACCESS-ESM1-5, BCC-CSM2-MR, CESM2, CNRM-ESM2-1,
E3SM-1-1-ECA, MIROC-ES2L, MRI-ESM2-0, IPSL-CM6A-LR, NorESM2)
constitute formulations broadly based on the Century model

and its

derivatives, for which the slowest-cycling (passive) pool is conceptually
aligned with physicochemical protection. Although parameteriza-
tions in implemented soil models may have diverged from published
values, protected C appears to be defined consistently across the nine
remaining CMIP6 ESMs (Supplementary Table 3). To facilitate future
benchmarking studies and reduce potential uncertainties, we highly
encourage more detailed documentation on parameterizations imple-
mented in soil models within CMIP ESMs.

The protected C pool in each model was then compared with the

total soil C for each gridcell globally to calculate the proportion of
protected C (Supplementary Fig. 4). The global mean (with 90% CIs)
was calculated for each model and in each indicated temperature
regime (Fig.

and Supplementary Fig. 5). We then calculated the

spatially explicit multi-model mean across the ensemble of soil C
models (Fig.

), as well as the latitudinal means of the observation-

ally derived data product and multi-model ensemble (Fig.

). We

also quantified changes in soil C pools for each model under an
RCP8.5 (or ssp585) scenario by comparing projections from the end
of the twenty-first century (2090–2100) with the present day (2005–
2015). Results are shown as absolute and relative changes in soil C
stocks globally (Supplementary Figs. 6 and 7 and Supplementary
Table 4) and across time and latitude (as Hovmöller diagrams;
Supplementary Fig. 8).

Model pool interpretability

The comparison of protected C across different model formulations,
following the above pool nomenclature, has been proposed in several
past studies with first-order models

. We recognize that there are

still inherent mismatches between operationally defined, measurable
pools of soil organic C fractions and the modelled states that are simu-
lated by different soil biogeochemical models

–

. Although progress

is being made to help bridge this divide

, the global-scale models

that are the focus of this work make assumptions that are broadly
consistent with characteristics of mineral-associated and particulate
soil C fractions

–

(Supplementary Fig. 12).

Page 10

Nature Geoscience

Article

https://doi.org/10.1038/s41561-024-01384-7

Indeed, in non-permafrost mineral soils, the formation and com-

position of the passive pool in select first-order models (for example,
in the Century model of Parton et al.

) was formulated to resemble the

mineral-associated fraction, and mirrors the mathematical representa-
tions of physicochemically protected C in several microbially explicit
models, including MIMICS and CORPSE

. Namely, the protected

(passive in the first-order models) pool is formed primarily from the
microbial (active) pool, where this transfer coefficient is modulated
by soil clay and silt content, such that proportionally more C is stored
within the protected pool in fine-textured soils

–

. More decom-

posable (metabolic) litter enters only the microbial (active) pool,
whereas structural litter (controlled by lignin content) primarily enters
the particulate (slow) pool

. A small amount of the protected pool

forms from direct sorption of the particulate pool

. Consequently,

the C:N ratio of the protected (passive) pool often resembles that of the
microbial (active) pool (that is, typically ~11–12 for the protected pool
in the Century model

), which is consistent with the stoichiometry of

mineral-associated organic matter

. By contrast, the C:N ratio of the

slow pool is more closely aligned with the stoichiometry of particulate
organic matter (that is, ~12–20 in the Century model

Finally, losses from the protected pool are represented as a

first-order process in most soil C models

, including the ones

presented here, with typical intrinsic turnover times ranging between
75 and 1,000 years, which are broadly consistent with radiocarbon
estimates

(see Supplementary Table 3 and Supplementary Fig. 13).

Whereas there is heterogeneity within all underlying soil fractions in
practice, it is well accepted that mineral-associated organic C is gener-
ally slower-cycling (decades to millennia) than particulate organic C
(days to decades)

–

. Indeed, the intrinsic turnover times for protected

(passive) C ranged from 200 to 450 years in all analysed CMIP6 ESMs,
except for MIROC-ES2L and MRI-ESM2-0 for which the turnover times
were 1,000 years (Supplementary Table 3). By contrast, intrinsic turno-
ver times of the particulate (slow) pool ranged from 5 to 10 years in all
analysed CMIP6 ESMs, except for MIROC-ES2L for which the turnover
time was 25 years. These intrinsic turnover times are baseline (mini-
mum) pool-specific values, such that environmental modifiers (based
on temperature and moisture) further increase the turnover time in
most gridcells; for example, by a factor of ~2.5 in temperate forests
up to a factor of ~10 in Arctic regions for the Century model

. These

environmental controls are important when comparing pool-specific
turnover times with global radiocarbon data (Supplementary Fig. 14),
and future work could estimate the distribution of C ages for each
modelled pool across gridcells globally. Leveraging Δ

C measurements

from the International Soil Radiocarbon Database

, and following

Shi et al.

to convert between turnover time and Δ

C, we show that

the protected (passive) pool turnover time within many first-order
models is broadly aligned with mean C ages of the mineral-associated
fraction (Supplementary Table 3 and Supplementary Figs. 13 and 14).
However, future model versions may further improve the correspond-
ence by reducing the intrinsic turnover time of protected C, especially
in MIROC-ES2L and MRI-ESM2-0, and strengthening its environmental
controls to better match the global distribution of Δ

C measurements

(Supplementary Fig. 14).

By contrast, the turnover time of protected C is fixed at 75 years in

CORPSE and ranges from 15 to 80 years in MIMICS based on a depend-
ence on soil texture

. We note that the turnover times of MIMICS and

CORPSE appear to be biased low compared with radiocarbon data
of soil fractions (Supplementary Figs. 13 and 14), which may explain
their strong responses under future conditions (Supplementary
Figs. 6–8) and suggests that bulk soil C ages may also be too young in
these models. In fact, subsequent versions of MIMICS that have yet to be
integrated into the global biogeochemical testbed have implemented
much longer turnover times of protected C—for instance, ranging from
180 to 1,000 years (that is, increased by a factor of 12.5) in MIMICS-CN

and in the C-only MIMICS independently fitted to site-level data at

the Reynolds Creek Experimental Watershed

. These modifications

in subsequent versions of MIMICS have also resulted in much higher
proportions of protected C (namely, 50–80%) that are better aligned
with observational constraints. Thus, we emphasize that all models,
including microbially explicit and first-order formulations, have the
opportunity to be falsified or refined in future iterations, but this does
not preclude the ability to benchmark them here.

The range of conceptual representations and parameterizations of

physicochemical protection across ecosystem to global-scale models
(for example, see the schematic in Sulman et al.

) presents an urgent

need to develop and apply observationally derived constraints. We
posit that both the microbially explicit and first-order soil C models in
this study (Supplementary Tables 1–3) constitute conceptualizations
that are compatible with benchmarking protected soil C to refine future
model parameterizations. Although we recognize that there are still
inherent uncertainties in the interpretability of modelled pools and
operationally defined fractions, we believe that this study represents
a valuable benchmark on model pools globally. Observational synthe-
ses show that mineral-associated C is typically an old (Supplementary
Fig. 13 and Supplementary Table 3) and large (Fig.

and Supplementary

Fig. 5) component of total soil C. Thus, the slowest-cycling passive
pool in first-order models—which broadly aligns with radiocarbon
measurements of mineral-associated C (Supplementary Figs. 13 and
14)—should similarly represent a large fraction of the total soil C pool
in mineral soils. Such a correspondence would also improve model
performance metrics based on bulk soil C ages

. Since pool-specific

turnover times appear to be broadly aligned with observations, future
studies may seek to refine other parameters that strongly affect the
proportion of protected C—for example, increasing the fraction of
the active pool that flows to the passive pool (Supplementary Figs. 12
and 15) and its mineralogical controls

. Such parameters may also be

informed by observational estimates of the relative flows of plant- or
microbial-derived C to mineral-associated C

. Ultimately, we con-

tend that benchmarking soil C pool distributions concurrently with
radiocarbon ages in future work will present a powerful constraint on
soil C dynamics across global models.

Data analysis

We investigated the relationship between soil C stocks (within unpro-
tected and protected C pools) and MAT across the global data product
and models. Soil C stocks for each pool were log-transformed and
a linear regression was used to reflect a hypothesized exponential
relationship with temperature (as in Hartley et al.

). To ensure that

the observed global patterns were not driven by confounding envi-
ronmental variables, we used a multiple linear regression to control
for net primary productivity, mean annual precipitation and clay and
silt content in our analyses. In all cases, the climate and productivity
covariates were sourced from their respective model output or data
product. Clay and silt content was obtained from the Harmonized
World Soil Database

as it was not reported for most CMIP6 models.

The climatological temperature sensitivity was then calculated from
the slope of the log-linear temperature relationship and, following
Hartley et al.

, summarized as the proportional decline in C stocks for

every 10 °C (analogous to a climatological Q

, but for C stocks) with

corresponding parametric 95% CIs (Supplementary Table 6).

This calculation was then repeated for protected and unprotected

soil C pools in the data and models, across all mean temperatures
above 0 °C (Fig.

) and across two underlying temperature regimes,

namely cool (<15 °C) and warm (≥15 °C) (Fig.

and Supplementary

Fig. 9). Soils with mean temperatures below 0 °C were not considered
in our analyses due to distinct controls on soil C decomposition (for
example, permafrost dynamics

) and consequent implications on

the interpretability of protected pools in first-order models in these
regions (Supplementary Fig. 3), as well as higher uncertainties in the
data product within this temperature regime

. Furthermore, the same

Page 11

Nature Geoscience

Article

https://doi.org/10.1038/s41561-024-01384-7

analyses were performed on an observational synthesis

(Supplemen-

tary Figs. 10 and 11). The general agreement between the observations
and data product provides further evidence for the robustness of
the reported temperature sensitivity trends between protected and
unprotected soil C across climate regimes (Fig.

). Differences in clima-

tological temperature sensitivity between protected and unprotected
soil C pools, and across the data and models, were deemed statistically
significant when the 95% CIs did not overlap (Supplementary Table 6).

Data availability

All data used in this manuscript are fully open-access and available at
the references provided. The CMIP6 model output is available from
the Earth System Grid Federation at

https://esgf-node.llnl.gov/search/

cmip6/

, the Biogeochemical Testbed model output is available on the

Climate Data Gateway at

https://doi.org/10.5065/d6nc600w

and the

observationally derived dataset is available on Zenodo at

https://doi.

org/10.5281/zenodo.6539765

Code availability

All code relating to this study is available from the corresponding author
on GitHub at

https://github.com/katerina-georgiou/global-soilCpools

and is archived on Zenodo at

https://doi.org/10.5281/zenodo.10515706

References

55. Harris, I., Jones, P. D., Osborn, T. J. & Lister, D. H. Updated

high-resolution grids of monthly climatic observations – the CRU
TS3.10 Dataset. Int. J. Climatol.

, 623–642 (2014).

56. Schneider, U. et al. GPCC Full Data Reanalysis Version 6.0 at

0.5°: Monthly Land-Surface Precipitation from Rain-Gauges
built on GTS-based and Historic Data (GPCC, 2011);

https://doi.

org/10.5676/DWD_GPCC/FD_M_V6_050

57. Friedl, M. A. et al. MODIS Collection 5 global land cover:

algorithm refinements and characterization of new datasets.
Remote Sens. Environ.

114

, 168–182 (2010).

58. Zhao, M., Heinsch, F. A., Nemani, R. R. & Running, S. W.

Improvements of the MODIS terrestrial gross and net primary
production global data set. Remote Sens. Environ.

, 164–176

(2005).

59. Wieder, W. R., Boehnert, J., Bonan, G. B. & Langseth, M. Regridded

Harmonized World Soil Database v1.2 (ORNL DAAC, 2014);

https://

doi.org/10.3334/ORNLDAAC/1247

60. Hengl, T. et al. SoilGrids250m: global gridded soil information

based on machine learning. PLoS ONE

https://doi.org/10.1371/

journal.pone.0169748

(2017).

61. Georgiou, K. et al. Globally-gridded data for manuscript: global

stocks and capacity of mineral-associated soil organic carbon.
Zenodo

https://doi.org/10.5281/zenodo.6539765

(2022)

62. Ploton, P. et al. Spatial validation reveals poor predictive

performance of large-scale ecological mapping models. Nat.
Commun.

, 4540 (2020).

63. Meyer, H. & Pebesma, E. Machine learning-based global maps of

ecological variables and the challenge of assessing them. Nat.
Commun.

, 2208 (2022).

64. Thornton, P. E. & Rosenbloom, N. A. Ecosystem model spin-up:

estimating steady state conditions in a coupled terrestrial carbon
and nitrogen cycle model. Ecol. Modell.

189

, 25–48 (2005).

65. Lawrence, D. M. et al. Parameterization improvements

and functional and structural advances in Version 4 of the
Community Land Model. J. Adv. Model. Earth Syst.

https://doi.

org/10.1029/2011MS00045

(2011).

66. Lee, W. L. et al. Taiwan Earth System Model Version 1: description

and evaluation of mean state. Geosci. Model Dev.

, 3887–3904

(2020).

67. Christensen, B. T. in Evaluation of Soil Organic Matter Models Vol.

38 (eds Powlson, D. S. et al.) 143–159 (Springer, 1996).

68. Elliot, E. T., Paustian, K. & Frey, S. D. in Evaluation of Soil Organic Matter

Models Vol. 38 (eds Powlson, D. S. et al.) 161–179 (Springer, 1996).

69. Lehmann, J. & Kleber, M. The contentious nature of soil organic

matter. Nature

528

, 60–68 (2015).

70. Zhang, Y. et al. Simulating measurable ecosystem carbon and

nitrogen dynamics with the mechanistically defined MEMS 2.0
model. Biogeosciences

, 3147–3171 (2021).

71. Dangal, S. R. S. et al. Improving soil carbon estimates by linking

conceptual pools against measurable carbon fractions in
the DAYCENT model version 4.5. J. Adv. Model. Earth Syst.

e2021MS002622 (2022).

72. Abramoff, R. et al. The Millennial model: in search of measurable

pools and transformations for modeling soil carbon in the new
century. Biogeochemistry

137

, 51–71 (2018).

73. Parton, W. J. & Rasmussen, P. E. Long-term effects of crop

management in wheat-fallow: II. CENTURY model simulations.
Soil Sci. Soc. Am. J.

, 530–536 (1994).

74. Cotrufo, M. F., Ranalli, M. G., Haddix, M. L., Six, J. & Lugato, E. Soil

carbon storage informed by particulate and mineral-associated
organic matter. Nat. Geosci.

, 989–994 (2019).

75. Rocci, K. S. et al. Aligning theoretical and empirical

representations of soil carbon-to-nitrogen stoichiometry with
process-based terrestrial biogeochemistry models. Soil Biol.
Biochem.

189

, 109272 (2024).

76. Lawrence, C. R. et al. An open-source database for the synthesis

of soil radiocarbon data: International Soil Radiocarbon Database
(ISRaD) version 1.0. Earth Syst. Sci. Data

, 61–76 (2020).

77. Kyker-Snowman, E., Wieder, W. R., Frey, S. D. & Grandy, A. S.

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

, 4413–4434 (2020).

78. Georgiou, K. et al. Divergent controls of soil organic

carbon between observations and process-based models.
Biogeochemistry

, 775–792 (2021).

79. Angst, G., Mueller, K. E., Nierop, K. G. J. & Simpson, M. J. Plant- or

microbial-derived? A review on the molecular composition of
stabilized soil organic matter. Soil Biol. Biochem.

156

, 108189 (2021).

Acknowledgements

K.G. was supported as a Lawrence Fellow at Lawrence Livermore
National Laboratory (LLNL) by the LLNL-LDRD Program under Project
No. 21-ERD-045 and 24-LW-053. J.P.-R. and E.W.S. were supported by the
US Department of Energy (DOE) Office of Science, Office of Biological
and Environmental Research, Genomic Science Program as part of
the LLNL Microbes Persist Scientific Focus Area, SCW1632. Work at
LLNL was conducted under the auspices of DOE Contract DE-AC52-
07NA27344. W.R.W. was supported by National Science Foundation
Grants 1926413 and 2031238 and by USDA NIFA-AFRI Grant 2020-
67019-31395. C.D.K., W.J.R. and Q.Z. were supported by the US DOE
Biological and Environmental Research (BER) Program at LBNL under
DOE Contract DE-AC02-05CH11231 through the Regional and Global
Model Analysis Program (RUBISCO SFA). N.J.B. was supported by the US
DOE BER Early Career Research Program under Contract FP00005182.
R.Z.A. and B.N.S. were supported by the US DOE BER at Oak Ridge
National Laboratory under DOE Contract DE-AC05-00OR22725. A.A.
acknowledges the research environment Biodiversity and Ecosystem
Services in a Changing Climate (BECC) at Lund University and funding
from the Swedish Research Council (2021-05344).

Author contributions

K.G. conceived and designed the study, based on ideas developed
through discussions with C.D.K., W.R.W., W.J.R. and R.B.J. All authors
contributed to the study development. K.G., W.R.W., M.D.H. and Q.Z.
collected and contributed to the model output. M.D.H., W.J.P., W.R.W.,
C.D.K., A.A., R.Z.A., D.P. and B.N.S. provided insights on the model

Georgiou et al. (2024). Emergent temperature sensitivity of soil organic carbon driven by mineral associations

Page 1

Page 2

Page 3

Page 4

Page 5

Page 6

Page 7

Page 8

Page 9

Page 10

Page 11

Page 12

Document Outline