# Process-oriented analysis of dominant sources of uncertainty in the land carbon sink

Michael O'sullivan, Pierre Friedlingstein, Stephen Sitch, Peter Anthoni, Almut Arneth, Vivek Arora, Vladislav Bastrikov, Christine Delire, Daniel Goll, Atul Jain, et al.

► To cite this version:

Michael O'sullivan, Pierre Friedlingstein, Stephen Sitch, Peter Anthoni, Almut Arneth, et al. Process-oriented analysis of dominant sources of uncertainty in the land carbon sink. Nature Communications, 2022, 13 (1), pp.4781. 10.1038/s41467-022-32416-8. hal-03763263

HAL Id: hal-03763263

https://hal.science/hal-03763263v1

Submitted on 29 Aug 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_2.png" /></figure>

# Process-oriented analysis of dominant sources of uncertainty in the land carbon sink

Received: 11 February 2022

Accepted: 28 July 2022

Published online: 15 August 2022

<figure><img src="image_3.png" /></figure>

Check for updates

Michael O'Sullivan<fnref n="1" />✉, Pierre Friedlingstein ①,², Stephen Sitch ③,
Peter Anthoni⁴, Almut Arneth⁴, Vivek K. Arora⁵, Vladislav Bastrikov⁶,
Christine Delire ⑦, Daniel S. Goll ⑥, Atul Jain ⑧, Etsushi Kato ⑨,
Daniel Kennedy¹⁰, Jürgen Knauer ¹¹,¹², Sebastian Lienert ⑬¹³,
Danica Lombardozzi ⑬¹⁰, Patrick C. McGuire¹⁴, Joe R. Melton ⑤,
Julia E. M. S. Nabel ⑬¹⁵,¹⁶, Julia Pongratz¹⁵,¹⁷, Benjamin Poulter ⑬¹⁸,
Roland Séférian⁷, Hanqin Tian¹⁹, Nicolas Vuichard⁶, Anthony P. Walker ⑬²⁰,
Wenping Yuan ⑬²¹, Xu Yue ⑬²² & Sönke Zaehle ⑬¹⁶

The observed global net land carbon sink is captured by current land models. All models agree that atmospheric CO₂ and nitrogen deposition driven gains in carbon stocks are partially offset by climate and land-use and land-cover change (LULCC) losses. However, there is a lack of consensus in the partitioning of the sink between vegetation and soil, where models do not even agree on the direction of change in carbon stocks over the past 60 years. This uncertainty is driven by plant productivity, allocation, and turnover response to atmospheric CO₂ (and to a smaller extent to LULCC), and the response of soil to LULCC (and to a lesser extent climate). Overall, differences in turnover explain ~70% of model spread in both vegetation and soil carbon changes. Further analysis of internal plant and soil (individual pools) cycling is needed to reduce uncertainty in the controlling processes behind the global land carbon sink.

<footnote n="1">¹College of Engineering, Mathematics and Physical Sciences, University of Exeter, Exeter EX4 4QF, UK. ²Laboratoire de Météorologie Dynamique, Institut Pierre-Simon Laplace, CNRS-ENS-UPMC-X, Paris, France. ³College of Life and Environmental Sciences, University of Exeter, Exeter EX4 4RJ, UK. ⁴Karlsruhe Institute of Technology, Institute of Meteorology and Climate Research/Atmospheric Environmental Research, 82467 Garmisch-Partenkirchen, Germany. ⁵Canadian Centre for Climate Modelling and Analysis, Climate Research Division, Environment and Climate Change Canada, Victoria, BC, Canada. ⁶Laboratoire des Sciences du Climat et de l'Environnement, LSCE/IPSL, CEA-CNRS-UVSQ, Université Paris-Saclay, F-91198 Gif-sur-Yvette, France. ⁷CNRM, Université de Toulouse, Météo-France, CNRS, Toulouse, France. ⁸Department of Atmospheric Sciences, University of Illinois, Urbana, IL 61821, USA. ⁹Institute of Applied Energy (IAE), Minato-ku, Tokyo 105-0003, Japan. ¹⁰National Center for Atmospheric Research, Climate and Global Dynamics, Terrestrial Sciences Section, Boulder, CO 80305, USA. ¹¹Hawkesbury Institute for the Environment, Western Sydney University, Penrith, NSW, Australia. ¹²CSIRO Oceans and Atmosphere, Canberra, ACT 2101, Australia. ¹³Climate and Environmental Physics, Physics Institute and Oeschger Centre for Climate Change Research, University of Bern, Bern, Switzerland. ¹⁴Department of Meteorology, University of Reading, Reading, UK. ¹⁵Max Planck Institute for Meteorology, Bundesstr. 53, 20146 Hamburg, Germany. ¹⁶Max Planck Institute for Biogeochemistry, Jena, Germany. ¹⁷Ludwig-Maximilians-Universität München, Luisenstr. 37, 80333 München, Germany. ¹⁸NASA Goddard Space Flight Center, Biospheric Sciences Laboratory, Greenbelt, MD 20771, USA. ¹⁹Schiller Institute for Integrated Science and Society, Department of Earth and Environmental Sciences, Boston College, Chestnut Hill, MA 02467, USA. ²⁰Climate Change Science Institute & Environmental Sciences Division, Oak Ridge National Lab, Oak Ridge, TN 37831, USA. ²¹School of Atmospheric Sciences, Sun Yat-sen University, Zhuhai, Guangdong 510245, China. ²²School of Environmental Science and Engineering, Nanjing University of Information Science and Technology (NUIST), Nanjing, China. ✉e-mail: m.osullivan@exeter.ac.uk</footnote>
Over the last 60 years, there has been a continuous rise in anthro-
pogenic CO₂ emissions. Around the equivalent of a quarter of these
emissions have been taken up by the land biosphere (known as the
natural land sink), acting as a strong negative feedback to mitigate
climate change¹. To be able to project the global carbon cycle (and
hence the climate response) in the future, we need to understand the
underlying processes and their timescales that drive the contemporary
land sink. There are many distinct but interdependent mechanisms
that regulate the flow of carbon into, through, and out of the land. The
timescales at which these processes (e.g. photosynthesis, allocation,
plant growth, litterfall, plant mortality and soil turnover) act on the
carbon cycle range from days to centuries, and the interplay between
these determines the changes in land carbon storage².

Recent work suggests the global net land sink is located in northern latitudes³. The attribution of drivers is highly uncertain, but a combination of increasing atmospheric CO₂ concentrations⁴, reactive nitrogen (N) deposition⁵ and atmospheric warming⁶ are likely responsible. Tropical lands are probably nearer net zero carbon sinks due to large land-use and land cover changes (LULCC) carbon losses counteracting the 'natural' sink⁷,⁸. Atmospheric observations also show that the Amazon forest is currently carbon neutral ('natural' sink = LULCC-driven source)⁹, and this provides further evidence of a limited net tropical sink. Yet, the largest gross fluxes between the land and atmosphere are often found in tropical regions¹⁰, and so changes in tropical ecosystem functioning can have significant global impacts.

Process-based dynamic global vegetation models (DGVMs) that simulate processes of carbon uptake and release can help to elucidate the roles of individual drivers (rising atmospheric CO₂, changes in climate, nutrient deposition, and LULCC), attribute the change to processes, and quantify regional sinks over timescales needed given the multitude of timescales over which terrestrial processes act (i.e. beyond the period of reliable empirical and remote-sensed data).

DGVMs are used to estimate the natural land sink (hereafter
simply land sink) as part of the global carbon budget (GCB⁸). Globally,
the DGVM multi-model mean estimate of the global land sink is con-
sistent with the global carbon budget residual land sink (the difference
between fossil and land-use emissions and atmospheric and ocean
sinks, see Table 5 in ref. 8), however, there is a significant spread across
models, and here we show the spread widens at regional scales, when
quantifying changes in vegetation and soil carbon, or attributing
changes in internal processes to external drivers. Therefore, while
model ensembles are helpful in analysing global-scale processes, they
must be interpreted in the context of their process representation and
unique and collective biases²,⁴.

