CO2 fertilization is spatially distinct from stomatal conductance reduction in controlling ecosystem water-use efficiency increase

It is well known that global ecosystem water-use efficiency (EWUE) has noticeably increased over the last several decades. However, it remains unclear how individual environmental drivers contribute to EWUE changes, particularly from CO2 fertilization and stomatal suppression effects. Using a satellite-driven water–carbon coupling model—Penman–Monteith–Leuning version 2 (PML-V2), we quantified individual contributions from the observational drivers (atmospheric CO2, climate forcing, leaf area index (LAI), albedo and emissivity) across the globe over 1982–2014. The PML-V2 was well-calibrated and showed a good performance for simulating EWUE (with a determination coefficient (R 2) of 0.56) compared to observational annual EWUE over 1982–2014 derived from global 95 eddy flux sites from the FLUXNET2015 dataset. Our results showed that global EWUE increasing trend (0.04 ± 0.004 gC mm−1 H2O decade−1) was largely contributed by increasing CO2 (51%) and LAI (20%), but counteracted by climate forcing (−26%). Globally, the CO2 fertilization effect on photosynthesis (23%) was similar to the CO2 suppression effect on stomatal conductance (28%). Spatially, the fertilization effect dominated EWUE trend over semi-arid regions while the stomatal suppression effect controlled over tropical forests. These findings improve understanding of how environmental factors affect the long-term change of EWUE, and can help policymakers for water use planning and ecosystem management.


