Thank you for visiting nature.com. You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

# Soil moisture–atmosphere feedback dominates land carbon uptake variability

## Abstract

Year-to-year changes in carbon uptake by terrestrial ecosystems have an essential role in determining atmospheric carbon dioxide concentrations1. It remains uncertain to what extent temperature and water availability can explain these variations at the global scale2,3,4,5. Here we use factorial climate model simulations6 and show that variability in soil moisture drives 90 per cent of the inter-annual variability in global land carbon uptake, mainly through its impact on photosynthesis. We find that most of this ecosystem response occurs indirectly as soil moisture–atmosphere feedback amplifies temperature and humidity anomalies and enhances the direct effects of soil water stress. The strength of this feedback mechanism explains why coupled climate models indicate that soil moisture has a dominant role4, which is not readily apparent from land surface model simulations and observational analyses2,5. These findings highlight the need to account for feedback between soil and atmospheric dryness when estimating the response of the carbon cycle to climatic change globally5,7, as well as when conducting field-scale investigations of the response of the ecosystem to droughts8,9. Our results show that most of the global variability in modelled land carbon uptake is driven by temperature and vapour pressure deficit effects that are controlled by soil moisture.

## Main

Improving the ability of Earth system models (ESMs) to correctly reproduce the observed variability in land carbon fluxes is essential for building confidence in projections of the long-term response of the carbon cycle to a warming and changing climate10. This research agenda has been evolving rapidly in the past decade thanks to coordinated model comparison experiments11,12, theoretical advances13, model developments14,15, as well as new observations from ground-based networks16,17 and satellite platforms18. Yet, the spread among ESMs remains substantial19,20 and highlights the need to better constrain the sensitivity of increasingly complex biogeochemical models to changes in atmospheric and hydrological drivers such as radiation21, temperature7, soil water availability3 and vapour pressure deficit (VPD; a measure of atmospheric dryness that depends on air temperature and humidity). In particular, it remains unclear whether temperature or soil moisture is the dominant driver of the inter-annual variability (IAV) in land carbon uptake at the global scale2,3,4,5. Here, we investigate the extent to which temperature, VPD and soil moisture effects co-vary as a result of soil moisture–atmosphere feedback, and reconcile conflicting assessments of the sensitivity of global carbon fluxes to these variables.

Soil moisture drought is one of the key prerequisites for the development of extremely high temperatures22,23,24, whereas atmospheric dynamics control the onset of such extremes25. During droughts, low soil moisture content limits evapotranspiration, which is the most efficient surface cooling flux26. This modification of the surface energy balance increases the air temperature, lowers the relative humidity and thus raises VPD. The importance of such soil moisture–atmosphere feedback, hereafter referred to as land–atmosphere coupling (LAC), is confirmed by both models and observations27,28,29. In current carbon cycle models, the impacts of soil moisture, temperature and VPD on ecosystem productivity and respiration are usually parameterized using stress functions. Typically, simulated photosynthesis rates are limited by low soil moisture content and extreme temperatures via a scaling of Vcmax (the maximum rate of Rubisco carboxylase activity)30 or through a downregulation of stomatal conductance (gs) in response to VPD, relative humidity or a soil water stress function31,32. Ecosystem respiration and fire occurrences are also controlled by soil moisture content, temperature or atmospheric dryness33,34. Because of this situation, the overall influence of soil moisture can potentially occur as (1) a direct impact on photosynthesis and respiration processes through the soil water stress regulation or (2) as an indirect response to extreme temperature and VPD anomalies resulting from LAC.

Here, we investigate the magnitude of these two different causal pathways (that is, direct and indirect) using coupled climate model simulations from the Global Land-Atmosphere Coupling Experiment, Coupled Model Intercomparison Project 5 (GLACE-CMIP5)6 (Methods). To identify the overall influence of soil moisture variability on carbon fluxes and atmospheric conditions, we use an experiment (experiment A) in which the non-seasonal variability in soil moisture is artificially removed. This is achieved by forcing the soil moisture in experiment A to follow the mean seasonal soil moisture cycle calculated from a reference control simulation (CTL) (Extended Data Figs. 1, 2). Experiment A thus simulates the temperature, VPD and carbon fluxes that would occur under climatologically normal soil moisture conditions. We note that sea surface temperatures (SSTs) are identical in experiments A and CTL. This ensures that the main differences between experiment A and CTL are due to the different soil moisture conditions and are not caused by differences in SST patterns (Methods). Using this framework, previous studies have shown that suppressing the non-seasonal soil moisture variability in experiment A strongly reduces the magnitude of temperature and VPD extremes compared to the control simulation6,27,35 (Extended Data Fig. 3). Here, by comparing the carbon flux anomalies of experiment A with those of the control simulation, we are able to estimate the overall magnitude of soil moisture effects (that is, direct and indirect effects) on the IAV of net biome production (NBP; which represents the net land carbon uptake). Because we focus on the IAV, all presented figures are based on anomalies (de-seasoned and de-trended data) from the period 1960–2005, unless otherwise noted.

Our results show that suppressing non-seasonal variability in soil moisture leads to a 91% (standard deviation of ±2.3%) decrease in the variance of global mean NBP, consistently across all of the four participating climate models (Fig. 1a, Supplementary Table 1). In other words, without soil moisture variability, the IAV of net land carbon uptake is almost eliminated. This primarily occurs because of a reduction in the IAV of gross primary production (GPP) (Fig. 1b, c, Supplementary Table 1) and to a lesser extent because of a reduction in the IAV of ecosystem respiration and disturbance fluxes (the sum of autotrophic and heterotrophic respiration, fires and any other modelled disturbance). As explained above, both direct soil moisture effects and indirect temperature and VPD effects related to LAC can be responsible for the widespread reduction of NBP variability occurring in experiment A (Fig. 2a).