These deficiencies are of high significance when trying to under-
stand past changes but also when DGVMs are used to make predictions
of future carbon cycling¹¹, as limited trust can be placed in projections
when certain fundamental processes and carbon-climate feedbacks

are not fully captured². Therefore, it is of high priority to identify the
leading sources of uncertainty in DGVMs and understand the rela-
tionship between drivers and processes on varying spatial scales. In
this study, we use a triple (3D-matrix) approach to identify where
models agree and, just as importantly, disagree, and thus guide future
modelling efforts:

* What are the (1) external drivers (concurrent rises in atmospheric CO₂ and N deposition, climate and LULCC), (2) main regions (tropics and extra tropics) and (3) processes (production vs turnover) primarily responsible for the changes in the global net land carbon sink?

We use the suite of 18 DGVMs from the GCB2021 (TRENDYv10; ref. 8) to quantify changes in net carbon exchange and carbon stocks over the period 1959–2020. TRENDYv10 provides a set of simulations to attribute these changes to drivers, and we use a process attribution framework to decompose changes in carbon stocks into those driven by productivity and turnover separately. This framework enables us to express productivity and turnover-induced changes in carbon stocks in units of PgC, which allows for a direct comparison between both processes, which have units of carbon per unit time and time, respectively (see Methods).

Results

Drivers of global and regional land sinks

The global net land sink is derived from the difference between the
fossil fuel emissions (E<sub>FOS</sub>) and the CO<sub>2</sub> accumulation in the atmo-
sphere (G<sub>ATM</sub>) and uptake by the oceans (S<sub>OCEAN</sub>). We refer to this
estimate of the net land sink as the 'observed' GCB budget constraint (=
E<sub>FOS</sub>·G<sub>ATM</sub>·S<sub>OCEAN</sub>; see methods and ref. 8) which grew from 0.2 ± 0.4
(mean ± std. dev) PgC yr<sup>-1</sup> in the 1960s to 1.7 ± 0.6 PgC yr<sup>-1</sup> in the
decade 2011–2020, with DGVMs (S3 simulation; see Methods) cap-
turing the increase (–0.1 ± 0.6 to 1.6 ± 0.5 PgC yr<sup>-1</sup>; Fig. 1). DGVMs
suggest the enhanced net sink over the past 60 years has mainly been
driven by rising atmospheric CO<sub>2</sub> concentrations and nitrogen
deposition (1.2 ± 0.2 PgC yr<sup>-1</sup> over 1960–1969 rising to 3.5 ± 0.8 PgC yr<sup>-1</sup>
over 2011–2020) (Fig. 1b), with the sink partially offset by relatively
constant net LULCC emissions of 1.3 ± 0.5 PgC yr<sup>-1</sup> (predominantly
arising from tropical deforestation and shifting cultivation<sup>12</sup>; Fig. 2 and
Supplementary Fig. 1). DGVMs indicate that climatic variability drives
the large year-to-year changes (±2 PgC yr<sup>-1</sup>) in the land sink (Fig. 1), with
long-term (multi-decadal) climate trends reducing the land sink since
the 1980s (–0.4 ± 0.5 PgC yr<sup>-1</sup> over 1980–2020) (Fig. 1b).

The CO₂ and N deposition-driven land carbon sink is located in northern and tropical forests (Fig. 2). LULCC-induced carbon losses occur across the globe but are most apparent in tropical latitudes and certain hotspot regions across the northern hemisphere (China, USA and West Eurasia), which recently have been (within our study period)

<figure><img src="image_4.png" /><figcaption>**Fig. 1 | Global net land sink and attribution to drivers.** a Net annual land carbon sink (PgC yr⁻¹) as estimated by dynamic global vegetation models (DGVMs) with all drivers varying (red) and the top-down Global Carbon Budget (GCB) constraint (black) and b the decomposition of the DGVM net sink into contributions from</figcaption></figure>

rising atmospheric CO₂ concentrations and N deposition (blue), changes in climate (yellow), and land-use and land cover change (grey). Thick lines show a locally fitted regression. Shading around DGVM estimates corresponds to 1σ, and the uncertainty on the GCB constraint is taken from ref. 8.
<figure><img src="image_5.png" /></figure>

**Fig. 2 | Driver attribution to spatial changes in total ecosystem carbon.** Maps show the multi-model mean **a** net change in ecosystem (vegetation and soil) carbon (kgC m⁻²) from 1959–2020 and **b** the contribution of each driver to overall change; CO₂ and N deposition (green), climate (red), and land-use and land cover change

(blue). Stippling in panel **a** indicates <80% of models agree on the direction of change. The colours in panel **b** are calculated by assigning a red-green-blue (RGB) value to each grid depending on the relative magnitude of change due to each driver. Transparency is determined by the magnitude of the net change in panel **a**.

densely forested (Fig. 2 and Supplementary Fig. 1). There are signs of European carbon sink in part driven by a/reforestation or natural land (re)establishment. The impact of changes in climate can also be detected across the globe, with vegetation and soil carbon losses in the Amazon, the Sahel and South Africa, and climate-induced carbon gains in east Brazil, Australia, and across the high northern latitudes (Fig. 2 and Supplementary Fig. 1). In summary, the DGVM ensemble agrees on the large sink located in the world's forests. However, the models are not in full agreement in the direction of change in ecosystem carbon across much of the globe, in particular in regions with competing CO₂ and LULCC effects (Fig. 2a).

## Key processes and uncertainties behind the sink

The top-down global GCB budget constraint combined with the bottom-up DGVM estimates give high confidence the land has been a net sink of carbon over the last 60 years, however, the relative contribution of each driver (e.g. CO₂, N deposition, climate, LULCC), and ecological attribution (vegetation, soil) of the sink remains elusive.

The DGVM multi-model mean suggests a net global increase in both vegetation ($ΔC_v = 28 \pm 26$ (mean ± std. dev) PgC) and soil ($ΔC_s = 21 \pm 32$ PgC) carbon stocks, although the uncertainty is large (Fig. 3a, b), with two models (LPX-Bern and YIBs) simulating a net loss of vegetation carbon and three models (CLASSIC-N, ISAM and LPJ-GUESS) simulating a net loss of soil carbon (Fig. 3c).

This relatively high uncertainty is in part due to large opposing fluxes (Fig. 3a) driven by increasing atmospheric CO₂ and N deposition (carbon increase) vs LULCC (carbon decrease at global scale), with soil carbon gains (from rising CO₂) also counteracted by negative fluxes driven by climate change (Fig. 3b). In general, there are no observational constraints on long-term (pre-satellite era) changes in global vegetation or soil carbon, and so reducing model uncertainty is challenging. Therefore, although the DGVMs generally capture the global net carbon sink, the attribution to vegetation or soil stocks is

completely unconstrained (Fig. 3c) and is driven by alternate model structures and parameterisations between DGVMs.

DGVMs have inherently different baseline productivity, allocation, and turnover rates as well as varying degrees of sensitivity of processes to environmental change. To gain a deeper understanding of the processes driving carbon sink changes, we use our attribution framework (see Methods) to associate the modelled changes in vegetation and soil carbon to changes in inputs (NPP for vegetation and litterfall for soil) and outputs (turnover rates for vegetation and soil). For each of these processes, we identify their changes due to CO₂, climate change or LULCC, according to the models.

The ensemble-simulated increase in global vegetation carbon over the last 60 years ($ΔC_v = 28 \pm 26$ PgC; Fig. 4a) was driven by enhanced plant productivity in response to rising CO₂ concentrations and N deposition, with small (but highly uncertain) losses due to changes in vegetation turnover ($ΔC_{v,input,CO2} = 87 \pm 26$ PgC, $ΔC_{v output,CO2} = -3 \pm 21$ PgC; Fig. 4a). Although the central estimate of the net CO₂ and nitrogen deposition biomass response appears relatively well constrained at global scale ($ΔC_{v,CO2} = 84 \pm 31$ PgC), model estimates range from a gain of 30 PgC (YIBs) to a gain of 150 PgC (LPJ) in biomass due to rising atmospheric CO₂ (Supplementary Fig. 2).

