Climate-driven vegetation greening further reduces water availability in drylands

2.1 Vegetation dynamics

Three long-term satellite-based LAI datasets were applied to analyze changes in global terrestrial vegetation during 1982–2016 derived from the Global Inventory Modeling and Mapping Studies (GIMMS3g v1) LAI (http://sites.bu.edu/cliveg/datacodes/; Zhu et al., 2013), the Global Land Surface Satellites (GLASS v3) LAI (http://www.glass.umd.edu/; Xiao et al., 2016), and the Global Mapping (GLOBMAP v3) LAI (https://zenodo.org/record/4700264; Liu et al., 2012). We first composited the original data (8-day or 15-day temporal resolution depending on the product) to the monthly temporal resolution by averaging the composites in the same month. Before the averaging process, we weighted each composite by the number of days that the composites fall in the given month. Subsequently, monthly LAI data were linearly resampled to the 0.5° spatial resolution. In this study, the growing season-averaged LAI was used as a proxy of vegetation growth. To this end, we identified the climatological growing season based on Vegetation Index and Phenology (VIP) satellite-based product (https://vip.arizona.edu/vipdata/V4/DATAPOOL/PHENOLOGY/; Zhang et al., 2014). The VIP product provides the date of onset and end of growing season for the period 1981–2016, on the basis of a piecewise logistic model-based land surface phenology (LSP) detection algorithm in combination with daily enhanced vegetation index generated from AVHRR LTDR and Moderate Resolution Imaging Spectroradiometer (MODIS) CMG records. We first resampled the original VIP dataset (0.05°) to 0.5° spatial resolution, and further calculated the multiyear mean onset date and end date of growing season at the grid-cell scale. Subsequently, we identified which months belong to the growing season for each vegetated grid cell (Figure S1), and used as the climatological reference period to derive a multiyear time series of growing season LAI, similar to Forzieri et al. (2017). Considering a fixed climatological growing season at pixel level over the entire observational period 1982–2016 preserves the growing season LAI retrievals from possible effects associated with inter-annual changes in the phenological cycle. We derived a global vegetated land mask from the MODIS land cover map (MCD12Q1, https://lpdaac.usgs.gov/products/mcd12q1v006/) for the year 2001 (Friedl & Sulla-Menashe, 2019), which is the starting year of this product and is the closest year to the start of our observational period (1982). All grid cells (0.5° × 0.5° resolution) for which the MCD12Q1 dominant plant functional type is vegetation class (including tree, shrub, grass, and cropland) were defined as vegetated areas and were included in our analyses (Figure S2). For each vegetated grid cell, we retrieved the temporal trend in growing season-averaged LAI (δLAI) for the three satellite-based LAI datasets, and evaluated the significance using the nonparametric Mann–Kendal test (Figure 1). Besides, considering that grid cells have different area extents depending on their latitude, we computed the global mean LAI and its trend by the weighted average method based on the grid-cell area (Figure 2a). The same methodology is applied consistently for global-scale aggregated metrics described in the following sections.

We recognized that changes in LAI can only partially reflect the vegetation control on water and energy exchanges between land and atmosphere. However, we argued that other important variables, such as surface albedo and emissivity, are strongly driven by vegetation density and therefore correlated with LAI, as demonstrated in Text S1 and Figure S3. Their inclusion in our modeling framework as additional vegetation predictors would bring issues of causality between predictors and target variables and limited independent information.

To characterize the spatial distribution of rainfed agricultural lands, we used the SPAM dataset (https://www.mapspam.info/data/; Yu et al., 2020), which provides cropland maps for different farming systems (rainfed and irrigation) for three different epochs, 2000, 2005, and 2010. Grid cells with <10% of the irrigated area in all three periods were considered as the rainfed cropland and were used for further analysis.

2.2 Climate drivers

To explore our results under different climatic conditions, the global vegetated areas were classified into humid, semihumid, semiarid, and arid by adopting the climate classification of the United Nations Environment Program (United Nations Environment Program, 1997). Such classification is based on the aridity index, which is defined as the ratio of potential evapotranspiration to precipitation (i.e., ET_p/P; Figure S4). The thresholds between humid and semihumid classes, between semihumid and semiarid classes, and between semiarid and arid classes were set to ET_p/P = 1, 2, 5, respectively (Ukkola et al., 2016). Drylands in this study are defined as regions with ET_p/P ≥ 2 (i.e., sum of semiarid and arid classes), slightly more strict compared with related literature (ET_p/P ≥ 1.54; Huang et al., 2016; Zhou et al., 2021).

2.3 Evapotranspiration products

Evapotranspiration estimates were derived at the global level for the period 1982–2016 from five different observation-based datasets including the Breathing Earth System Simulator (BESS; http://environment.snu.ac.kr/bess_flux/; Jiang & Ryu, 2016), the FLUXCOM Model Tree Ensemble (https://www.bgc-jena.mpg.de/geodb/projects/Home.php; Jung et al., 2019), the Global Land Data Assimilation System (GLDAS v2.0 + v2.1; https://ldas.gsfc.nasa.gov/gldas; Rodell et al., 2004), the Global Land Evaporation and Amsterdam Model (GLEAM v3.3a; https://www.gleam.eu/; Martens et al., 2017; Miralles et al., 2011), and the Process-based Land Surface Evapotranspiration/Heat Fluxes (PLSH; http://files.ntsg.umt.edu/data/; Zhang et al., 2015). These ET products are provided at the monthly scale, but with distinct spatial resolutions, ranging from 0.0833° (PLSH) to 0.5° (BESS and FLUXCOM), and thus were linearly resampled to the common spatial resolution 0.5°. To cover the whole period, FLUXCOM dataset used here is the ensemble mean that encompasses members of different machine learning methods and energy balance corrections within the “RS + METEO” setup. Besides, the GLDAS dataset used here was generated by merging two separate product versions v2.0 (time coverage 1948–2014) and v2.1 (time coverage 2000–present). In order to derive a GLDAS ET product covering our entire reference period (1982–2016), we built a regression relationship between the two datasets over their overlapping period (2000–2014) at the grid scale, and subsequently, extended the GLDAS v2.1 time series back to 1982 using such relationship with GLDAS v2.0. We assessed the consistency of the integrated dataset (i.e., GLDAS v2.1 and its extension) with GLDAS v2.0 for the period 1982–2014, as shown in Figure S5.

The above-mentioned ET products have been largely validated in previous studies (Bai & Liu, 2018; Jiang & Ryu, 2016; Jung et al., 2010; Miralles et al., 2011; Zhang et al., 2015). To further increase the confidence in the use of such datasets in our global modeling framework, we additionally evaluated their performance against surface observations in a dedicated experiment (for details, see Text S2). Results show high consistency between simulated and observed ET values and confirm the validity of the ensemble of global ET products used in this study (Figure S6).

Moreover, GLEAM provides in addition to ET estimates, its key components including transpiration (E_c), canopy interception loss (E_i), and soil evaporation (E_s). Estimates of ET components from GLEAM have been well validated in previous works against ground observations (Miralles et al., 2010; Talsma et al., 2018; Wei et al., 2017) and are therefore used in this study to explore the ET-related physical and biological responses to environmental changes.

2.4 Simulations of the LSMs

We used monthly LAI simulations for the period 1982–2016 from an ensemble of 14 LSMs completed within the TRENDY v7 project (Sitch et al., 2015), including CABLE-POP (Haverd et al., 2014), CLASS-CTEM (Melton et al., 2019), CLM5.0 (Lawrence et al., 2019), DLEM (Tian et al., 2010), ISAM (Jain & Yang, 2005), JSBACH (Raddatz et al., 2007), JULES (Best et al., 2011), LPJ-GUESS (Smith et al., 2001), LPX (Keller et al., 2017), OCN (Zaehle et al., 2010), ORCHIDEE (Krinner et al., 2005), SDGVM (Woodward & Lomas, 2004), SURFEX (Le Moigne et al., 2020), and VISIT (Ito, 2008). All model outputs were resampled to the common 0.5° spatial resolution. In order to disentangle the effect of climate change to LAI, we analyzed outputs from three sets of factorial simulation experiments: (S1) varying CO₂ only; (S2) varying CO₂ and climate; and (S3) varying CO₂, climate, and land use. For models that include the nitrogen cycle, time-varying nitrogen inputs were used for all three experiments.

Multiple observational datasets were used to force TRENDY LSM simulations, including time series of atmospheric CO₂ concentration (Keeling & Whorf, 2005), historical climate fields from the climatic research unit (CRU) observational dataset and its derivative (Harris et al., 2014), and land cover change data from the land-use harmonization product (LUH2; Chini et al., 2021). Additional details on the forcing input data used in the TRENDY v7 project were reported in Text S3. Furthermore, simulated LAI under the S3 scenario were compared with satellite observations to evaluate model performance. The overall high consistency (Figures S7 and S8) supports the use of the ensemble of LSMs to explore long-term changes in LAI and identify its key drivers, as also performed in previous studies over shorter time periods (Piao et al., 2015; Zhu et al., 2016).

2.5 Retrievals of climate-driven greening through data model integration

In total, the approach described above produced a 42-member ensemble of δLAI^CLI estimates, each one derived from a different combination of satellite LAI product and LSM (Figure S10; Table S1). The term δLAI^CLI excludes the effect of CO₂ fertilization, land cover change, and nitrogen deposition, as these components have been removed from the factorial simulations performed on the TRENDY runs (term δLAI^S2-δLAI^S1). Consequently, CO₂ fertilization, land cover change, and nitrogen deposition disentangled in Equations (2 and 3) do not influence the following analyses.

2.6 Sensitivity of ET to LAI changes

We applied Equation (4) to each unique combination of ET and LAI dataset, therefore resulting in a 15-member ensemble of ET sensitivity estimates (Figure S12). We also obtained the corresponding sensitivities of ET to changes in climate drivers (Figures S13–S15), which have been used to quantify the direct effect of climate change (see next section). Furthermore, based on the estimates of ET components provided by GLEAM v3.3a, we applied Equation (4) separately to ∆E_c, ∆E_i, and ∆E_s, and thus calculated the sensitivity of different ET components to changes in growing season-averaged LAI (Figures S16–S18). A Student's t-test was performed at the grid scale to determine the statistical significance of sensitivity of ET and its components. In order to properly represent the generality of relationships between sensitivities of ET and local environment conditions, we performed binned average analysis to minimize the uncertainty induced by spatial heterogeneity (difference in topography, soil property, and so on) that randomly affects sensitivities of ET and its components to LAI changes (Figure 3c,d; Figure S19).

2.7 Attribution of ET variations to climate-induced changes in LAI

The use of three satellite-observed LAI records, five ET products, and 14 model simulations leads to a resulting ensemble of 210 members for the term ($$; Figure S20), which makes possible to quantify the different sources of uncertainties in the derived signal. Following an analogous approach, trend in different ET components (i.e., E_c, E_i, and E_s) due to climate-induced LAI changes was calculated (Figure S21). The significance of such indirect effect of climate change is determined by whether sensitivity to LAI and climate-induced trend in LAI are both significant (p < .05) at the grid-cell scale.

2.8 Direct and indirect effects of climate change on WA

The fraction of the overall trend in climate-induced WA explained by a given term (P, direct climate effect on ET, indirect climate effect on ET) is then quantified as the ratio between the term-specific trend and the overall trend (Figure S25). The residual term (δWA − δWA^CLI) represents the shift in WA induced by other factors (δWA^OF), that is, changes in vegetation driven by elevated CO₂ concentration, land cover change, and nitrogen deposition (Figure S24e,f), which have been excluded in the previous analyses (see Section 2.5). Direction and magnitude of climate change effect on WA through multiple pathways were investigated across various climatic conditions (Figure 5). Furthermore, changes in WA were analyzed with a particular focus on global rainfed croplands (especially in drylands) to offer more detailed insights into the indirect effect of climate change (through LAI) on WA for food production, based on the cropland mapping derived from SPAM dataset (Figure 6).

2.9 Uncertainty analysis

We performed three experiments (i.e., Exp1, Exp2, and Exp3) to identify the sources of uncertainties, based on a total 210-member ensemble of indirect climatic effect on ET (i.e., $$, Figure S20). More specifically, Exp1, Exp2, and Exp3 enabled to quantify the uncertainty originated from the satellite LAI product, the ET product, and the LSM simulation, respectively (Figures S26 and S27). As an example for Exp1, separately for each satellite LAI product, we computed averaged $$ across multiple LSMs and ET products. We thus obtained three estimates of $$ (one for each satellite LAI product) and computed the average and standard deviation over the resulting three-member ensemble. Such standard deviation reflects the uncertainty originated from the satellite LAI product. Similar to the calculation steps in Exp1, but for Exp2, we computed averaged $$ across multiple LAI products and LSMs, separately for each ET product, and thus obtained five estimates of $$; meanwhile, for Exp3, we computed averaged $$ across multiple LAI and ET products, separately for each LSM, and thus obtained 14 estimates of $$. For the total uncertainty, the standard deviation was computed over the whole 210-member ensemble (Figure S26d).

3.1 Climate-induced variations in LAI

During the period 1982–2016, global vegetated land has experienced a widespread positive trend in LAI with a multiple product mean (MPM) of 4.03 ± 0.88 10⁻³ m² m⁻² year⁻¹ (δLAI, Figures 1d and 2a). Such greening signal is statistically significant (Mann–Kendall test, p < .05) over ~44% of the vegetated lands and consistent across single LAI products (47.3%, 44.9%, and 36.7% for GIMMS3g, GLOBMAP, and GLASS, respectively, Figure 1). Regions with the most notable increase in vegetation density (>20 10⁻³ m² m⁻² year⁻¹) are observed in southeastern China, Europe, Sahel, and eastern United States.

Factorial experiments based on 14 LSMs from the TRENDY project (see Section 2) show that all models, except CABLE-POP, agree that climate change has had an overall positive effect on vegetation growth, even though the inter-model spread is considerable (global average signal ranging from 0.62 ± 0.12 to 2.71 ± 0.48 10⁻³ m² m⁻² year⁻¹, Figure 2a; Figure S10). Results of multimodel ensemble mean (MMEM) show that the trend in LAI attributed to climate change (δLAI^CLI) is 1.35 ± 1.50 10⁻³ m² m⁻² year⁻¹, which is about 33.5% of the averaged observed trend in LAI (MPM δLAI). The contribution of climate change to the greening signal seems stronger in cold regions where low temperatures are a key limiting factor for vegetation growth (McManus et al., 2012; Myers-Smith et al., 2020). For instance, when focusing on northern latitudes (60°N–90°N), vegetation trend induced by climate change reaches 2.25 ± 2.17 10⁻³ m² m⁻² year⁻¹ (MMEM), equivalent to 68.2% of the observed overall trend (Figure 2b), indicating that climate change dominates the boreal vegetation dynamics over the past few decades. The strong positive effects of climate change occur in northern Canada, Scandinavia, and Siberia (Figure 2c). It is worth noting that, despite the negative impact of climate change in southeastern United States, Amazon, India, and southern China (Figure 2c), an increasing overall trend in LAI can still be observed in these regions (Figure 1), proving that the positive effects of CO₂ fertilization, together with changes in land cover, have effectively compensated the climate-induced decline of LAI.

3.2 Sensitivity of ET to variations in LAI

In order to disentangle the indirect effect of climate change on water resources, we first quantified the overall sensitivity of ET to changes in vegetation density (∂ET/∂LAI) within a multiple linear regression framework by ruling out the contribution of climate factors (Equation 3, see Section 2). At the interannual timescale, we found that ET is positively related to variations in LAI (∂ET/∂LAI > 0) across most of the global vegetated land (88.2%) with an average magnitude of 64.80 ± 41.94 mm/m⁻² m⁻² (Figure 3a). In particular, ∂ET/∂LAI reaches the highest value (>200 mm/m⁻² m⁻²) in desert margins, but is relatively low (even negative) in rainforests. Results obtained separately from five different ET products (BESS, FLUXCOM, GLDAS, GLEAM, and PLSH) show consistent patterns of ∂ET/∂LAI (Figure S12). A comprehensive set of modeling experiments was additionally performed to further test the validity of our approach and the robustness of ∂ET/∂LAI retrievals with respect to the potential impacts of nonlinear relationships, collinearity, and discrepancy in input datasets, which was carefully illustrated in Text S4 and S5, and meanwhile briefly presented in Section 4.

To better clarify the climate control on different ET-related biological and physical processes (Jasechko et al., 2013; Lian et al., 2018), we analyzed the sensitivity of ET to LAI changes separately for E_c, E_i, and E_s. We found that E_c has the highest sensitivity to LAI changes at almost all latitudinal bands, especially in the subtropics where it is about 30-fold more sensitive than E_i and E_s (Figure 3b). Not only that, as LAI increases, the resulting increase in E_c at the expense of the soil component makes E_c more sensitive to LAI changes than the total ET, except for the tropics where the sensitivity of E_i is dominant because of the high frequency of rainfalls combined with the high LAI which leads to relevant fraction rain interception. Under this context, E_c component exhibits the largest average contribution to global positive ∂ET/∂LAI (100.1%), followed by E_i (4.6%) and E_s (−5.5%). These results suggest that the positive sensitivity of ET to vegetation change is mostly driven from the strong feedback of its biological processes (i.e., stomatal control) rather than the physical ones (i.e., evaporation of intercepted rain or soil moisture), highlighting the critical role of LAI and E_c in land–atmosphere coupling.

The sensitivity of ET to LAI, and in particular the E_c component, shows a clear dependence on local environment gradients of aridity (ET_p/P: ratio of potential evapotranspiration to precipitation) and mean annual LAI (Figure 3c,d; Figure S19). The significant linear positive relationship between the ∂ET/∂LAI and the ET_p/P (R² = .86, p < .01) reveals an enhanced vegetation control on the hydrological cycle as climate becomes drier (Figure 3c). ∂ET/∂LAI declines exponentially with the increase of LAI (R² = .98, p < .01) and gradually approaches zero, suggesting that change in vegetation density has a fairly limited impact on ET in regions with high LAI (Figure 3d). As a consequence, large-scale vegetation changes over drylands with low-density vegetation have the potential to affect available water via ET to a much greater extent compared to regions with wetter climate or denser vegetation where the water balance is not so strictly under biological control. This finding highlights the relevance of the land–climate interplay and therefore of vegetation management on the availability of water resources and the sustainability of societies in arid regions.

3.3 Indirect effect of climate change on ET

The indirect effect of climate change on ET through adjustments in terrestrial vegetation ($$) was quantified by multiplying the observed sensitivity of ET to variations in LAI by the climate-induced trend in LAI (Equation 5, see Section 2). Based on the combination of three satellite-based LAI retrievals, five observation-driven ET products, and 14 LSMs, we derived a 210-member ensemble of $$ estimates that were used to account for the structural uncertainty originating from the different datasets and land surface schemes and to quantify their marginal effect (see Section 2). Globally, climate change has caused substantial variations in terrestrial vegetation over the past 35 years and increased indirectly ET by an average rate of 0.051 ± 0.067 mm year⁻² (Figure 4a) which corresponds to 8.9% of the overall positive trend in ET. Such indirect effect of climate change mediated by vegetation is of comparable magnitude to the direct effect originating from temperature, precipitation, and radiation ($$) which has led to an increase in ET of 0.129 ± 0.128 mm year⁻² (Figure S22). This highlights the relevance of the indirect effect via associated vegetation change in the long-term climate influence on the water cycle and the complex feedback mechanisms in the Earth system among climate change, greening trends, and variation of the water cycle. The increase in ET due to the climate-induced vegetation changes is larger in low vegetated lands (e.g., inland western United States, Mediterranean regions, Sahel, southern Africa, northern China), thanks to the combination of high increase in LAI (Figure 2c) and of the high ET sensitivity to vegetation changes (Figure 3a,d). On the contrary, decrease in ET through climate impacts on vegetation is found in some limited regions in central Asia as a consequence of the vegetation degradation (Figure 1) following prolonged drought stress (Piao et al., 2011; Zhang et al., 2016).

We found that the indirect climate effect is prominently driven by the climate-induced change in transpiration (E_c) which appears the dominant factor over 70.0% of the global vegetated land (Figure 4b; Figures S21 and S28) regardless of the climatic conditions (Figure S29) and clearly reflects a higher sensitivity of biological processes to changes in evaporative surface (Figure 3b). Nevertheless, exceptions occur in specific regions such as in tropical forests (e.g., Amazon, Congo, and Southeast Asia), where changes in canopy interception loss (E_i) associated with climate-induced variations in vegetation explain the largest fraction of the trend in ET (Figure 4b; Figures S21 and S28), largely driven by dominant sensitivity of E_i (Figure 3b). This occurs in regions with frequent rainfalls and high LAI like the humid tropics (Figure S19), where vegetation typically shows higher water interception efficiency and storage capacity (Wallace et al., 2013).

Considering that the direct and indirect effects of climate change on the water cycle are directional, their integrated effects on variations in ET can be additive (i.e., leading to an amplification of the signal) or offsetting (i.e., leading to a dampening of the signal). Overall, 53.4% of the global vegetated land exhibits additive effects (i.e., “+ +” and “– –” in Figure 4c) and 46.6% shows offsetting effects (i.e., “– +” and “+ −”). By considering the indirect effect from climate-induced vegetation changes, we consequently quantified the net effect of climate change on ET. In this respect, our study found that in 23.9% of the global vegetated land, despite an increase in ET due to the direct climate effect, this positive signal can be partly mitigated by the opposite effect of climate-induced vegetation change (i.e., “+ −,” Figure 4c), with an average decrease of 36.4% (Figure 4d). Climate warming and the associated increase in atmospheric moisture demand are the primary causes of the initial positive effect on ET in these regions (Figures S22a and S23a), while precipitation remains stable or even decreases (Figure S23b). The resulting vegetation browning mitigates the initial positive effect of higher temperature on ET (Figure 2c). Over a larger fraction of the globe (“+ +”: 36.8%, Figure 4c), the indirect effect tends to exacerbate the direct signal, leading to a further increase of 39.6% (Figure 4d). This amplification effect on ET occurs in western China, Europe, and Siberia, as well as in drylands (ET_p/P ≥ 2, see Figure S4) of southern Africa and northwestern India.

3.4 Implications for WA

Climate change, by altering both water supply and loss, influences regional WA and a range of water-dependent socio-economic sectors whose sustainability may be at risk in the decades to come (Greve et al., 2018; Zhou et al., 2021). In order to provide a comprehensive assessment on the hydrological implications of climate change, we disentangled the climate-induced trend in WA (δWA^CLI) originating from direct ($$) and indirect ($$) pathways (Equation 8, see Section 2). We found that in the past three and a half decades, climate change has led to a decline in WA with an average rate of −1.440 ± 0.178 mm year⁻², explaining 81.9% of total reduction (Figure S24). The residual term (i.e., 18.1%) is attributable to changes in vegetation driven by other factors than climate change, such as elevated CO₂ concentration, land cover change, and nitrogen deposition. We noticed that the climate-induced reduction in WA appears primarily caused by a decrease in precipitation (average contribution: 87.5%), and moderately driven by the ET response to climate change through direct and indirect pathways (average contributions of 9.0% and 3.6%, respectively; Figure S25). However, for a considerable part of global vegetated land (~20.0%) that has experienced low trend in precipitation (between −0.4 and 0.8 mm year⁻²), both climate-induced direct and indirect effects on ET play a far larger role in influencing WA (Figure 5a,b). For example, in regions with trend in precipitation ranging from −0.4 to 0 mm year⁻², variations in ET directly and indirectly related to climate change provide higher contributions to the climate-induced decrease in WA (29.8% and 11.5%, respectively), further amplifying the negative effect of precipitation decrease (gray bars in Figure 5a,b). These results reinforce the importance of accounting for the climate–vegetation interactions to properly assess the overall effect of climate change on WA, particularly in areas characterized by limited shifts in precipitation regime, such as arid (e.g., central Asia and Australia) and boreal regions (e.g., northern Canada and Siberia; Figures S23 and S25).

To further illustrate the underlying mechanisms which modulate the impact of climate change on WA, we averaged the marginal contributions of direct and indirect effects of climate change across gradients of precipitation change (δP) and aridity (ET_p/P). We found that the contribution of the indirect climate effect on WA becomes larger under drier climate (Figure 5d)—where both sensitivity of ET to vegetation alterations (Figure 3c) and resulting impact (Figure 4a) are larger—whereas the direct climate effect exhibits an opposite relationship with ET_p/P (Figure 5c). These patterns suggest that the local background climate strongly regulates the effect of ET change on WA. More importantly, the high positive contribution of the indirect climate effect on WA under high ET_p/P (Figure 5d) indicates that climate-induced vegetation change tends to exacerbate the declining WA in drylands, accelerating the water depletion of ecosystems that already suffer from water stress. Specifically, across the drylands with slight decrease in precipitation (δP: −0.4–0 mm year⁻²), climate-induced vegetation change amplifies the decline in WA by 22.0% (Figure 5d).

Such exacerbation in declining WA is evident in 42.9% of global land dominated by rainfed croplands, and a considerable fraction of them are located in drylands (e.g., Sahel, dryland belt of northern Eurasia, and Australia; Figure 6a,c). Obviously, climate-induced vegetation change and associated ET change are more menace to water resources for rainfed croplands in drylands with the average magnitude of −0.514 ± 0.467 mm year⁻², resulting in a further 12.6% decrease in WA for local agricultural production (Figure 6b,d). Especially for those with slight decrease in precipitation (δP: −0.4–0 mm year⁻²), indirect effect of climate change on ET contributes 69.9% of decrease in regional WA, and therefore poses an additional and tremendous pressure to agricultural water supply (Figure 6d).