Using a sensitivity analysis (equations (1), (2), Supplementary Figs. 13) of the local model response to anomalies in soil moisture, temperature, VPD and shortwave solar radiation in CTL versus experiment A, we isolate the contributions of direct soil moisture effects (Fig. 2b) versus indirect effects (Fig. 2c) to the overall reduction in NBP variability (Fig. 2a). Regionally, direct soil moisture effects are found in both temperate and tropical biomes, whereas indirect effects occurring through the feedback on temperature and VPD are mostly concentrated in semi-arid and tropical regions. Our sensitivity analysis also shows that most of the reduction in NBP variability found in experiment A occurs because of a reduction in the variance of the climatological drivers, rather than because of a change in the sensitivity of NBP to these drivers (Extended Data Fig. 4). These findings demonstrate that soil moisture can affect carbon uptake variability in two different and equally important ways. First, soil moisture variability has direct effects on NBP, mostly because plant photosynthesis is reduced when soils become dry below a certain threshold (Fig. 2b); second, it enhances temperature and VPD anomalies through LAC, thus leading to indirect effects on NBP (Fig. 2c, Extended Data Fig. 5). Importantly, some regions can be more sensitive to indirect effects (that is, soil moisture feedback mechanisms on temperature and VPD) than to direct soil moisture effects (Extended Data Fig. 6). We note that because disentangling the individual contributions of temperature and VPD to NBP variability is not straightforward, only their joint contribution is reported here (see Methods for a discussion).

When aggregating these results to the global scale (Fig. 3a), we find that indirect effects alone are on average (across models) responsible for most (60%) of the global NBP IAV, whereas direct soil moisture effects account for only 20%. Suppressing direct and indirect effects together leads to a net decrease in NBP variance of about 90% (consistent with Fig. 1) as a result of the positive covariance between the direct and indirect effects (Supplementary Tables 2, 3). Finally, the temperature (T) and VPD effects that are independent of soil moisture conditions and still persist in experiment A ($${{\rm{NBP}}}_{{\rm{nonLAC}}}^{{\rm{T}}\& {\rm{VPD}}}$$) account for only 9% of the overall global NBP variability, whereas radiation effects account for the remaining 11%. As a result of spatial aggregation (Fig. 3b), indirect effects also tend to increase in relative importance as they are spatially more coherent (probably owing to atmospheric mixing) and do not average out as fast as the direct effects2. In summary, the largest fraction of the global mean NBP IAV is driven by anomalies in temperature and VPD that represent an indirect response to soil moisture variability (given that they do not occur in its absence, as demonstrated by the experiment). This finding reconciles opposing perspectives on the roles of temperature versus water availability2,3,4,5, because the apparent importance of either driver actually depends on whether the indirect (feedback) effects are attributed to temperature or soil moisture (see Extended Data Fig. 7, Supplementary Fig. 5). Although it is not possible to replicate the factorial experiment with observations (this would require manipulating soil moisture everywhere on the planet), we assess the degree to which the reference simulations reflect real observations. Evaluating the control simulations against observational estimates, we find that the modelled sensitivity of global NBP IAV to the different meteorological drivers (Fig. 3) agrees well with two independent observational products (Extended Data Fig. 8). Taking into account the uncertainty of these observations, the spatial patterns of NBP IAV simulated by the models are also in reasonable agreement with real-world variability (Supplementary Fig. 6; see discussion in Methods).

More generally, our results show that the areas where NBP IAV is the largest overall (Fig. 4a) often correspond to those where the reduction of temperature and VPD variability due to prescribing soil moisture is the strongest (Fig. 4b, c). In other words, NBP variability tends to be larger where LAC is stronger (Fig. 4d). These known hotspots of LAC36 match well with earlier studies that suggested that semi-arid regions dominate global NBP IAV37,38, even though our analysis refines these previous findings (Extended Data Fig. 9) by also including regions that are usually classified as temperate or humid but that are affected by LAC for only a few dry months during the year (for example, eastern Europe22, Amazon basin39).

These results also improve our understanding of the sensitivity of land carbon uptake IAV to tropical mean temperature40,41, which has been used to constrain coupled climate model projections7,42. Here, we find that the IAV of mean tropical land temperature is barely changed in the experiment with prescribed soil moisture (Extended Data Fig. 10). This is because suppressing soil moisture anomalies reduces temperature extremes only in a couple of hotspot regions (Fig. 4b, Extended Data Fig. 3) with little impact on the overall tropical mean. Thus, although the IAV in global land carbon uptake has been empirically found to be sensitive to tropical mean temperature in numerous studies5,41, our results suggest that this sensitivity does not represent a strong mechanistic link, and thus might not necessarily represent the most adequate model constraint. In fact, the El Niño Southern Oscillation and SST in general may be the confounding driver of both tropical mean temperature and the precipitation patterns that cause the soil moisture anomalies leading to NBP variability.

In conclusion, we show that the IAV in land carbon uptake simulated by ESMs is primarily driven by anomalies in temperature and VPD, which are themselves controlled by soil moisture variability. These indirect soil moisture effects occur through LAC and account for 60% (±18%) of the simulated global land carbon uptake IAV. They explain why the simulated global NBP variability (1) mainly arises from tropical and semi-arid regions37,38, which are known hotspots of LAC6,36,43, (2) is predominantly a temperature and VPD response (at the global scale) according to land surface models and empirical sensitivity analyses2,5 and (3) is also largely dependent on soil moisture variability according to coupled climate simulations4. Our results reveal that soil moisture–atmosphere feedback mechanisms represent a dominant source of variability in global carbon uptake and thus reconcile previous conflicting assessments2,3,4,5. To some extent, we note that these findings might be symptomatic of how land surface models were developed in the first place. Parameterizing a strong sensitivity of carbon uptake to observed VPD or temperature can constitute a simpler way for a land-surface model to achieve good skill, especially when soil water stress and soil moisture dynamics are only represented approximately. As a result, even though models strongly agree that direct and indirect soil moisture effects together dominate land carbon uptake variability, the actual partitioning between direct and indirect effects may be more dependent on modelling approaches. More generally, our results illustrate the importance of differentiating estimates of ecosystem sensitivity to natural droughts, as opposed to artificial droughts (for example, rainfall exclusion experiments), given that only the former incorporates LAC and its impact on temperature and humidity. Because soil and atmospheric dryness do not equally respond to climate change27,44, the direct and indirect soil moisture effects identified here might affect future NBP in different ways. Since current climate models have a large spread in their representation of vegetation response to dryness45 and of LAC strength46, this could introduce uncertainties in the feedback that are difficult to diagnose from offline land surface model evaluation efforts47, with potentially large impacts on carbon fluxes, as demonstrated here. We also note that long-term changes in vegetation structure and composition might alter the ecosystem’s future response4 to and control9,48,49 of soil moisture–atmosphere feedback. Thus, more physical and holistic representations of the response of vegetation to soil and atmospheric dryness might have a strong potential to reduce key uncertainties in current projections of future terrestrial carbon fluxes.