Model spread is driven by a combination of uncertainty in the response of both NPP (e.g. driven by uncertainty in leaf-level photosynthetic response, scaling to canopy and landscape scales, carbon allocation and nutrient limitation) and vegetation turnover ($τ_v$) (e.g. driven by the treatment of nutrient limitation coupled with differences in allocation and in some cases forest-stand packing constraints¹³) to rising atmospheric CO₂ and N deposition (Fig. 5 and Supplementary Fig. 3). There is some confidence in the direction of change in global and regional NPP, albeit not in the magnitude, whereas models disagree on the direction of $Δτ_v$ due to rising CO₂ (Fig. 4). CO₂ and N deposition-driven $Δτ_v$ depends on changes in stand dynamics (in some models) and plant allocation. For example, an increase in short-lived
<figure><img src="image_6.png" /></figure>

**Fig. 3 | Temporal changes in global vegetation and soil carbon stocks.** Time-series show the change in global a vegetation (ΔC<sub>v</sub>) and b soil (ΔC<sub>s</sub>) carbon stocks from 1959–2020 due to each of the three external drivers (CO<sub>2</sub> and N deposition, climate, land-use, and land cover change). Lines represent the mean of the dynamic global vegetation models (DGVMs) and shading the ±1σ of the DGVMs. The DGVM

output is first smoothed using a fourth-order spline. The cumulative net (sum of three drivers) change in global carbon stocks by 2020 is shown (red crosses show each model and red circle shows the model mean). c Shows the change in vegetation and soil stocks for each of the 18 models and the grey region is the Global Carbon Budget net land sink constraint (see Methods).

root production, at the expense of wood growth, to alleviate nutrient limitations on plant growth can reduce τ<sub>v</sub>. This increase in root allocation also changes vegetation nitrogen demand (because of the substantially higher C:N in wood compared to roots), which feeds back onto the whole-plant productivity response to rising CO<sub>2</sub>. Conversely, some models increase their wood allocation fraction when NPP (or production) increases, which will increase τ<sub>v</sub><sup>14</sup>.

This increase in vegetation stocks has been partially offset by emissions associated with LULCC from vegetation (ΔC<sub>v,LULCC</sub> = -58 ± 21 PgC; Fig. 4a), predominantly in the tropics (Supplementary Fig. 1), which increases vegetation turnover (loss of biomass to the atmosphere or wood products). About half (ΔĈ<sub>v输入,LULCC</sub> = -23 ± 22 PgC, -40%; Fig. 4a) of the net LULCC-driven global vegetation carbon losses are from a reduction in inputs (NPP) to the land, as crops or pastures often have lower productivity than the forests that they replace. The other half (ΔC<sub>v, output,LULCC</sub> = -29 ± 18 PgC; Fig. 4a) is attributable to changes in turnover related to forests containing woody biomass, which has a slower turnover than leaves and fine roots. The turnover effect is further aggravated through the phenomenon that due to their longer-lived biomass, forests take up more carbon under rising CO2 levels, but this sink is lost by clearing for agricultural use. This is known as the “loss of additional sink capacity” (LASC)<sup>15,16</sup> and accounts for the impact of environmental (CO<sub>2</sub> and climate) changes on carbon uptake of deforested land (found in S3 simulation) compared to potential vegetation (found in S2 simulation). Ref. 16 attributes ~40% of DGVM estimated LULCC vegetation and soil carbon losses to the LASC, albeit with a different methodology and study period.

All models simulate a net loss of global biomass following LULCC, although the magnitude is not well constrained. Model spread is driven

by differences in LULCC process representation in models<sup>7</sup>. We find that models that include wood harvest (routine harvest of established managed forests), grazing harvest, or shifting cultivation simulate larger C<sub>v</sub> losses than models without these land management processes (Supplementary Fig. 4), in line with earlier studies showing the high importance of land management on vegetation carbon stocks<sup>17</sup>. Further, in northern ecosystems, models do not agree on the direction of ΔNPP or Δτ<sub>v</sub> following LULCC (Supplementary Fig. 5a). In addition to the land management processes mentioned above, uncertainty in large-scale regional changes in C<sub>v</sub> is in part driven by smaller-scale regional compensations, where losses in Russia and USA are somewhat countered by a forest regrowth sink in Eurasia (Supplementary Fig. 1). Therefore, the balance between these two opposing fluxes determines net northern changes in C<sub>v</sub> from LULCC, with uncertainties in both carbon loss fluxes and regrowth uptake<sup>18,19</sup>. In addition, spatial patterns of northern LULCC losses and gains are not entirely consistent between models (Supplementary Fig. 6), which adds to uncertainty in local, regional, and hemispheric net LULCC-driven changes in vegetation carbon.

Changes in climate have a much smaller but similarly uncertain influence on vegetation processes at a global scale, with no agreement among models on the direction of change in biomass driven by shifts in NPP (ΔĈ<sub>v输入,CLIM</sub> = -1 ± 14 PgC over 1959–2020) or vegetation turnover (ΔC<sub>v, output,CLIM</sub> = 3 ± 9 PgC) at global scale (Fig. 4a). However, models align closer at regional scales, where northern warming has stimulated productivity (ΔĈ<sub>v, input,CLIM,North</sub> = 11 ± 6 PgC; Supplementary Fig. 5a), and in the tropics, reductions in productivity (due to warming and/or changes in precipitation) lead to losses of vegetation (ΔĈ<sub>v, input,CLIM,Tropics</sub> = -12 ± 12 PgC; Supplementary Fig. 5b). Overall,
<figure><img src="image_7.png" /><figcaption>**Fig. 4** | Process and driver attribution of changes in global vegetation and soil carbon stocks. Change in global a vegetation (ΔC<sub>v</sub>) and b soil (ΔC<sub>s</sub>) carbon stocks over 1959–2020 (PgC). The contribution to net changes in carbon stocks (green bars) from changes in inputs (net primary productivity for vegetation (ΔNPPr<sub>v,1959</sub>) and vegetation to soil flux for soil (Δf<sub>vs</sub>τ<sub>s,1959</sub>), red bars), outputs/turnover (NPP<sub>1959</sub>Δτ<sub>v</sub> for vegetation and f<sub>vs,1959</sub>Δτ<sub>s</sub> for soil, orange bars), and the interaction</figcaption></figure>

term (ΔNPPΔτ<sub>v</sub> for vegetation and Δf<sub>vs</sub>Δτ<sub>s</sub> for soil, blue bars) are shown. The bars depict the multi-model mean with the range as ±1σ of the models. The arrows show the direction of change in carbon stocks due to each process. The panels from left to right show the changes due to all drivers varying (ALL), changes in atmospheric CO<sub>2</sub> and N deposition, climate (CLIM), and land-use and land cover change (LULCC).

<figure><img src="image_8.png" /><figcaption>**Fig. 5** | Attribution of uncertainty to processes in modelled changes in carbon stocks. The relative uncertainty (defined as the standard deviation among model estimates) in the change in global vegetation (ΔC<sub>v</sub>) and soil (ΔC<sub>s</sub>) carbon stocks resulting from each of the driving terms (Eqs. 9 and 10 in Methods). The four terms are baseline input (baseline productivity NPP<sub>1959</sub>) for vegetation and baseline litterfall/mortality (f<sub>vs,1959</sub>) for soil), change in inputs (change in productivity ΔNPP for vegetation and change in litterfall/mortality (Δf<sub>vs</sub>) for soil), baseline turnover (τ<sub>v,1959</sub> for vegetation and τ<sub>s,1959</sub> for soil), and change in turnover (Δτ<sub>v</sub> for vegetation and Δτ<sub>s</sub> for soil). For each of these terms in Eqs. 9 and 10, we calculate the standard deviation in ΔC<sub>v</sub> and ΔC<sub>s</sub> using the multi-model mean values of all other terms in the equations and the individual model values for that term.</figcaption></figure>

these opposing regional trends cancel and explain the small global signal and the model uncertainty at the global scale.