Introduction
The arising atmospheric CO 2 due to anthropogenic CO 2 and other greenhouse gases emissions since 1860s (IPCC 2013) have induced multiple changes in biogeophysical and biochemical processes including surface energy balance, terrestrial vegetation dynamics, and coupled carbon-water cycle (Ciais et al 2013, Green et al 2017, Jung et al 2017, Bonan and Doney 2018. The most remarkable aspect is that the greenhouse effect-induced global warming and related environmental changes (e.g. geographical patterns of global precipitation) are expected to increase drought (Dai 2012) and evapotranspiration (ET) (Jung et al 2010, Zhang et al 2012, 2021, Ma et al 2020, accelerate global water cycle (Gedney et al 2006, Huntington 2006, Betts et al 2007, Mankin et al 2019, and significantly influence vegetation water-use efficiency (e.g. Franks et al 2013, Keenan et al 2013, gross primary productivity (GPP) and land carbon sink (e.g. Ciais et al 2005, Sitch et al 2015. Owing to the increased atmospheric CO 2 , surface climate change, and vegetation structure change (measured by leaf area index (LAI) in this study), global ecosystem water-use efficiency (EWUE) have experienced a significant long-term increase, especially since 1980s. The EWUE, defined as the ratio of GPP to ET (Law et al 2002), is a critical metric of potential water costs to maintain a given rate of photosynthesis, providing important information of how terrestrial ecosystem interacts with global carbon and water cycles (Franks et al 2013, Lavergne et al 2019. Within global terrestrial ecosystems, forest EWUE plays the most important role. Eddy covariance data-based observations (Keenan et al 2013), tree-ring records (van der Sleen et al 2014) and observation-driven modeling studies (Huang et al 2015, Cheng et al 2017, Zhou et al 2017 demonstrated that the global forest EWUE has increased significantly in response to elevated CO 2 over the last three decades (Lavergne et al 2019). In particular, much longer tree-ring records also indicated that the increased CO 2 has induced substantial increases in forest EWUE over the last century (Frank et al 2015, Dekker et al 2016, which was about two times larger than that simulated by the Coupled Model Intercomparison Project Phase 5 (CMIP5) Earth system models (Dekker et al 2016). In addition, studies also indicated that the EWUE in global grassland ecosystems also have changed significantly. Based on the Moderate Resolution Imaging Spectroradiometer (MODIS) products, Gang et al (2016) found that all grassland types has globally experienced an increase in 2000-2013, but regionally about 40% of grassland ecosystems showed a decreased EWUE. Yang et al (2022) found that annual EWUE in Northwest China (mainly covered by desert, steppe and meadow etc) showed a significant increasing trend over 1982-2015, which is primarily dominated by GPP increase (88%). Another recent analysis showed that the alpine grass EWUE also has an increasing trend in the Qinghai-Tibet plateau in China since 2000s (Ji et al 2022).
However, it remains poorly understood that how the long-term increase in EWUE was contributed individually or interactively by major forcing factors (e.g. atmospheric CO 2 , climate change, LAI) since 1980s. Particularly, how CO 2 fertilization-induced change in GPP and CO 2 stomatal suppressioninduced change in ET contributed separately to the change in EWUE remains unknown. On one hand, the CO 2 fertilization effect stimulated plant leaf photosynthesis and enhanced vegetation productivity (i.e. GPP) at an ecosystem scale (Ainsworth and Long 2005, Schimel et al 2015, Wang et al 2020. On the other hand, increased CO 2 also induced the suppression effect on plant/canopy transpiration (E c ) by reducing leaf/canopy stomatal conductance which regulates exchanges of both CO 2 and H 2 O molecules between the plant and the atmosphere (Field et al 1995, Katul et al 2009, Damour et al 2010, Medlyn et al 2011, Swann et al 2016, Zhang et al 2021. Yang et al (2018) attributed ET changes using CMIP6 models and showed that an increase in ET caused by a warming-induced vapor pressure deficit increase is almost entirely offset by a decrease in ET caused by stomatal conductance reduction (or resistance increase) driven by rising CO 2 concentration. Using a remote sensing-based ET method, Zhang et al (2021) showed that the increased CO 2 -induced canopy stomatal conductance reduction can largely offset the strong greening-induced increase in ET and E c over Eurasia since 1980s. To identify such an effect on EWUE, change in the plant water-use efficiency (PWUE), i.e. the rate of carbon uptake via GPP per unit of water loss through E c , should also be attributed.
The LAI, due to the combined effect from increased CO 2 fertilization effect, climate change, land use change, and nitrogen deposition, has increased globally and was quite evident by satellitebased observations since 1980s (Norby and Zak 2011, Zhu et al 2016, Chen et al 2019. The LAI is expected to increase E c and precipitation interception by vegetation (E i ), thereby boosting the rate of terrestrial ET (total ET is defined as the sum of soil evaporation E s , E c and E i ) (Zhang et al 2015, 2021, Lemordant et al 2018, Zeng et al 2018, which further reduces PWUE and EWUE by the increased E c and ET. In this study, to isolate the individual contributions of environmental drivers (i.e. climate forcing, atmospheric CO 2 , LAI, surface albedo and emissivity) to the long-term changes in EWUE and PWUE at global and regional scales, we conducted a group of sensitivity modeling experiments by applying a satellite data-driven water-carbon coupling process model. The experiments produced global simulations at a high spatial and temporal resolution (0.05 • × 0.05 • , half-monthly) over 1982-2014 under a series of scenarios that allow for one or more drivers vary over the historical period as forcing inputs to the model. Based on these modeling experiments, we aimed to (a) identify global/regional dominant drivers of the long-term trends in annual EWUE andPWUE over 1982-2014 and quantify the relative contributions of the increased CO 2 fertilization and stomatal suppression effects to the global/regional trends in EWUE and PWUE.

Satellite data-driven water-carbon coupling process model (PML-V2)
The version 2 of the Penman-Monteith-Leuning model (PML-V2) was used in this study. Compared to the older version PML-V1 that does not estimate GPP and CO 2 's impact on ET, the PML-V2 considers both effects from rising CO 2 , and increasing LAI and climate forcing on GPP and ET. Specifically, the PML-V2 was developed by coupling the widelyused photosynthesis model (Farquhar et al 1980) and a canopy stomatal conductance model (Yu et al 2004) with the Penman-Monteith energy balance equation (Monteith 1965) to jointly estimate GPP and ET (Leuning et al 2008, Zhang et al 2010, Gan et al 2018. The PML-V2 model has been applied to successfully produce the MODIS LAI-based global GPP and ET products at a 500 m and 8 d resolution since 2002, which were noticeably better than most widely used GPP and ET products (Zhang et al 2019). The PML-V2 was also used to explore the long-term impacts from increased CO 2 and strong greening on ET trends over Eurasia since 1980s (Zhang et al 2021).
In the PML-V2 model, transient GPP and canopy stomatal conductance to water vapor (G c ) are calculated by the photosynthesis process (A gs ) with the constraint of vapor pressure deficit (D) at surface, where f D is the D constraint function, D min , D max , D 0 are parameters to identify the constraint of atmospheric vapor pressure deficit, P 1 = A m βI 0 η, P 2 = A m βI 0 , P 3 = A m ηC a , P 4 = βI 0 ηC a , and A m = 0.5V m . I 0 is the photosynthetically active radiation (PAR, in mol) from shortwave downward radiation.C a is the atmospheric CO 2 concentration (in ppm or mol mol −1 ). V m,25 , β, η, m are key parameters for photosynthesis process, i.e. maximum photosynthetic carboxylation efficiency (V cmax ), initial photochemical efficiency (β), initial value of the slope of CO 2 response curve (η), and stomatal conductance coefficient (m) in PML-V2 model (see table S1 available online at stacks.iop.org/ERL/ 17/054048/mmedia). The PML-V2 estimates ET by separately calculating the ET's three components: canopy transpiration (E c ), soil evaporation (E s ) and rainfall interception (E i ), where E c is estimated using the Penman-Monteith energy balance equation but the canopy stomatal conductance (G c , in a unit of m s −1 ) is not a constant and changes through coupling with the photosynthesis process (equation (2)), The soil evaporation (E s ) depends on absorbed energy flux and soil water deficit, where the surface available energy (A = R n − G) is divided into canopy absorbed energy (A c ) and soil absorbed energy (A s ), ; R n is the net radiation and G is the ground heat flux (W m −2 ). ε = ∆ γ , and ∆ is the slope of the curve relating saturation water vapor pressure to temperature (kPa • C −1 ). ρ is the density of air (g m −3 ); c p is the specific heat of air at constant pressure (J g −1 • C −1 ); G a is the aerodynamic conductance (m s −1 ); f is a dimensionless variable that determines the water availability for soil evaporation. The PML-V2 model also simulates the E i based on a revised Gash-model scheme (van Dijk and Bruijnzeel 2001), which has been shown in supplementary text 1.
In PML-V2, 11 key parameters (e.g. maximum photosynthetic carboxylation efficiency (V cmax ), stomatal conductance coefficient (m), and specific water storage capacity (S leaf ), etc, see table S1) were well calibrated across ten plant functional types (PFTs) based on the NASA MCD12Q1 IGBP land cover classes (Friedl et al 2010) against global 95-site eddy covariance data-derived GPP and ET based on the FLUXNET2015 dataset and AVHRR GIMMS3g-based LAI over 2002-2014 (see table S1). The key parameters for each PFT were optimized by minimizing the objective function (F): F = 2 − (NSE ET + NSE GPP ), where NSE ET and NSE GPP are the Nash-Sutcliffe efficiency (NSE) for ET and GPP, respectively obtained by comparing the PML-V2 simulations and the flux tower observations (Zhang et al 2019). The performance of PML-V2 simulated GPP and ET was shown in figure S1 and discussed in the section 3.

Observational forcing datasets
In this study, we applied the PML-V2 model to produce historical coupled estimations of global GPP, ET and its components E c , E s and E i over 1982-2014 at a 0.05 • and half-monthly resolution. The climate forcing inputs for the PML-V2 model were land surface climate variables (i.e. precipitation, nearsurface air temperature, pressure and relative humidity, wind speed at 10 m, downward longwave and shortwave radiations) from the GLDAS V2.0 observational 3 hourly climate dataset (Rodell et al 2004, Beaudoing andRodell 2019). The LAI forcing was the global satellite-based half-monthly LAI from the satellite AVHRR GIMMS3g V3 product (Zhu et al 2013). The global satellite-based 8 d surface albedo for shortwave and surface emissivity for longwave radiation from the GLASS V4 product (Liu et al 2013, Cheng andLiang 2014) were also used. The CO 2 forcing for PML-V2 was the global monthly mean atmospheric CO 2 records from NOAA (ftp://aftp. cmdl.noaa.gov/products/trends/co 2 /co 2 _mm_gl.txt).  (Huntzinger et al 2013). It is noted that the PML model does not include nitrogen-cycle process and therefore direct impact of nitrogen deposition was not included. However, the relative contribution (∼9%) from global nitrogen deposition was largely reflected in changes in LAI (Zhu et al 2016), which is considered as a driver in the PML-V2 model. Using the PML-V2 model with observational climate data, atmospheric CO 2 concentration, satellitederived LAI and albedo and emissivity datasets as forcings, we conducted five scenarios of sensitivity simulations (S1-S5, see table 1) by setting all timevarying drivers (S1), or one or more drivers using fixed climate at 1982 (constant-1982) with others time-varying (S2-S4), or all constant-1982 drivers (i.e., S5 as the control run) as inputs over 1982-2014. Therefore, differences between pairs of these sensitivity simulations allow for a robust assessment of model-derived carbon and water fluxes (e.g. GPP, ET, E c , E s , E i , EWUE, and PWUE) sensitivity (including variability and trend) to individual drivers (S1-S2 for albedo and emissivity's contribution, S2-S3 for LAI 's contribution, S3-S4 for CO 2 's contribution, and S4-S5 for climate's contribution) or interactions between two or more drivers (S2-S4 for LAI and CO 2 's contribution, S3-S5 for climate and CO 2 's contribution, S2-S5 for LAI, CO 2 and climate's contribution, S1-S5 for all drivers (i.e., historical change)). For example, estimating the increased atmospheric CO 2 's contributions to historical changes in annual GPP, ET, Ec and corresponding EWUE and PWUE is given by equations (9)-(11): where the subscripts S3 and S4 represent output variables from simulation scenarios S3 and S4 respectively (table 1). The increased CO 2 's impacts on EWUE and PWUE can be further divided into the CO 2 fertilization effect and the CO 2 -induced suppression effect on plant stomata (equations (14) and (15)): ∆PWUE CO2 = ∆PWUE CO2,fertilization We then also estimated the elevated CO 2 -induced fertilization and suppression effects' impacts on EWUE and PWUE by the following equations (16)- (19):

Relative contribution of individual drivers
We calculated the relative contributions of individual drivers (table 1) to the global or regional long-term changes or trends in the studying variables (e.g. annual EWUE and PWUE) by (Zhang et al 2021): where X i is one of the four major model drivers (climate forcing, CO 2 , LAI, surface albedo and emissivity) or one of the CO 2 fertilization and stomatal suppression effects, and Y is one of the studying variables (GPP, ET, E c , E s , E i , EWUE, and PWUE). From equation (2), a negative (positive) value of the relative contribution of X i means that driver X i has a negative (positive) contribution to the overall change or trend in studying variable Y.

Results
We validated the mean annual EWUE simulated by PML-V2 against observation-based annual EWUE over 2002-2014 derived from 95 eddy flux sites from the FLUXNET2015 dataset. Result shows that the linear correlation between observation-based and PML-simulated EWUE values across sites was 0.75 (R 2 = 0.56), with a relative bias of 24% and a root mean squared error (RMSE) of 0.96 gC mm −1 H 2 O, respectively ( figure 1(a)). In addition, the PML-V2simulated half-monthly ET and GPP also showed a good performance with the FLUXNET2015-based dataset, e.g. the NSE was 0.70 for ET and 0.74 for GPP, and the RMSE was 0.7 mm d −1 for ET and 2.05 gC m −2 d −1 for GPP, and bias were less than 5% (figures S1(a) and (b)). Over 1982-2014, the mean ± standard deviation (s.d.) values of historical global total GPP and ET were 150.66 ± 3.69 GtC yr −1 and 66.56 ± 0.66 × 10 3 km 3 yr −1 respectively. The multiyear mean global GPP falls into the full range (119-169 GtC yr −1 ) of multiple process-based models (Anav et al 2015) Figure 1(b) shows that the annual global total GPP had a significant growing trend of 3.33 ± 0.06 GtC yr −1 decade −1 (R 2 = 0.76) over 1982-2014, while the annual global total ET had no significant trend (R 2 = 0.21). As a result, we found that the annual EWUE (=GPP/ET) estimated by the PML-V2 model had a mean ± s.d. of 2.26 ± 0.05 gC mm −1 H 2 O and an increasing trend of 0.040 ± 0.004 gC mm −1 H 2 O decade −1 (R 2 = 0.71) for 1982-2014. We also compared the PML-V2 estimated EWUE with those from the CMIP6 models (Eyring et al 2016) and the MsT-MIP models (Huntzinger et al 2013). The CMIP6 models provide the state-of-the-art simulations of historical changes in annual GPP and ET forced by observational external greenhouse gases (e.g. CO 2 , CH 4 , NO x etc) emissions and human-induced landuse change. The terrestrial ecosystem models of MsTMIP were run offline and driven by observational surface climate, atmospheric CO 2 , land-use change and nitrogen deposition. We found that the global EWUE from an ensemble of 12 CMIP6 models increased continuously from 1982 to 2014, with a linear trend of 0.036 ± 0.015 gC mm −1 H 2 O decade −1 (R 2 = 0.79) which is very close to the PML-V2 simulated EWUE trend ( figure 1(c)). However, the simulated global EWUE over 1982-2010 from an ensemble of 13 MsTMIP models showed a much larger trend of 0.058 ± 0.032 gC mm −1 H 2 O decade −1 (R 2 = 0.75) and a larger uncertainty range across models compared to the CMIP6 models ( figure 1(c)).
With the sensitivity modeling experiments using PML-V2, we were able to isolate the individual contributions of environmental drivers, i.e. climate, CO 2 , LAI, albedo and emissivity, and their interactions (see section 2). We found that at the global scale, the summed effect of all individual drivers' contribution can explain >99% of changes in global annual GPP, ET, E c , E s and E i , indicating that the interactions between drivers (e.g. contribution of LAI and CO 2 minus sum of contribution of LAI and contribution of CO 2 ) are negligible in this study ( figure S2). Noticeably, the increased CO 2 had resulted in a substantial increase in global annual GPP and a decrease in global annual ET and E c over 1982-2014 (figures S2(a)-(c)). The direct contributions from increased CO 2 -induced fertilization effect on plant photosynthesis (43.94%) and increased LAI-induced growth of plant production (54.74%) together explain >98% of the increasing trend in global annual GPP over 1982-2014 (figures S2(a) and (f)). However, the global annual total ET with no significant trend is result from the positive contributions of LAI (23.30%) and climate change (31.33%) and the negative contributions of increased CO 2 (−44.39%) and surface albedo and emissivity (−0.98%). The negative contribution of increased CO 2 to global annual ET almost all comes from the CO 2 -induced suppression effect (−49.45%) on canopy/plant transpiration (E c ) through reducing stomatal conductance, counteracting almost all of the positive contributions from LAI (23.98%) and climate change (26.48%) to E c (figure S2(f)). As the other two components of ET, the E s and E i both have larger contributions to ET increase than E c ( figure  S2(f)). The larger contribution of E i (68.52%) than E s (28.03%) is mainly because that the larger vegetation cover from increased LAI leads to smaller E s and larger E i .
We further isolated the individual contributions of environmental drivers to the trends in annual EWUE and PWUE, and found that the interactive contributions between drivers are negligible (figure 2). We found that the increasing trend in global PWUE (=GPP/E c ) for 1982-2014 was 0.090 ± 0.005 gC mm −1 H 2 O decade −1 (R 2 = 0.90), which is about two times larger than global EWUE trend ( figure 2(c)). For both the global annual EWUE and PWUE, increased CO 2 is the largest contributor as it induced 51.36% and 71.93% of the increasing trends in EWUE and PWUE, respectively (figure 2). The CO 2 -driven increase in EWUE is the combined contribution of the CO 2 fertilization effect-induced increase in GPP and the CO 2 stomatal suppression effect-induced reduction in ET, which respectively contributed 23.01% and 28.35% of the global EWUE trend. For the CO 2 -driven increase in PWUE, the CO 2 fertilization effect and stomatal suppression effect respectively contributed 22.67% and 49.26% of the global PWUE trend. In contrast, climate change as the second largest driver contributed −26.35% and −24.19% of the EWUE and PWUE trends respectively. Additionally, the LAI alone induced 20% of the EWUE increase but only a small impact (3.84%) on the PWUE increase (figure 2(c)), indicating that the LAI-induced increase in E i (66.76%) dominated the LAI-induced increase in EWUE ( figure S2(f)).
At regional scales, significant positive trends (⩾0.02 gC mm −1 H 2 O decade −1 ) in annual EWUE for 1982-2014 covered nearly half (47.42%) of global land area, including the western USA, western Europe, central and southern Africa, India, central South America, central North China, and eastern Australia ( figure 3(a)). These regions are dominated by semi-arid grassland or grasslands, with the mean annual EWUE being mostly smaller than that in forest-dominated regions (see figure S3(a)). In these regions, trends in annual PWUE for 1982-2014 were consistently positive (⩾0.05 gC mm −1 H 2 O decade −1 over 45.78% of global land area, figure 3(b)). Whereas the significant negative trends in EWUE and PWUE respectively covered only 11.90% and 7.44% of global land area mainly over the >50 • N latitudes boreal forest or tundra-dominated regions (figure 3) where the mean annual PWUE are mostly larger than 6 gC mm −1 H 2 O (see figure S3(b)).
We further investigated the dominant drivers which have the largest relative contributions to regional trends in EWUE and PWUE. Most regions with positive trends in EWUE and PWUE were dominated by the increased CO 2 and LAI, while regions with negative trends in EWUE and PWUE were dominated by climate change (figure 4). Except for the unvegetated regions (29.23% of global land area), the CO 2 , LAI and climate change are the largest three contributors of regional trends in EWUE, covering 29.37%, 21.45% and 19.26% respectively of global land area ( figure 4(a)). For trends in EWUE, most of tropical forests were dominated by the increased CO 2 , while the semi-arid regions and boreal forests were dominated by the LAI or climate change ( figure 4(a)). For trends in PWUE, the CO 2 is the largest contributor, accounting for 52.50% of global land area including most tropical forests and semi-arid regions, while the climate change covering 15.65% of global land area most of which are over eastern Europe and eastern Asia ( figure 4(b)).
Due to the dominant contributions of the increased CO 2 to the global/regional trends in EWUE and PWUE (figures 2-4), we further isolated the relative contributions of increased CO 2 fertilization effect and stomatal suppression effect-induced changes in EWUE and PWUE to the regional trends in EWUE andPWUE over 1982-2014 (figure 5). The CO 2 -dominanted positive trends in EWUE and PWUE over most tropical forests (figure 4) were controlled by the increased CO 2 stomatal suppression effect which contributes 30% to more than 50% of the regional trends in EWUE and PWUE (figure 5). The CO 2 stomatal suppression effect also contributed about 10%-30% of the overall trends in EWUE and PWUE over temperate and boreal forests where the CO 2 , LAI or climate change were the dominant drivers (figure 5). Whereas most of global semi-arid regions with positive trends in EWUE and PWUE were controlled by the increased CO 2 fertilization effect, which particularly contributes more than 40% of the regional trends in EWUE (figure 5).

Discussion
Our results support previous studies (Cheng et al 2017) that the historical long-term growth of global GPP since 1980s was mainly driven by the increased ecosystem-scale water-use efficiency (EWUE). We estimated that the increased EWUE was largely contributed by increased CO 2 (51.36%) and LAI (20%) but was counteracted by climate change (−26.35%). The PML-V2 model estimated global EWUE trend (0.040 ± 0.004 gC mm −1 H 2 O decade −1 ) is also smaller than the previous estimate (0.064 ± 0.017 gC mm −1 H 2 O decade −1 ) using four process-based ecosystem models driven by observational climate (Huang et al 2015), and is much smaller than the previous estimate (0.137 ± 0.043 gC mm −1 H 2 O decade −1 ) using an analytical leaf-to-ecosystem scaling WUE model forced by observational CO 2 , climate and satellitebased LAI since 1980s (Cheng et al 2017). Most previous modeling studies agreed that the global EWUE increases over last three decades were mainly driven by rising CO 2 (Huang et al 2015, Cheng et al 2017, El Masri et al 2019. The WUE model used by Cheng et al (2017) showed a relative contribution of 47.83% from increased CO 2 (that was estimated with equation (12)  . This large discrepancy of increased CO 2 contributions between ecosystem models and our PML model-based or the WUE model-based estimates is because that the LAI was not considered as an external driver but being a predicted variable in ecosystem models driven by CO 2 and climate change. This study divided the CO 2 impact on EWUE into the direct effects (i.e. the fertilization effect and stomatal suppression effect), and the indirect effect that was reflected by LAI changes due to increased CO 2 (with 70% contribution), climate change and nitrogen deposition (Zhu et al 2016).
Our new findings suggested that the increased CO 2 dominant impact (51%) on global EWUE increase includes contributions from fertilization effect on photosynthesis (23%) and suppression effect on stomatal conductance (28%). For global PWUE increase, the CO 2 's contribution increased to 71.93%, with 22.67% from fertilization effect and 49.26% from stomatal suppression effect. This result, unexpectedly, indicates that stomatal suppression effect on plant transpiration almost dominated the long-term PWUE increase over last four decades. The CO 2 dominant contributions were largely reflected by very strong plant physiological responses to arising CO 2 , which expresses the relationship between stomatal conductance and leaflevel water-use efficiency by quantifying the intrinsic water-use efficiency (iWUE) from the diffusion equation [iWUE = A/g = (C a − C i )/1.6], linking photosynthesis rate (A) to stomatal conductance (g) (Ehleringer et al 1993). The CO 2 fertilization effect induces an increase in A, and the CO 2 suppression effect induces a reduction in g, which maintain constant intercellular CO 2 concentration (C i ) in response to increased atmospheric CO 2 , C a (Keenan et al 2013).
Our new findings indicated that the CO 2dominated positive trends in PWUE over tropical and European forests were dominantly controlled by the increased CO 2 stomatal suppression effect. Whereas most of global semi-arid regions with positive trends in PWUE, were dominated by the increased CO 2 fertilization effect. The increasing trends in PWUE over tropical and European forests are compatible with observation-based iWUE derived from eddy covariance data (Keenan et al 2013) and tree-ring records (van der Sleen et al 2014, Frank et al 2015, Dekker et al 2016. The CO 2 stomatal suppression effect-induced dominant impact on tropical forests can largely offset the increase in ET due to climate warming, leading to a little change of tropical ET, which may further change regional-to-global water cycles, such as increase river runoff (Gedney et al 2006, Piao et al 2007. While the strong control of CO 2 fertilization effect over semi-arid regions may cause larger LAI and more ET, which then induce soil drying and enhance severity of droughts. In summary, this study highlights the important role of increased CO 2 -induced fertilization effect and stomatal conductance reduction in enhancing worldwide plant-to-ecosystem WUE and the potential of regulating global terrestrial carbon and water cycles.
These findings could help policymakers for water use planning and ecosystem management. For example, in near future, vegetation expansion over semi-arid regions may induce risks of regional water resource safety, as CO 2 fertilization-induced WUE increase does not reduce water loss. While for forest regions, the increased WUE could be beneficial to afforestation over tropical or temperature wet regions, as stomatal suppression-induced WUE increase for forests leads to lower loss of water from ecosystems.

Data availability statement
The data that support the findings of this study are openly available at the following URL/DOI: https://doi.org/10.6084/m9.figshare.14184740. Data will be available from 18 November 2021.