## Methods

### Model experiment

The presented results are based on the Global Land–Atmosphere Coupling Experiment – Coupled Model Intercomparison Project phase 5 (GLACE-CMIP5) numerical experiment6. This model experiment was originally designed to investigate soil moisture–climate feedback mechanisms under historical and future scenarios, and notably their impact on extreme heat events6. Its experimental design is inspired from the original GLACE experiment43, which focused on the role of soil moisture in seasonal weather predictability. Six ESMs were used for global climate simulations: the Community Climate System Model 4 (CCSM4), the European community ESM (EC-Earth), the European Centre/Hamburg Model 6 (ECHAM6), the Geophysical Fluid Dynamics Laboratory model (GFDL), the Institut Pierre‐Simon Laplace model (IPSL) and the Australian Community Climate and Earth System Simulator (ACCESS). Model outputs for carbon fluxes are available only for four models (CCSM4, ECHAM6, GFDL and IPSL), and the availability of certain variables is limited in some cases (Supplementary Table 4), which explains why some analyses cannot be conducted with all models (for example, Fig. 1c).

The control (CTL) and soil moisture (experiment A) experiments consist of coupled atmosphere–land simulations (Extended Data Fig. 2) using prescribed SSTs, sea ice, land use and atmospheric CO2 concentrations from each of the model’s fully coupled reference CMIP5 runs (except for CCSM4, where the reference CMIP5 run itself is used as the control simulation). Unlike so-called ‘offline’ simulations, in which a land surface model is driven by a fixed meteorological forcing, a coupled simulation resolves water and energy exchanges between the land and the atmosphere, allowing land processes to feed back to the atmosphere and influence it locally. The model simulations cover the historical period since 1950 and the 21st century (RCP8.5 scenario). Further details documenting the control experiment, including the description of the atmospheric and land model components, can be found in Seneviratne et al.6. The only forced difference between the CTL and experiment A simulations is the soil moisture variability. In experiment A, soil moisture is prescribed to a reference climatology (seasonal cycle) calculated from the control run over the period 1971–2000 (Extended Data Fig. 1). Thus, the main difference (on a climatological timescale) between the two simulations is related to the change in soil moisture. It is worth noting that at finer, meteorological, timescales (for example, daily time series), the internal variability inherent to general circulation models will also lead to differences between the two simulations.

Prescribing soil moisture implies that the water balance is not necessarily conserved. An investigation of this imbalance with the Community Earth System Model (CESM)50 showed a positive net imbalance (that is, the sum of all water additions and subtractions) of the order of +8% globally (relative to the annual mean precipitation), associated with an overall increase in land evapotranspiration. We note that in some specific regions, less water may be added than is removed (negative imbalance) and that temperature extremes are found to be reduced in both cases (positive or negative imbalance) as a result of the suppressed LAC. Although there is no apparent impact on global mean precipitation50, there are some changes in the distribution of precipitation (for example, an increase in extreme events over the tropics51). We do not expect changes in precipitation between CTL and experiment A to have any impact on carbon fluxes (because soil moisture is prescribed).

To enable a consistent comparison, we re-grid all model outputs to a common resolution of 2° using conservative re-gridding and compute monthly averages. The entire analysis presented here is focused on the IAV over the period 1960–2005. We note that VPD is first calculated from daily averages of temperature and relative humidity, and is only then averaged to monthly means. The IAV corresponds to the signal remaining after removing the seasonal cycle, as well as any long-term linear trend on a monthly basis (the long-term trend of each month is subtracted). For the ECHAM6 model, two grid cells located in the Tibetan plateau are discarded from the whole analysis, because spurious spikes are present in heterotrophic respiration for experiment A. We also discard Greenland and Antarctica to maintain a comparable spatial coverage for all models. Although this work focuses on anomalies (that is, deviations from the seasonal cycle), we also illustrate the seasonal cycles of NBP, GPP and respiration and disturbance simulated in CTL and experiment A in Supplementary Fig. 7. For completeness, we also provide time series of global mean soil moisture, temperature, VPD and radiation IAV (similar to Fig. 1) in Supplementary Fig. 8.

### Comparison of the control simulations with observational estimates