Increased production and loss of biomass enhance soil inputs (ΔĈ<sub>s,Input,ALL</sub> PgC) and drives the growth in global soil carbon stocks (ΔC<sub>s</sub> PgC), albeit with low confidence in the magnitude (Fig. 4b). Our analysis also indicates that enhanced CO<sub>2</sub> concentrations increase soil turnover, which is a direct result of enhanced inputs to the fast-turnover litter and surface soil pools. Our methodology simplifies the TRENDY model structures into a single soil pool to represent the entire soil system, and so a relative increase in the *faster* pools moves the aggregate soil pool turnover (τ<sub>s</sub>) towards that of the litter and surface

soil pools—a phenomenon known as false-priming<sup>20</sup>. In addition, changes to soil conditions (e.g. increased soil moisture<sup>13</sup>) following increased atmospheric CO<sub>2</sub> can impact τ<sub>s</sub>, although it is difficult to separate the false-priming and actual changes in τ<sub>s</sub> with the set of simulations used in this study. This simplification leads to low confidence in the partition of net soil carbon changes into the input and turnover-driven changes, whereas the net change is well constrained (ΔC<sub>s,CO2</sub>PgC).

Reduced litterfall (input to soil) due to the conversion of forests to agricultural land is the predominant pathway of soil carbon loss, in line with previous studies<sup>21</sup> (ΔĈ<sub>s,Input,LULCC</sub> = -36 ± 76 PgC, ΔC<sub>s,LULCC</sub> = -25 ± 29 PgC). This is particularly evident in the tropics (Supplementary Fig. 7b) with large-scale deforestation over the last several decades, while deforestation and reforestation are more in balance in the extra tropics. However, there is relatively high uncertainty in the magnitude of LULCC impacts on soil carbon inputs (ΔĈ<sub>s,Input,LULCC</sub>) (Fig. 4b), with models ranging from losses of over 100 PgC to gains of 45 PgC over the period 1959–2020. Land management processes such as crop and wood harvesting (in particular, the frequency of harvest and fraction of biomass removed and respired elsewhere) are leading uncertainties in simulating LULCC impacts on soil carbon (Supplementary Fig. 4)<sup>22,23</sup>. Turnover-driven changes play less of a role in soil C changes in our models. This could be due to offsetting effects such as a decrease in soil C due to the altered quality of the litter input (a higher fraction of faster decomposing material for agriculture as compared to forests) or a reduction in fire-related losses with the transformation of natural ecosystems<sup>21</sup>.

Increased heterotrophic respiration rates and subsequent reductions in the turnover time of the soil pool due to global warming entirely offset the increase in soil carbon from climate-enhanced NPP and litterfall (Fig. 4b). Overall, changes in climate caused a net loss (ΔC<sub>s,CLIM</sub> = -13 ± 12 PgC) of global soil carbon. We find the impact of climate change on carbon storage is accelerating, in particular for northern soil and tropical biomass, where carbon stock changes from increased/decreased productivity (for northern soil and tropical
biomass, respectively) have accelerated by 19 TgC yr⁻² and 9 TgC yr⁻², since the turn of this century (Supplementary Fig. 8). This equates to annual mean carbon gain/losses due to changes in productivity of 1.6 PgC yr⁻¹ and -0.3 PgC yr⁻¹ over 2001–2020 compared to 0.6 PgC yr⁻¹ and -0.2 PgC yr⁻¹ over 1961–2000, for northern soil and tropical biomass, respectively (Supplementary Fig. 8b, e). This acceleration is likely driven by warming in the north (2001–2020 April-August mean 0.9 °C warmer than 1961–2000) and the tropics (2001–2020 annual mean 0.6 °C warmer than 1961–2000)²⁴,²⁵, as well as potential changes in sensitivity of the biosphere to warming²⁶.

## Discussion

The 'observed' GCB global net land sink is captured by the TRENDYv10 DGVM ensemble. We find a robust agreement that rising atmospheric CO₂ concentrations and N deposition drive the land carbon sink, and changes in climate along with LULCC lead to a source of carbon in the atmosphere. Further, the DGVMs corroborate that the net land sink is located in northern latitudes, with the tropics more carbon neutral, in line with recent independent top-down estimates³,⁸.

However, the models do not entirely agree on the partition of the net land sink between vegetation and soil. Thirteen of 18 models indicate an increase in global vegetation and soil carbon stocks, although the magnitude of change is highly unconstrained (ranging from -0 to +80 PgC over the past 60 years). The remaining five models simulate a reduction in either vegetation or soil stocks, highlighting the disparity between models regarding internal carbon cycling processes.

There are no observations of global carbon stocks covering the entire study period, but a recent synthesis of observation-based biomass change estimates in global forests suggests a net biomass sink of 0.3–2.1 PgC yr⁻¹ since 2000²⁷, giving confidence to the models that simulate net gains. However, these observational studies do not constrain the magnitude of the vegetation sink any more than the DGVM mean biomass sink of -0.1 to 1.7 PgC yr⁻¹ (averaged over 2001–2020). Furthermore, we cannot even be certain in the direction of change in global soil carbon as direct observations of large-scale changes in soil carbon stocks do not exist²,⁴.

Our analysis indicates that baseline vegetation turnover rates and the turnover response to rising atmospheric CO₂ and N deposition are the key uncertainties in modelled Δ*C*<sub>v</sub>, with a smaller but still significant contribution from changes in NPP. Differences in modelled vegetation turnover are driven by variations in whole-plant mortality rates, simulated tissue lifespan, and the allocation of NPP to plant components with inherently different turnover rates and their responses to changing environmental conditions¹⁴,²⁸. For example, some models may reduce plant mortality and turnover as atmospheric CO₂ concentrations rise. Increased carbon stores could supply maintenance respiration in periods of photosynthetic stress²⁹,³⁰, and increases in water-use efficiency may alleviate the impact of drought³¹,³². In contrast, shifts in plant community composition or self-thinning dynamics may increase stand-level mortality³³. These structural differences lead to no agreement on the direction of Δτ<sub>v</sub> in the models presented here.

In addition to turnover-related uncertainties, the increase in NPP and biomass due to rising CO₂ is not well constrained. The choice of photosynthesis and stomatal model determines the leaf-level response, with a wide range of implementations and outcomes in DGVMs³⁴. Scaling from leaf to canopy level causes further discrepancies due to different structural assumptions between models as to the vertical distribution of light, nitrogen, photosynthetic capacity, and the treatment of sunlit and shaded leaves³⁴. Further, NPP and biomass production is constrained by stoichiometric nutrient requirements. Nitrogen availability mediates the biomass response to rising CO₂³⁵, and therefore it is likely that differences in the inclusion and representation of carbon-nitrogen coupling between models

contribute to uncertainty in Δ*C*<sub>v</sub><sup>36,37</sup>. Overall, these structural uncertainties manifest as a four-fold range in estimates of the CO₂ effect on NPP and biomass (gains of 41 PgC (CLASSIC-N) up to 157 PgC (ISAM) over the previous 60 years).

Furthermore, LULCC-driven reductions in NPP and woody carbon stocks have removed a substantial portion of global vegetation carbon. However, the exact magnitude of carbon loss is not well constrained by the DGVMs. Biases in the carbon density of deforested land³⁸ and model-specific choices of whether to replace forest or grassland for agricultural expansion are leading sources of error when calculating LULCC losses³⁹. Our results highlight that, firstly, variation in the representation of relevant processes, in particular, forest management (e.g. treatment of wood harvest) and shifting cultivation, leads to systematic differences between modelled Δ*C*<sub>v</sub><sup>7,40</sup>. Second, DGVMs are missing carbon losses associated with degraded forests, which can exceed deforestation losses, indicating simulated LULCC losses could be underestimated⁴⁰⁻⁴². Finally, uncertainty exists with LULCC maps used to drive DGVMs as historical land-use is not perfectly known and products differ in the land cover types and transitions included⁴⁰.

The models suggest that soil has been a net carbon sink globally over the past six decades, although there is over 100% uncertainty on the magnitude (21 ± 32 PgC). This relatively large range is due to the opposing impacts of rising CO₂ and climate/LULCC effects which partly cancel out. Further, the impact of LULCC on soil carbon is also difficult to quantify robustly with the suite of DGVMs used here, as models do not agree on the sign of Δf<sub>vs</sub> or Δτ<sub>s</sub> following LULCC. In general, changes in soil carbon stocks are difficult to constrain due to a lack of large-scale observations. It is likely historical LULCC caused additional carbon fluxes as a result of land erosion⁴³,⁴⁴, degradation of agricultural soils⁴⁵,⁴⁶, and losses from drained peatlands⁴⁷. These management processes are not included in models and most do not simulate peatlands at all, potentially causing errors in simulated soil carbon fluxes.

Overall, the soil carbon sink is caused by increased litter inputs from CO₂-driven vegetation growth. There is some evidence additional litter and root exudates enhance soil carbon stocks when nitrogen availability is high⁴⁸. It is important to stress that the model simulations include widespread nitrogen deposition, which may have helped to sustain a strong CO₂ response of biomass and soil carbon sinks⁵,⁴⁹.

However, additional soil inputs can accelerate organic matter decomposition to release plant-available nitrogen (via priming effects⁸⁰), which reduces soil carbon storage⁵¹⁻⁵³. DGVMs do not simulate priming effects or explicitly account for microbial activity, mineral association, aggregation, or mycorrhizal fungi interactions, all of which regulate soil carbon turnover and accumulation rates⁵⁴. DGVMs represent soil carbon with a set of cascading pools, where decomposition losses are determined by first-order decay rates (although this is an active area of model development⁵⁵). Therefore, the modelled soil sink has to be interpreted in the context of model structure and the limited ability of DGVMs to capture governing soil processes.

In general, it is difficult to attribute processes (soil inputs vs turnover) to changes in soil carbon for two main reasons. First, soil model structure leads to false-priming effects²⁰. The apparent reduction in soil turnover is a consequence of increased surface soil carbon relative to deep soil, as opposed to a change in actual turnover rates. Second, our process attribution methodology treats the soil as a single pool. Therefore, the calculation of τ<sub>s</sub> (C<sub>s</sub>/R<sub>s</sub>) may include deep, inactive (on decadal timescales) soil carbon, and underestimate the turnover rate of the active soil. Hence, our process attribution of soil carbon changes is not well constrained. For example, we estimate a 300 PgC increase in soil carbon from increased inputs for CLM5.0, a result which is in part driven by the large carbon stocks in high northern latitudes leading to a large estimate of τ<sub>s</sub>⁵⁶. We do account for this bias by scaling output to actual changes in modelled carbon stocks (see
methodology), but our results highlight the urgent need for evaluation of individual carbon pools and fluxes (leaf, root, wood and each soil pool) separately, and not bulk terms.

In addition to the CO₂-driven sink, evidence exists for a substantial impact of forest regrowth on the northern sink over the previous 60 years¹⁸,¹⁹. DGVMs have some capacity to capture regrowth effects (predominantly following land abandonment), however, most do not simulate forest demography or detailed forest management (e.g. enhanced stock densities⁵⁷) and so could underestimate the actual regrowth carbon sink³⁸. Further, not all DGVMs simulate disturbance in unmanaged land (e.g. wind and pests), and therefore, natural disturbances beyond fire and the subsequent regrowth will not be captured by the DGVMs here¹⁸. Moreover, not all models simulate fire mortality, and fire models do not always capture observed spatiotemporal patterns and associated carbon emissions⁵⁹, as well as large uncertainties due to issues with forcing data⁶⁰.

In addition to direct human-caused carbon losses, alternative drivers of mortality can reduce carbon stocks. For example, climate-induced mortality, combined with a shorter tree lifespan resulting from faster growth⁶¹,⁶², currently weaken the Amazon forest sink⁶³. However, DGVMs show no sign of increased carbon turnover due to changes in climate or faster growth, although this discrepancy is not unexpected as DGVMs do not include detailed mortality processes, e.g. drought-mortality⁶⁴, although this is an area of current model development⁶⁵,⁶⁶. The DGVMs do indicate a decline in tropical productivity due to climate change, in line with remote-sensing and upscaled in-situ observations⁶⁷, indicating current temperatures may exceed tropical plant thresholds²⁴.

Here we have shown that the long-term evolution of the global net land sink is well estimated by an ensemble of state-of-the-art DGVMs. However, a complete process-based understanding is still lacking due to several model shortcomings. Specifically, above/below-ground partitioning of the sink is highly uncertain due to too simplistic representation of internal carbon cycling. Model improvement into plant allocation, tissue lifespan and mortality, as well as the inclusion of process-based soil carbon and nutrient cycling, should be a high priority moving forward.

## Methods

### TRENDYv10 models

We analyse output from 18 DGVMs that are part of a recent model intercomparison project, TRENDYv10 and the Global Carbon Budget 2021 (GCB2021)⁸. The models included in the analysis here are CABLE-POP, CLASSIC, CLASSIC-N, CLM5.0, DLEM, IBIS, ISAM, ISBA-CTRIP, JSBACH, LPJ-GUESS, LPJ, LPX-Bern, OCN, ORCHIDEE, ORCHIDEEv3, SDGVM, VISIT, and YIBs (see ref. 8 for model descriptions and setup). One model (JULES-ES-1.1) from TRENDYv10 is not included due to incomplete data. Note, CLASSIC-N⁶⁸ and ORCHIDEE⁶⁹ were not part of GCB2021 due to the inclusion of alternate model versions (CLASSIC and ORCHIDEEv3) but are included in TRENDYv10 and this study.

The models are forced with a merged monthly Climate Research Unit (CRU)⁷⁰ and 6-hourly Japanese 55-year Reanalysis (JRA-55)⁷¹ data set. The models are also forced with atmospheric CO₂⁷², gridded nitrogen deposition⁷³ and nitrogen fertiliser⁷⁴. DGVMs use the HYDE (v3.3) land-use change data set⁷⁵,⁷⁶, which provides annual pasture and cropland areas at a global scale, and includes improvements in the spatial distribution of agricultural regions⁷⁷. Several models (CABLE-POP, CLM5.0, JSBACH, LPJ-GUESS, LPJ, and VISIT) also use harmonised land-use change data (LUH2-GCB2021⁸), which provides information on sub-grid-scale land-use transitions.

To isolate the response of the land to each driver (CO₂, climate, LULCC), each model performs four simulations: S0 (fixed pre-industrial atmospheric CO₂ and land-use, recycled 1901–1920 climate), S1 (transient atmospheric CO₂, recycled 1901–1920 climate, and fixed pre-industrial land-use), S2 (transient atmospheric CO₂, transient

climate, and fixed pre-industrial land-use), and S3 (transient atmospheric CO₂, transient climate, and transient industrial land-use). Nitrogen deposition varies temporally in simulations S1–S3. Therefore, the transient CO₂ (+N deposition) effect on the ‘natural’ land sink is calculated by S1–S0, S2–S1+S0 is the climate effect on the ‘natural’ land sink, S3–S2 is the LULCC effect, and S3 is the net effect. There exists an artefact in the HYDE3.3 data causing a large land-use transition and emission peak around 1960. To correct this, we replace the LULCC estimates for 1959–1961 with the average of 1958 and 1962 in each DGVM.