We evaluate the simulated IAV of NBP, soil moisture, temperature and VPD against available observations in Supplementary Figs. 6, 911. For NBP IAV (Supplementary Fig. 6), we note that although observational estimates of NBP variability exist, they do not agree well with each other, reflecting our limited knowledge of net carbon fluxes globally52,53 (Supplementary Fig. 6g, ‘obs vs obs’). To focus on time periods in which these observational datasets are more reliable globally, we use the period 1980–2010 for the FLUXCOM RS+METEO dataset and the period 2000–2018 for the CAMS atmospheric CO2 inversion. We show that models correlate with these observational estimates as much as the observations themselves correlate with each other (Supplementary Fig. 6g, ‘models vs obs’). We also find that there is little consensus on the overall (de-trended) NBP IAV amplitude. The global mean NBP standard deviation of the different models ranges from 0.86 petagrams of carbon per year (Pg C yr−1) for CCSM4 to 2.76 Pg C yr−1 for GFDL. When comparing with observational products (Supplementary Fig. 6h), we find that—excluding FLUXCOM RS+METEO, which is known to underestimate the global NBP IAV52—the CAMS atmospheric CO2 inversion53 suggests a value of 0.68 Pg C yr−1, whereas dynamic vegetation models used for the Global Carbon Project1 suggest a range of 0.53 to 1.50 Pg C yr−1. Thus, some models (GFDL in particular) seem to overestimate the overall NBP variability. However, regardless of how close they are to observations or other estimates, all models are unanimous that the global NBP variance is reduced by about 90% when prescribing soil moisture and that indirect effects dominate this response (Figs. 1, 3).

We evaluate spatial patterns of IAV for soil moisture, temperature and VPD against available observational datasets in Supplementary Figs. 911. The simulated soil moisture IAV patterns agree reasonably well with the total soil moisture from the ERA5-Land reanalysis54 and with satellite observations of shallow soil moisture (5–10 cm depth) from the ESA CCI Combined product v4.555 (Supplementary Fig. 9). Regarding temperature and VPD IAV, we find that models and observational sources56,57 are in reasonable agreement (Supplementary Figs. 10, 11). Finally, we also evaluate spatial patterns of global long-term mean GPP, which is arguably better constrained by observations than long-term mean NBP. We find that the models agree very well with the observational data52,58 in terms of spatial patterns (Supplementary Fig. 12). For global mean GPP, two models produce a relatively high global mean GPP (of about 150 Pg C yr−1). However, such values are not entirely unrealistic according to other satellite-based estimates (for example, Joiner et al.59 report 140 Pg C yr−1).

### Sensitivity analysis

In Figs. 2, 3 we reproduce the approach by Jung et al.2, which consists of a local month-wise linear regression of the NBP model output against the main meteorological drivers (which are also deseasonalized and detrended):

$${\bf{N}}{\bf{B}}{{\bf{P}}}_{s,m}^{* }={\beta }_{s,m}^{{\rm{SM}}}{\bf{S}}{{\bf{M}}}_{s,m}+{\beta }_{s,m}^{{\rm{T}}}{{\bf{T}}}_{s,m}+{\beta }_{s,m}^{{\rm{VPD}}}{\bf{V}}{\bf{P}}{{\bf{D}}}_{s,m}+{\beta }_{s,m}^{{\rm{R}}}{{\bf{R}}}_{s,m},$$
(1)

where s is the spatial index (grid point), m is the month index (1 to 12) and β are regression coefficients. NBP, SM, T, VPD and R are N × 1 vectors, where N is the number of years; NBP denotes the net biome production anomaly, SM represents the total soil moisture anomaly, T denotes the 2-m air temperature anomaly, VPD represents the vapour pressure deficit anomaly and R is the surface downward solar radiation anomaly. In the main text, the four components of equation (1) are referred to using the more compact notation:

$${{\rm{NBP}}}^{* }={{\rm{NBP}}}^{{\rm{SM}}}+{{\rm{NBP}}}^{{\rm{T}}}+{{\rm{NBP}}}^{{\rm{VPD}}}+{{\rm{NBP}}}^{{\rm{R}}},$$
(2)

where NBPSM, NBPT, NBPVPD and NBPR correspond to the soil-moisture-driven, temperature-driven, VPD-driven and radiation-driven NBP, respectively, and NBP* is the overall result of the regression. This regression is applied to CTL and experiment A simulations separately (each regression is referred to using the appropriate notation $${{\rm{NBP}}}_{{\rm{CTL}}}^{\ast }$$ or $${{\rm{NBP}}}_{{\rm{ExpA}}}^{\ast }$$). In Fig. 2b, c, the difference in annual NBP variability is calculated by subtracting the standard deviation of the components of equation (2) from both experiments (for example, $$\Delta \sigma ({{\rm{NBP}}}^{{\rm{SM}}})\,=$$ $$\sigma ({{\rm{NBP}}}_{{\rm{ExpA}}}^{{\rm{SM}}})-\sigma ({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{SM}}})$$).

Because this statistical approach does not incorporate other potential sources of NBP variability as explanatory variables (ecosystem memory in particular, but also fires) and can only capture linear relationships within a given month, it should not be expected to capture the full complexity of ESM outputs. Our evaluation shows that this approach is able to reproduce a correct NBP IAV at the global (Supplementary Figs. 1, 2) and local (Supplementary Fig. 3) scales, although the overall NBP variability is generally underestimated because of the reasons mentioned above. We also apply this statistical approach to two fully independent observational estimates of NBP fluxes. We use the FLUXCOM RS+METEO dataset (GSWP3 version) over the period 1981–201052, which is a machine learning-based upscaling of flux tower measurements, and the CAMS v18r3 dataset53, which is an atmospheric CO2 inversion, over the period 2000–2018. We find that the overall partitioning of global NBP IAV between the different drivers is similar to what models are suggesting (Extended Data Fig. 8). The ability of the regression to reproduce these observational estimates is shown in Supplementary Fig. 13. For FLUXCOM, the sensitivity analysis is able to capture the variability almost perfectly. This is only possible because we use the same predictors as the ones used by the machine learning algorithms (that is, the GSWP3 meteorological forcing60). As a result, there is perfect internal consistency between FLUXCOM NEE and its predictors. For the CAMS inversion, however, such internal consistency does not exist. Using ERA5-Land54 soil moisture, temperature, VPD and radiation as predictors, we find that the sensitivity analysis agrees relatively well with the models, even though it underestimates the magnitude of CAMS NBP anomalies at the global scale. Locally, this regression performs moderately well (Supplementary Fig. 13f), which is nonetheless a reasonable result when considering the very high uncertainty of regional NBP anomalies when they are derived from CO2 inversions at the sub-continental scale53.

Of particular interest to this work is the difference in NBP variance between CTL and experiment A (Fig. 2a). We find that this difference can be reproduced very well by the sensitivity analysis for three out of the four models (Supplementary Fig. 4). Differences are underestimated for the CCSM4 model, but this seems to occur rather uniformly and most spatial patterns are preserved (thus the ratio of NBP variance between CTL and experiment A estimated from the regression is close to the actual one; see Supplementary Table 3). Closer inspection of the regression residuals suggests that ecosystem memory and lag effects (which cannot be captured by equation (1)) might be particularly important for this model. It is interesting to note that for some models (for example, GFDL), the NBP variance can also increase locally when seasonal soil moisture is prescribed (Supplementary Fig. 4). This occurs only in a few arid regions that have almost no NBP variability in the control simulation and where soil moisture is extremely low except during occasional wet years. Prescribing a mean seasonal soil moisture in those regions causes small amounts of soil water to be available every year (instead of every few years), which increases the overall NBP variability.

Finally, several alternative formulations to equation (1) were tested. The chosen formulation (equation (1)) is the one that best reproduces the model NBP outputs. Potential alternative formulations may consist in (i) using only soil moisture, temperature and radiation, as in Jung et al.2; (ii) including an interaction term between temperature and soil moisture instead of VPD; (iii) replacing VPD by relative humidity. Using any of these three alternative formulations does not affect the main finding of the study, that is, that most of the global NBP variability is driven by indirect soil moisture effects (see Supplementary Figs. 5, 14, 15).

### Joint analysis of temperature and VPD effects

In Figs. 2, 3 the contributions of temperature and VPD are represented as a sum (NBPT&VPD = NBPT + NBPVPD). This is because temperature and VPD are correlated to some extent (VPD is calculated from the temperature and the relative humidity), so that the ability of the sensitivity analysis to attribute NBP anomalies to either one of these two variables (that is, temperature versus VPD) might be limited in some cases. We recognize this potential limitation by analysing the joint contribution of these two variables. For completeness, individual contributions are also illustrated in Extended Data Figs. 4, 5. With the caveats mentioned above, Extended Data Fig. 4 shows that VPD has a much larger role than temperature in the reduction of NBP variability occurring between CTL and experiment A. However, this does not mean that temperature is less sensitive than VPD to prescribing soil moisture. Rather, Extended Data Fig. 5 shows that the sensitivity analysis attributes more NBP variability to VPD to begin with, but that both the VPD-driven and temperature-driven NBP variability are reduced in experiment A.

### Variance contributions at different levels of aggregation

In Fig. 3, Extended Data Fig. 7 and Supplementary Figs. 5, 1416, the contribution of different drivers to NBPCTL variance is computed at different levels of spatial aggregation. The following different levels of aggregation are used: 2°, 3°, 4°, 5°, 6°, 7.5°, 9°, 10°, 12°, 15°, 18°, 20°, 22.5°, 30°, 36°, 45°, 60°, 90°, 180° and 360° (that is, global). Contributions are calculated as follows. Similarly to Jung et al.2, the different NBP time series (NBPSM, NBPT&VPD and NBPR) are first aggregated to the given spatial resolution. After aggregation, the variance of the time series (that is, $${\sigma }^{2}({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{SM}}\,})$$ and so on) are computed at each grid point. Then, the variance of the T&VPD contribution $${\sigma }^{2}({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{T}}\& {\rm{VPD}}\,})$$ is decomposed at each grid point into an LAC-dependent and non LAC-dependent contribution, as explained in Supplementary Information section 2. After that, and similar to Jung et al.2, the global spatial average of the variances is calculated for each of the four contributions (for example, $$\overline{{\sigma }^{2}({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{SM}}\,})}$$). The relative contribution of a component at a given level of spatial aggregation (as shown in Fig. 3b) is then calculated by normalizing that global spatial average against the sum of all components:

$$\begin{array}{c}{\rm{Contribution}}\,({{\rm{NBP}}}^{{\rm{SM}}})\,=\\ \frac{\overline{{\sigma }^{2}({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{SM}}\,})}}{\overline{{\sigma }^{2}({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{SM}}\,})}+\overline{{\sigma }^{2}({{\rm{NBP}}}_{{\rm{LAC}}}^{{\rm{T}}\& {\rm{VPD}}\,})}+\overline{{\sigma }^{2}({{\rm{NBP}}}_{{\rm{nonLAC}}}^{{\rm{T}}\& {\rm{VPD}}\,})}+\overline{{\sigma }^{2}({{\rm{NBP}}}_{{\rm{CTL}}}^{{\rm{R}}\,})}}.\end{array}$$
(3)

Identically to Jung et al.2, the spread in the contributions estimated by the four different models shown in Extended Data Fig. 7 is reported in two different ways. The outer uncertainty bounds represent the standard deviation of the contribution estimated by the four models. The inner uncertainty bounds represent the standard deviation between the four estimates, but after removing the mean contribution of each model across all levels of aggregation. Thus, the inner uncertainty bounds show the uncertainty in the tendency of the contribution (its change from regional to global scale).

## Data availability

GLACE-CMIP5 model outputs can be obtained from S.I.S. (sonia.seneviratne@ethz.ch). FluxCom data are available at http://www.fluxcom.org/CF-Download/. CAMS data are available from the Atmosphere Data Store at https://atmosphere.copernicus.eu/data. ERA5 and ERA5Land data are available from the Climate Data Store at https://cds.climate.copernicus.eu. VPM-GPP is available at https://doi.org/10.6084/m9.figshare.c.3789814. ESA CCI Soil Moisture is available at https://www.esa-soilmoisture-cci.org. CRU TS data are available at https://crudata.uea.ac.uk/cru/data/hrg/. GSWP3 data are available at https://doi.org/10.20783/DIAS.501. The corresponding author can also be contacted at vincent.humphrey@bluewin.ch. Source data are provided with this paper.