Net Biome Productivity (NBP) from the S3 simulation represents the net land sink and can therefore be compared to the ‘observed’ Global Carbon Budget land constraint, which is calculated in the GCB as fossil fuel emissions–atmospheric carbon growth rate–ocean carbon sink (see Fig. 1 in the main text and Table 5 in ref. 8 and data taken from https://doi.org/10.18160/gcp-2021). Atmospheric CO₂ concentration measurements began in 1959⁷⁸, and so this independent constraint on the land sink covers the period 1959–2020, which defines the study period used throughout our analysis.

### Data processing

First, we calculate for each model the global and regional (two regions: north of 30°N and south of 30°N) annual mean values for NBP, net primary productivity (NPP), heterotrophic respiration ($R_h$), vegetation carbon ($C_v$) and soil carbon ($C_s$), which includes litter ($C_{litter}$) and coarse woody debris ($C_{CWD}$) for 1959–2020. We exclude the product pool from our analysis as we do not have carbon emission data from wood products available.

### Process attribution framework

To attribute changes in C pools to different processes, we quantify C pool changes attributable to changes in productivity and changes in turnover time. We do this using an analytical approximation of C pool dynamics corrected for non-steady state and more complicated behaviour of the models. First, we approximate the steady-state C pool size given the inputs and turnover time of a given pool in a given year. We then quantify the change in the steady-state pool sizes estimated for 1959 and the year in question and partition this change proportionally among inputs and turnover time. Finally, to adjust to non-steady state conditions and more complicated model structures (i.e. multiple pool dynamics, not a single pool as in this approximation), we adjust the attributed proportions by adding the difference in actual modelled change in carbon stock and the steady-state approximation, scaled by the relative magnitude of productivity and turnover-driven changes.

The change in carbon pools can be expressed as:

$$ \frac{dC}{dt} = \text{INPUT} - \text{OUTPUT}, \qquad (1) $$

where we define $\frac{dC}{dt}$ as the change in carbon pool per year. To quantify the internal processes impacting the output flux, we can define $OUTPUT = \frac{C}{\tau}$, where $\tau$ is the mean turnover time (the inverse of the mean turnover rate) of the carbon pool, C. Therefore, we can rewrite Eq. 1 as:

$$ \frac{dC}{dt} = \text{INPUT} - \frac{C}{\tau} \qquad (2) $$

We separate the land into two main pools, vegetation ($C_v$) and soil ($C_s$). $C_v$ is a direct model output and includes leaf, root and wood biomass. Similarly, $C_{litter}$, $C_{CWD}$ are direct model outputs, which we sum together for our analysis. The rest of the methodology, $C_s$ refers to the sum of these three pools. Not all models simulate litter and CWD pools. ISAM does not simulate litter and only CABLE-POP, CLM5.0,
DLEM and LPJ-GUESS simulate CWD. Further, for CLM5.0, the litter pool is incorporated in the soil pool output, and so is also excluded from the sum to avoid double counting.

We can now specify Eq. 2 for vegetation and soil as:

$$ \frac{dC_v}{dt} = \text{NPP} - \frac{C_v}{\tau_v} \qquad (3) $$

$$ \frac{dC_s}{dt} = f_{vs} - \frac{C_s}{\tau_s}, \qquad (4) $$

where $f_{vs}$ is the carbon flux from vegetation to soil, from either litterfall, mortality or direct transfer from roots to soil.

We calculate annual $\tau_v$ (using Eq. 3), $\tau_s$ and $f_{vs}$ using model output as follows:

$$ \tau_v = \frac{C_v}{\text{NPP} - \Delta C_v}, \\ \tau_s = \frac{C_s}{R_h}, $$

$$ f_{vs} = \Delta C_s + R_h, $$

where $\Delta C_v$ and $\Delta C_s$ are the annual changes in vegetation and soil carbon. We want to be able to define $\Delta C$ in terms of input fluxes and turnover times only, as these are the two aggregated mechanisms driving changes in the pools. To do this, we use the fact that the pools and the input and output fluxes are generally much larger than the change in carbon pools each year to make the assumption that $\Delta C \simeq 0$. Rearranging Eqs. 3 and 4, we can estimate the steady-state carbon pools in any particular year as:

$$ C_v \approx \hat{C}_v = \text{NPP}\tau_v \qquad (5) $$

$$ C_s \approx \hat{C}_s = f_{vs}\tau_s, \qquad (6) $$

where $\hat{C}_v$ and $\hat{C}_s$ are steady-state approximations of the actual carbon pools $C_v$ and $C_s$. We can now define the change in carbon pools over our 60-year study period, 1959–2020, in terms of inputs and turnover time:

$$ \Delta C_v \approx \Delta \hat{C}_v = \text{NPP}_{2020}\tau_{v,2020} - \text{NPP}_{1959}\tau_{v,1959} \qquad (7) $$

$$ \Delta C_s \approx \Delta \hat{C}_s = f_{vs,2020}\tau_{s,2020} - f_{vs,1959}\tau_{s,1959}. \qquad (8) $$

In order to separate the impacts of changes in inputs and turnover times, we next define NPP$_{2020}$ = NPP$_{1959}$ + $\Delta$NPP, $f_{vs,2020} = f_{vs,1959} + \Delta f_{vs}$, and $\tau_{2020} = \tau_{1959} + \Delta \tau$. We now substitute these into Eqs. 7 and 8:

$$ \Delta \hat{C}_v = \Delta \text{NPP} \tau_{v,1959} + \text{NPP}_{1959} \Delta \tau_v + \Delta \text{NPP} \Delta \tau_v \qquad (9) $$

$$ \Delta \hat{C}_s = \Delta f_{vs} \tau_{s,1959} + f_{vs,1959} \Delta \tau_s + \Delta f_{vs} \Delta \tau_s \qquad (10) $$

The right-hand side of the equation contains three terms corresponding to the change in carbon storage due to changes in inputs (NPP for vegetation and litterfall for soil), outputs (bulk turnover), and an interaction term. Our approach, therefore, extends the factor separation approach by ref. 21, which applied it to attribute simulated soil carbon changes into input-driven, turnover-driven change and a synergy term to also cover vegetation carbon changes.

The change in vegetation carbon due to each term can be written as: $\Delta \hat{C}_{v,\text{input}} = \Delta \text{NPP} \tau_{v,1959}$, $\Delta \hat{C}_{v,\text{output}} = \text{NPP}_{1959} \Delta \tau_v$, $\Delta \hat{C}_{v,\text{interaction}} = \Delta \text{NPP} \Delta \tau$ and $\Delta \hat{C}_v = \Delta \hat{C}_{v,\text{input}} + \Delta \hat{C}_{v,\text{output}} + \Delta \hat{C}_{v,\text{interaction}}$, with equivalent terms for soil carbon.

Finally, we adjust the one-pool steady-state estimates for the difference in the model simulations that arise due to the combination of steady-state assumption (Eqs. 5–8), the grouping of all pools into single vegetation and soil pools (as the distribution between component pools changes over time), and the linearisation not capturing the full dynamics of land carbon cycling (Eqs. 9 and 10). Therefore, the estimates of changes in carbon pools ($\Delta \hat{C}_v$ and $\Delta \hat{C}_s$) do not equal the actual change in C pools ($\Delta C_v$ and $\Delta C_s$), which are a direct model output (see Supplementary Fig. 9).

To address this issue, we calculate for each model the difference between the actual change in carbon stocks and our approximation ($\delta_v = \Delta C_v - \Delta \hat{C}_v$ and $\delta_s = \Delta C_s - \Delta \hat{C}_s$). We make the simple assumption that productivity and turnover are proportionally responsible for the mismatch between steady state and the actual carbon stock, and so adjust each term on the right-hand side of Eqs. 9 and 10 by the relative magnitude of each term multiplied by $\delta_v$ or $\delta_s$ for vegetation and soil, respectively. For example, $\Delta \hat{C}_{v,\text{input}}$ is adjusted by adding:

$$ \frac{|\Delta \hat{C}_{v,\text{input}}|}{|\Delta \hat{C}_{v,\text{input}}| + |\Delta \hat{C}_{v,\text{output}}| + |\Delta \hat{C}_{v,\text{interaction}}|} * \delta_v. \qquad (11) $$

We calculate Eqs. 9 and 10 for each of the eighteen models, apply the appropriate scaling and use them to produce Fig. 4.

## Uncertainty analysis

The spread of model estimates of changes in carbon stocks and fluxes is used to determine the uncertainty of our results. For land flux estimates (NBP) and $\Delta C_v$ and $\Delta C_s$, we show ‘mean ± std. dev.’. We decompose the uncertainty (defined as the model spread) in $\Delta \hat{C}_v$ and $\Delta \hat{C}_s$ into four components: baseline input (NPP$_{1959}$ or $f_{vs,1959}$), change in inputs ($\Delta$NPP or $\Delta f_{vs}$), baseline turnover ($\tau_{v,1959}$ or $\tau_{s,1959}$), and change in turnover ($\Delta \tau_v$ or $\Delta \tau_s$). For each of these terms in Eqs. 9 and 10, we calculate the standard deviation in $\Delta \hat{C}_v$ and $\Delta \hat{C}_s$ using the multi-model mean values of all other terms in the equations and the individual model values for that term.

## Data availability

All code and post-processed data generated in this study are available at https://doi.org/10.5281/zenodo.6884342. The raw model output is available at the following sftp site: trendy-v10@trendy.ex.ac.uk. Access will be granted by first contacting Stephen Sitch (S.A.Sitch@exeter.ac.uk).

## References

1. Friedlingstein, P. et al. Climate–carbon cycle feedback analysis: Results from the C4MIP model intercomparison. *J. Clim.* **19**, 3337–3353 (2006).

2. Canadell, J. G., Monteiro, P. M. S., Costa, M. H. & Da Cunha, L. C. Global carbon and other biogeochemical cycles and feedbacks. (2021).

3. Ciais, P. et al. Five decades of northern land carbon uptake revealed by the interhemispheric CO gradient. *Nature* **568**, 221–225 (2019).

4. Walker, A. P. et al. Integrating the evidence for a terrestrial carbon sink caused by increasing atmospheric CO2. *N. Phytol.* **229**, 2413–2445 (2021).

5. O’Sullivan, M. et al. Have synergies between nitrogen deposition and atmospheric CO₂ driven the recent enhancement of the terrestrial carbon sink? *Glob. Biogeochem. Cycles* **33**, 163–180 (2019).

6. Forkel, M. et al. Enhanced seasonal CO₂ exchange caused by amplified plant productivity in northern ecosystems. *Science* **351**, 696–699 (2016).
7. Arneth, A. et al. Historical carbon dioxide emissions caused by land-use changes are possibly larger than assumed. *Nat. Geosci.* **10**, 79–84 (2017).

8. Friedlingstein, P. et al. Global carbon budget 2021. *Earth Syst. Sci. Data* **14**, 1917–2005 (2022).

9. Gatti, L. V. et al. Amazonia as a carbon source linked to deforestation and climate change. *Nature* **595**, 388–393 (2021).

10. Pan, Y. et al. A large and persistent carbon sink in the world's forests. *Science* **333**, 988–993 (2011).

11. Lovenduski, N. S. & Bonan, G. B. Reducing uncertainty in projections of terrestrial carbon uptake. *Environ. Res. Lett.* **12**, 044020(2017).

12. Chini, L. et al. Land-use harmonization datasets for annual global carbon budgets. *Earth Syst. Sci. Data* **13**, 4175–4189 (2021).

13. Walker, A. P. et al. Predicting long-term carbon sequestration in response to CO2 enrichment: how and why do current ecosystem models differ? *Glob. Biogeochem. Cycles* **29**, 476–495 (2015).

14. De Kauwe, M. G. et al. Where does the carbon go? A model-data intercomparison of vegetation carbon allocation and turnover processes at two temperate forest free-air CO2 enrichment sites. *N. Phytol.* **203**, 883–899 (2014).

15. Pongratz, J., Reick, C. H., Houghton, R. A. & House, J. I. Terminology as a key uncertainty in net land use and land cover change carbon flux estimates. *Earth Syst. Dyn.* **5**, 177–195 (2014).

16. Obermeier, W. A. et al. Modelled land use and land cover change emissions – a spatio-temporal comparison of different approaches. *Earth Syst. Dyn.* **12**, 635–670 (2021).

17. Erb, K.-H. et al. Unexpectedly large impact of forest management and grazing on global vegetation biomass. *Nature* **553**, 73–76 (2018).

18. Pugh, T. A. M. et al. Role of forest regrowth in global carbon sink dynamics. *Proc. Natl Acad. Sci. USA* **116**, 4382–4387 (2019).

19. Kondo, M. et al. Plant regrowth as a driver of recent enhancement of terrestrial CO2 uptake. *Geophys. Res. Lett.* **45**, 4820–4830 (2018).

20. Koven, C. D. et al. Controls on terrestrial carbon feedbacks by productivity versus turnover in the CMIP5 Earth System Models. *Biogeosciences* **12**, 5211–5228 (2015).

21. Nyawira, S. S., Nabel, J. E. M., Brovkin, V. & Pongratz, J. Input-driven versus turnover-driven controls of simulated changes in soil carbon due to land-use change. *Environ. Res. Lett.* **12**, 084015 (2017).

22. Nyawira, S. S., Nabel, J. E. M., Don, A., Brovkin, V. & Pongratz, J. Soil carbon response to land-use change: evaluation of a global vegetation model using meta-data. *Biogeosci. Discuss.* https://doi.org/10.5194/bg-2016-161 (2016).

23. Boysen, L. R. et al. Evaluation of soil carbon dynamics after forest cover change in CMIP6 land models using chronosequences. *Environ. Res. Lett.* **16**, 074030 (2021).

24. Sullivan, M. J. P. et al. Long-term thermal sensitivity of Earth's tropical forests. *Science* **368**, 869–874 (2020).

25. Nemani, R. R. et al. Climate-driven increases in global terrestrial net primary production from 1982 to 1999. *Science* **300**, 1560–1563 (2003).

26. Wang, X. et al. A two-fold increase of carbon cycle sensitivity to tropical temperature variations. *Nature* **506**, 212–215 (2014).

27. Xu, L. et al. Changes in global terrestrial live biomass over the 21st century. *Sci. Adv.* **7**, eabe9829 (2021).

28. Friend, A. D. et al. Carbon residence time dominates uncertainty in terrestrial vegetation responses to future climate and atmospheric CO2. *Proc. Natl Acad. Sci. USA* **111**, 3280–3285 (2014).

29. Signori-Müller, C. et al. Non-structural carbohydrates mediate seasonal water stress across Amazon forests. *Nat. Commun.* **12**, 2310 (2021).

30. Martinez-Vilalta, J. et al. Dynamics of non-structural carbohydrates in terrestrial plants: a global synthesis. *Ecol. Monogr.* **86**, 495–516 (2016).

31. De Kauwe, M. G. et al. Forest water use and water use efficiency at elevated CO2: a model-data intercomparison at two contrasting temperate forest FACE sites. *Glob. Chang. Biol.* **19**, 1759–1779 (2013).

32. Jiang, M., Kelly, J. W. G., Atwell, B. J., Tissue, D. T. & Medlyn, B. E. Drought by CO2 interactions in trees: a test of the water savings mechanism. *N. Phytol.* **230**, 1421–1434 (2021).

33. Forrester, D. I. et al. Self-thinning tree mortality models that account for vertical stand structure, species mixing and climate. *Ecol. Manag.* **487**, 118936 (2021).

34. Rogers, A. et al. A roadmap for improving the representation of photosynthesis in Earth system models. *N. Phytol.* **213**, 22–42 (2017).

35. Reich, P. B., Hungate, B. A. & Luo, Y. Carbon-nitrogen interactions in terrestrial ecosystems in response to rising atmospheric carbon dioxide. *Annu. Rev. Ecol. Evol. Syst.* https://doi.org/10.1146/annurev.ecolsys.37.091305.110039 (2006).

36. Zaehle, S. et al. Evaluation of 11 terrestrial carbon-nitrogen cycle models against observations from two temperate free-air CO₂ enrichment studies. *N. Phytol.* **202**, 803–822 (2014).

37. Meyerholt, J., Sickel, K. & Zaehle, S. Ensemble projections elucidate effects of uncertainty in terrestrial nitrogen limitation on future carbon uptake. *Glob. Chang. Biol.* **26**, 3978–3996 (2020).

38. Yang, H. et al. Comparison of forest above-ground biomass from dynamic global vegetation models with spatially explicit remotely sensed observation-based estimates. *Glob. Chang. Biol.* **26**, 3997–4012 (2020).

39. Gasser, T. & Ciais, P. A theoretical framework for the net land-to-atmosphere CO₂ flux and its implications in the definition of 'emissions from land-use change'. *Earth Syst. Dyn.* **4**, 171–186 (2013).

40. Pongratz, J. et al. Land use effects on climate: current state, recent progress, and emerging topics. *Curr. Clim. Change Rep.* https://doi.org/10.1007/s40641-021-00178-y (2021).

41. Matricardi, E. A. T. et al. Long-term forest degradation surpasses deforestation in the Brazilian Amazon. *Science* **369**, 1378–1382 (2020).

42. Baccini, A. et al. Tropical forests are a net carbon source based on aboveground measurements of gain and loss. *Science* **358**, 230–234 (2017).

43. Doetterl, S. et al. Erosion, deposition and soil carbon: a review of process-level controls, experimental tools and models to address C cycling in dynamic landscapes. *Earth Sci. Rev.* **154**, 102–122 (2016).

44. Van Oost, K. et al. The impact of agricultural soil erosion on the global carbon cycle. *Science* **318**, 626–629 (2007).

45. Sanderman, J., Hengl, T. & Fiske, G. J. Soil carbon debt of 12,000 years of human land use. *Proc. Natl Acad. Sci. USA* **114**, 9575–9580 (2017).

46. Don, A., Schumacher, J. & Freibauer, A. Impact of tropical land-use change on soil organic carbon stocks - a meta-analysis. *Glob. Chang. Biol.* **17**, 1658–1670 (2011).

47. Qiu, C. et al. Large historical carbon emissions from cultivated northern peatlands. *Sci. Adv.* **7**, eabf1332 (2021).

48. Hungate, B. A. et al. Assessing the effect of elevated carbon dioxide on soil carbon: a comparison of four meta-analyses. *Glob. Chang. Biol.* **15**, 2020–2034 (2009).

49. Churkina, G. et al. Synergy of rising nitrogen depositions and atmospheric CO₂on land carbon uptake moderately offsets global warming. *Glob. Biogeochem. Cycles* **23**, (2009).

50. Cheng, W. et al. Synthesis and modeling perspectives of rhizosphere priming. *N. Phytol.* **201**, 31–44 (2014).

51. van Groenigen, K. J., Qi, X., Osenberg, C. W., Luo, Y. & Hungate, B. A. Faster decomposition under increased atmospheric CO₂ limits soil carbon storage. *Science* **344**, 508–509 (2014).

52. Kuzyakov, Y., Horwath, W. R., Dorodnikov, M. & Blagodatskaya, E. Review and synthesis of the effects of elevated atmospheric CO₂
on soil processes: no changes in pools, but increased fluxes and accelerated cycles. *Soil Biol. Biochem.* **128**, 66–78 (2019).

53. Fontaine, S., Bardoux, G., Abbadie, L. & Mariotti, A. Carbon input to soil may decrease soil carbon content. *Ecol. Lett.* **7**, 314–320 (2004).

54. Averill, C., Turner, B. L. & Finzi, A. C. Mycorrhiza-mediated competition between plants and decomposers drives soil carbon storage. *Nature* **505**, 543–545 (2014).

55. Abramoff, R. Z. et al. Improved global-scale predictions of soil carbon stocks with Millennial Version 2. *Soil Biol. Biochem.* **164**, 108466 (2022).

56. Lawrence, D. M. et al. The community land model version 5: description of new features, benchmarking, and impact of forcing uncertainty. *J. Adv. Model. Earth Syst.* **11**, 4245–4287 (2019).

57. Kauppi, P. E. et al. Carbon benefits from Forest Transitions promoting biomass expansions and thickening. *Glob. Chang. Biol.* **26**, 5365–5370 (2020).

58. Zaehle, S. et al. The importance of age-related decline in forest NPP for modeling regional carbon balances. *Ecol. Appl.* **16**, 1555–1574 (2006).

59. Li, F. et al. Historical (1700–2012) global multi-model estimates of the fire emissions from the Fire Modeling Intercomparison Project (FireMIP). *Atmos. Chem. Phys.* **19**, 12545–12567 (2019).

60. Ramo, R. et al. African burned area and fire carbon emissions are strongly impacted by small fires undetected by coarse resolution satellite data. *Proc. Natl. Acad. Sci. USA* **118**, e2011160118 (2021).

61. Brienen, R. J. W. et al. Forest carbon sink neutralized by pervasive growth-lifespan trade-offs. *Nat. Commun.* **11**, 4241 (2020).

62. Bugmann, H. & Bigler, C. Will the CO2 fertilization effect in forests be offset by reduced tree longevity? *Oecologia* **165**, 533–544 (2011).

63. Brienen, R. J. W. et al. Long-term decline of the Amazon carbon sink. *Nature* **519**, 344–348 (2015).

64. McDowell, N. et al. Drivers and mechanisms of tree mortality in moist tropical forests. *N. Phytol.* **219**, 851–869 (2018).

65. Anderegg, W. R. L. & Venturas, M. D. Plant hydraulics play a critical role in Earth system fluxes. *New Phytol.* **226**, 1535–1538 (2020).

66. Eller, C. B. et al. Stomatal optimization based on xylem hydraulics (SOX) improves land surface model simulation of vegetation responses to climate. *New Phytol.* **226**, 1622–1637 (2020).

67. O'Sullivan, M. et al. Climate-driven variability and trends in plant productivity over recent decades based on three global products. *Glob. Biogeochem. Cycles* **34**, e2020GB006613 (2020).

68. Asaadi, A. & Arora, V. K. Implementation of nitrogen cycle in the CLASSIC land model. *Biogeosciences* **18**, 669–706 (2021).

69. Lurton, T. et al. Implementation of the CMIP6 forcing data in the IPSL-CM6A-LR model. *J. Adv. Model. Earth Syst.* **12**, e2019MS001940 (2020).

70. Harris, I., Osborn, T. J., Jones, P. & Lister, D. Version 4 of the CRU TS monthly high-resolution gridded multivariate climate dataset. *Sci. Data* **7**, 109 (2020).

71. Kobayashi, S. et al. The JRA-55 reanalysis: general specifications and basic characteristics. *J. Meteorol. Soc. Japan* **93**, 5–48 (2015).

72. Dlugokencky, E. & Tans, P. Trends in atmospheric carbon dioxide, National Oceanic and Atmospheric Administration, EarthSystem Research Laboratory (NOAA/ESRL). http://www.esrl.noaa.gov/gmd/ccgg/trends/global.html (2020).

73. Hegglin, M., Kinnison, D. & Lamarque, J.-F. CCMI nitrogen surface fluxes in support of CMIP6 - version 2.0. https://doi.org/10.22033/ESGF/INPUT4MIPS.1125 (2016).

74. Lu, C. & Tian, H. Global nitrogen and phosphorus fertilizer use for agriculture production in the past half century: Shifted hot spots and nutrient imbalance. *Earth Syst. Sci. Data Discuss.* **9**, 181–192 (2017).

75. Klein Goldewijk, K., Beusen, A., Doelman, J. & Stehfest, E. Anthropogenic land use estimates for the Holocene - HYDE 3.2. *Earth Syst. Sci. Data* **9**, 927–953 (2017).

76. Klein Goldewijk, K., Dekker, S. C. & van Zanden, J. L. Per-capita estimations of long-term historical land use and the consequences for global change research. *J. Land Use Sci.* **12**, 313–337 (2017).

77. Rosan, T. M. et al. A multi-data assessment of land use and land cover emissions from Brazil during 2000–2019. *Environ. Res. Lett.* **16**, 074004 (2021).

78. Keeling, C. D. et al. Atmospheric carbon dioxide variations at Mauna Loa Observatory, Hawaii. *Tellus* **28**, 538–551 (1976).

## Acknowledgements

M.O.S., P.F. and S.S. have received funding from the European Union's Horizon 2020 research and innovation programme under Grant Agreement No. 821003 (project 4 C).

## Author contributions

M.O.S., P.F. and S.S. designed the study. M.O.S. conducted the analysis. M.O.S. wrote the manuscript with input from all authors. P.A., A.A., V.K.A., V.B., C.D., D.S.G., A.J., E.K., D.K., J.K., S.L., D.L., P.C.M., J.R.M., J.E.M.S.N., B.P., R.S., H.T., N.V., A.P.W., W.Y., X.Y. and S.Z. performed the land model simulations.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary information** The online version contains supplementary material available at https://doi.org/10.1038/s41467-022-32416-8.

**Correspondence** and requests for materials should be addressed to Michael O'Sullivan.

**Peer review information** *Nature Communications* thanks the other anonymous reviewer(s) for their contribution to the peer review of this work.

**Reprints and permission information** is available at http://www.nature.com/reprints

**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 license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license 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 license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2022