## Code availability

Code and documentation for CCSM4 is publicly available at https://www.cesm.ucar.edu/models/ccsm4.0/. Code and documentation for ECHAM6 (MPI-ESM) is available for scientific users at https://mpimet.mpg.de/en/science/modeling-with-icon/code-availability. Code and documentation for the GFDL model is publicly available at https://www.gfdl.noaa.gov/modeling-systems-group-public-releases/. Code and documentation for the IPSL model is publicly available at https://cmc.ipsl.fr/ipsl-climate-models/ipsl-cm5/. Model outputs were processed using the software Matlab 2019a.

## References

1. 1.

Friedlingstein, P. et al. Global carbon budget 2019. Earth Syst. Sci. Data 11, 1783–1838 (2019).

2. 2.

Jung, M. et al. Compensatory water effects link yearly global land CO2 sink changes to temperature. Nature 541, 516–520 (2017).

3. 3.

Humphrey, V. et al. Sensitivity of atmospheric CO2 growth rate to observed changes in terrestrial water storage. Nature 560, 628–631 (2018).

4. 4.

Green, J. K. et al. Large influence of soil moisture on long-term terrestrial carbon uptake. Nature 565, 476–479 (2019).

5. 5.

Piao, S. et al. Interannual variation of terrestrial carbon cycle: issues and perspectives. Glob. Change Biol. 26, 300–318 (2019).

6. 6.

Seneviratne, S. I. et al. Impact of soil moisture-climate feedbacks on CMIP5 projections: first results from the GLACE-CMIP5 experiment. Geophys. Res. Lett. 40, 5212–5217 (2013).

7. 7.

Cox, P. M. et al. Sensitivity of tropical carbon to climate change constrained by carbon dioxide variability. Nature 494, 341–344 (2013).

8. 8.

Novick, K. A. et al. The increasing importance of atmospheric demand for ecosystem water and carbon fluxes. Nat. Clim. Chang. 6, 1023–1027 (2016).

9. 9.

Anderegg, W. R. L., Trugman, A. T., Bowling, D. R., Salvucci, G. & Tuttle, S. E. Plant functional traits and climate influence drought intensification and land–atmosphere feedbacks. Proc. Natl Acad. Sci. USA 116, 14071–14076 (2019).

10. 10.

Ciais, P. et al. in Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change 465–570 (Cambridge Univ. Press, 2014).

11. 11.

Sitch, S. et al. Recent trends and drivers of regional sources and sinks of carbon dioxide. Biogeosciences 12, 653–679 (2015).

12. 12.

Jones, C. D. et al. C4MIP – The Coupled Climate–Carbon Cycle Model Intercomparison Project: experimental protocol for CMIP6. Geosci. Model Dev. 9, 2853–2880 (2016).

13. 13.

Tramontana, G. et al. Predicting carbon dioxide and energy fluxes across global FLUXNET sites with regression algorithms. Biogeosciences 13, 4291–4313 (2016).

14. 14.

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. 3, M03001 (2011).

15. 15.

Kennedy, D. et al. Implementing plant hydraulics in the Community Land Model, Version 5. J. Adv. Model. Earth Syst. 11, 485–513 (2019).

16. 16.

Baldocchi, D. et al. FLUXNET: a new tool to study the temporal and spatial variability of ecosystem-scale carbon dioxide, water vapor, and energy flux densities. Bull. Am. Meteorol. Soc. 82, 2415–2434 (2001).

17. 17.

Wunch, D. et al. The total carbon column observing network. Philos. Trans. R. Soc. A 369, 2087–2112 (2011).

18. 18.

Schimel, D. et al. Observing terrestrial ecosystems and the carbon cycle from space. Glob. Change Biol. 21, 1762–1776 (2015).

19. 19.

Friedlingstein, P. et al. Uncertainties in CMIP5 climate projections due to carbon cycle feedbacks. J. Clim. 27, 511–526 (2014).

20. 20.

Arora, V. K. et al. Carbon–concentration and carbon–climate feedbacks in CMIP6 models and their comparison to CMIP5 models. Biogeosciences 17, 4173–4222 (2020).

21. 21.

Mercado, L. M. et al. Impact of changes in diffuse radiation on the global land carbon sink. Nature 458, 1014–1017 (2009).

22. 22.

Seneviratne, S. I., Luthi, D., Litschi, M. & Schar, C. Land-atmosphere coupling and climate change in Europe. Nature 443, 205–209 (2006).

23. 23.

Hirschi, M., Mueller, B., Dorigo, W. & Seneviratne, S. I. Using remotely sensed soil moisture for land-atmosphere coupling diagnostics: the role of surface vs. root-zone soil moisture variability. Remote Sens. Environ. 154, 246–252 (2014).

24. 24.

Yin, D., Roderick, M. L., Leech, G., Sun, F. & Huang, Y. The contribution of reduction in evaporative cooling to higher surface air temperatures during drought. Geophys. Res. Lett. 41, 7891–7897 (2014).

25. 25.

Miralles, D. G., Teuling, A. J., van Heerwaarden, C. C. & Vila-Guerau de Arellano, J. Mega-heatwave temperatures due to combined soil desiccation and atmospheric heat accumulation. Nat. Geosci. 7, 345–349 (2014).

26. 26.

Bateni, S. M. & Entekhabi, D. Relative efficiency of land surface energy balance components. Wat. Resour. Res. 48, (2012).

27. 27.

Zhou, S. et al. Land–atmosphere feedbacks exacerbate concurrent soil drought and atmospheric aridity. Proc. Natl Acad. Sci. USA 116, 18848–18853 (2019).

28. 28.

Hirschi, M. et al. Observational evidence for soil-moisture impact on hot extremes in southeastern Europe. Nat. Geosci. 4, 17–21 (2011).

29. 29.

Wehrli, K., Guillod, B. P., Hauser, M., Leclair, M. & Seneviratne, S. I. Identifying key driving processes of major recent heat waves. J. Geophys. Res. D 124, 11746–11765 (2019).

30. 30.

Sellers, P. J. et al. A revised land surface parameterization (SiB2) for atmospheric GCMS. Part I: model formulation. J. Clim. 9, 676–705 (1996).

31. 31.

Leuning, R. A critical appraisal of a combined stomatal-photosynthesis model for C3 plants. Plant Cell Environ. 18, 339–355 (1995).

32. 32.

Medlyn, B. E. et al. Reconciling the optimal and empirical approaches to modelling stomatal conductance. Glob. Change Biol. 17, 2134–2144 (2011).

33. 33.

Yan, Z. et al. A moisture function of soil heterotrophic respiration that incorporates microscale processes. Nat. Commun. 9, 2562 (2018).

34. 34.

Metcalfe, D. B. et al. Shifts in plant respiration and carbon use efficiency at a large-scale drought experiment in the eastern Amazon. New Phytol. 187, 608–621 (2010).

35. 35.

Berg, A. et al. Interannual coupling between summertime surface temperature and precipitation over land: processes and implications for climate change. J. Clim. 28, 1308–1328 (2015).

36. 36.

Dirmeyer, P. A. The terrestrial segment of soil moisture–climate coupling. Geophys. Res. Lett. 38, L16702 (2011).

37. 37.

Ahlstrom, A. et al. The dominant role of semi-arid ecosystems in the trend and variability of the land CO2 sink. Science 348, 895–899 (2015).

38. 38.

Poulter, B. et al. Contribution of semi-arid ecosystems to interannual variability of the global carbon cycle. Nature 509, 600–603 (2014).

39. 39.

Levine, P. A. et al. Soil moisture variability intensifies and prolongs eastern Amazon temperature and carbon cycle response to El Niño–Southern Oscillation. J. Clim. 32, 1273–1292 (2019).

40. 40.

Keeling, C. D., Whorf, T. P., Wahlen, M. & Vanderplicht, J. Interannual extremes in the rate of rise of atmospheric carbon-dioxide since 1980. Nature 375, 666–670 (1995).

41. 41.

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

42. 42.

Huntzinger, D. N. et al. Uncertainty in the response of terrestrial carbon sink to environmental drivers undermines carbon-climate feedback predictions. Sci. Rep. 7, 4765 (2017).

43. 43.

Koster, R. D. et al. GLACE: the Global Land–Atmosphere Coupling Experiment. Part I: overview. J. Hydrometeorol. 7, 590–610 (2006).

44. 44.

Cook, B. I. et al. Twenty‐first century drought projections in the CMIP6 forcing scenarios. Earths Futur. 8, (2020).

45. 45.

Trugman, A. T., Medvigy, D., Mankin, J. S. & Anderegg, W. R. L. Soil moisture stress as a major driver of carbon cycle uncertainty. Geophys. Res. Lett. 45, 6495–6503 (2018).

46. 46.

Berg, A. & Sheffield, J. Soil moisture–evapotranspiration coupling in CMIP5 models: relationship with simulated climate and projections. J. Clim. 31, 4865–4878 (2018).

47. 47.

Collier, N. et al. The International Land Model Benchmarking (ILAMB) system: design, theory, and implementation. J. Adv. Model. Earth Syst. 10, 2731–2754 (2018).

48. 48.

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. Change Biol. 19, 1759–1779 (2013).

49. 49.

Swann, A. L. S., Hoffman, F. M., Koven, C. D. & Randerson, J. T. Plant responses to increasing CO2 reduce estimates of climate impacts on drought severity. Proc. Natl Acad. Sci. USA 113, 10019–10024 (2016).

50. 50.

Hauser, M., Orth, R. & Seneviratne, S. I. Investigating soil moisture–climate interactions with prescribed soil moisture experiments: an assessment with the Community Earth System Model (version 1.2). Geosci. Model Dev. 10, 1665–1677 (2017).

51. 51.

Lorenz, R. et al. Influence of land–atmosphere feedbacks on temperature and precipitation extremes in the GLACE-CMIP5 ensemble. J. Geophys. Res. D 121, 607–623 (2016).

52. 52.

Jung, M. et al. Scaling carbon fluxes from eddy covariance sites to globe: synthesis and evaluation of the FLUXCOM approach. Biogeosciences 17, 1343–1365 (2020).

53. 53.

Chevallier, F. et al. Objective evaluation of surface- and satellite-driven carbon dioxide atmospheric inversions. Atmos. Chem. Phys. 19, 14233–14251 (2019).

54. 54.

C3S. ERA5-Land Hourly Data from 1981 to Present C3S ERA5-Land reanalysis (2019); https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-land?tab=overview

55. 55.

Gruber, A., Scanlon, T., van der Schalie, R., Wagner, W. & Dorigo, W. Evolution of the ESA CCI Soil Moisture climate data records and their underlying merging methodology. Earth Syst. Sci. Data 11, 717–739 (2019).

56. 56.

Hersbach, H. et al. The ERA5 global reanalysis. Q. J. R. Meteorol. Soc. 146, 1999–2049 (2020).

57. 57.

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

58. 58.

Zhang, Y. et al. A global moderate resolution dataset of gross primary production of vegetation for 2000–2016. Sci. Data 4, 170165 (2017).

59. 59.

Joiner, J. et al. Estimation of terrestrial global gross primary production (GPP) with satellite data-driven models and eddy covariance flux data. Remote Sens. 10, 1346 (2018).

60. 60.

Kim, H. J. Global Soil Wetness Project Phase 3 Atmospheric Boundary Conditions (Experiment 1) (2017); https://doi.org/10.20783/DIAS.501

## Acknowledgements

This research was funded by a Postdoc.Mobility fellowship of the Swiss National Science Foundation (P400P2_180784). C.F. acknowledges funding through NASA IDS grant 80NSSC17K0687. P.G. acknowledges funding from NASA 80NSSC18K0998 and European Research Council synergy grant USMILE ERC CU18-3746. P.C. acknowledges funding from the ANR CLAND convergence institute. S.I.S. acknowledges partial support from the European Union’s Horizon 2020 Research and Innovation Program (grant agreement 821003 (4C)). We thank all modelling groups who participated in the GLACE-CMIP5 experiments and conducted the model runs, in particular F. Cheruy, S. Hagemann and D. Lawrence. We also thank G. Bonan, J. K. Green, M. Hirschi, D. Lawrence, D. Miralles, U. Weber and Y. Yin for comments on the analyses, data availability or the manuscript.

## Funding

Open access funding provided by Max Planck Society.

## Author information

Authors

### Contributions

V.H. designed and conducted the study. S.I.S. designed and coordinated the GLACE-CMIP5 climate model experiment. A.B., P.C., P.G., M.J., M.R., S.I.S. and C.F. provided feedback on the analyses, the figures and the manuscript.

### Corresponding author

Correspondence to Vincent Humphrey.

## Ethics declarations

### Competing interests

The authors declare no competing interests.

Peer review information Nature thanks Belinda Medlyn and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

## Extended data figures and tables

### Extended Data Fig. 1 Soil moisture treatments in CTL and experiment A simulations.

At each grid point, the seasonal cycle calculated from the control experiment (CTL) is prescribed into the factorial experiment (ExpA). These example times series are taken from the CCSM4 model at 2° N and 58° W (northeast Amazon region).

### Extended Data Fig. 2 Concept of GLACE.

Setup of the control simulation (left) and the experiment with prescribed seasonal soil moisture (right). Sea ice, land use, and atmospheric CO2 concentration also prescribed from CMIP5 in both experiments.

### Extended Data Fig. 3 Temperature and VPD extremes influenced by LAC.

a, b, Change in the 95th percentile between the distributions of de-seasoned and de-trended temperature (a) and VPD (b) between CTL and experiment A ($$\varDelta {Q}_{95}={Q}_{95}^{{\rm{ExpA}}}-{Q}_{95}^{{\rm{CTL}}}$$). The median ΔQ95 of all four models is reported. Suppressing non-seasonal soil moisture variability in experiment A reduces temperature and VPD extremes, demonstrating the role of LAC.

### Extended Data Fig. 4 Change in annual NBP variability between CTL and experiment A.

Evaluation of the latitudinal change in NBP standard deviation (Δσ) between CTL and experiment A, decomposed by meteorological driver according to the sensitivity analysis. Negative values indicate a decrease of the NBP variability in experiment A compared to CTL. The middle and right columns indicate how much of this change is due to a change in the variance of the meteorological driver between experiment A and CTL, or due to a change in the sensitivity of NBP to that driver respectively (also see equation (1)). The results for each model are normalized by the model’s NBP standard deviation (calculated across the entire space–time domain) and the median across models is depicted. Black dots indicate that at least one model disagrees on the sign of the change.

### Extended Data Fig. 5 NBP anomalies in CTL and experiment A.

Distributions (all grid points, at all time steps) of modelled NBP anomalies (left column) and their decomposition into meteorological drivers with the sensitivity analysis (other columns) for the control experiment (CTL) and the experiment with only seasonal soil moisture (ExpA). Rows correspond to each of the four climate models. Note the logarithmic scale on the vertical axis. By construction, there are no soil moisture-driven NBP anomalies in experiment A in the second column (because seasonal soil moisture is prescribed in this experiment). The magnitude of the temperature-driven and VPD-driven NBP anomalies (third and fourth columns) is substantially reduced in experiment A (as a result of soil moisture–atmosphere feedback mechanisms).

### Extended Data Fig. 6 Comparison of direct versus indirect effects.

Difference between the magnitudes of direct effects (Fig. 2b) versus indirect (feedback) effects occurring through temperature and VPD (Fig. 2c).

### Extended Data Fig. 7 Opposing perspectives on drivers of NBP IAV reconciled by soil moisture–atmosphere feedback.

a, Relative magnitude of individual NBP components across spatial scales (same as Fig. 3b). b, c, The apparent relative importance of the meteorological drivers depends on how the indirect effects of soil moisture on temperature and VPD are viewed. Outer uncertainty bounds indicate the model spread (ensemble mean ±1σ), inner uncertainty bounds indicate the spread (±1σ) in the tendency (that is, the relative change from local to global scale; see Methods).

### Extended Data Fig. 8 Sensitivity analysis compared to observational estimates.

a, Contribution of different meteorological drivers to global NBP IAV as estimated from the control simulations and from two independent observational products. Here, NBPT&VPD is not separated into a LAC and non-LAC contribution as in Fig. 3b (because this cannot be done with the observational datasets). b, Same as Fig. 3b, but based on FLUXCOM. c, Same as Fig. 3b, but based on CAMS.

### Extended Data Fig. 9 Contribution of LAC hotspots to global NBP IAV.

Global NBP IAV from the control experiment (CTL) calculated over all land grid cells versus only over the LAC hotspots identified in Fig. 4.

### Extended Data Fig. 10 Tropical temperature in CTL versus experiment A.

a, IAV in tropical (24° N to 24° S) mean land temperature in model experiments with and without variability in soil moisture (similar to Fig. 1a for NBP). b, Apparent sensitivity of global mean NBP to tropical mean temperature in CTL and experiment A.

## Supplementary information

### Supplementary Information

This file contains Supplementary Methods, Supplementary Tables 1–4 and Supplementary Figures 1–16.

## Rights and permissions

Reprints and Permissions

Humphrey, V., Berg, A., Ciais, P. et al. Soil moisture–atmosphere feedback dominates land carbon uptake variability. Nature 592, 65–69 (2021). https://doi.org/10.1038/s41586-021-03325-5

• Accepted:

• Published:

• Issue Date: