Quantifying Soil Moisture Impacts on Water Use Efﬁciency in Terrestrial Ecosystems of China

: Soil moisture ( SM ) signiﬁcantly affects the exchange of land surface energy and the stability of terrestrial ecosystems. Although some conclusions have been drawn about the effects of SM on the ecosystem water use efﬁciency ( WUE ), the inﬂuence mechanism and the quantitative assessment framework of SM on WUE are still unclear. This study provides an analysis framework for the feedback relationship between SM and WUE based on the dependence of the evaporation fraction on SM and output datasets from remote sensing and the Global Land Data Assimilation System. The results show that the range of WUE of terrestrial ecosystems of China was 0.02–19.26 g C/kg H 2 O in the growing season with an average value of 1.05 g C/kg H 2 O. They also show a downward trend in 43.99% of the total area. In the evapotranspiration ( ET ) pathway, SM negatively affected WUE , and the sensitivity coefﬁcient ranged from − 18.49 to − 0.04. In the net primary production ( NPP ) pathway, the sensitivity coefﬁcient ranged from − 68.66 to 43.19. Under the dual effects of the ET and NPP pathways, the inﬂuence of SM on WUE was negative in 84.62% of the area. Variation in SM led to signiﬁcant WUE variability. Generally, the percentage change in WUE ( ∆ WUE ) ranged from 0% to 190.86%, with an average value of 28.02%. The maximum ∆ WUE ranged from 0% to 758.78%, with an average value of 109.29%. The WUE of forest ecosystems showed strong resistance to SM variation, whereas that of non-forest vegetation was more sensitive to SM variation. This analytical framework provides a new perspective on the feedback relationship between WUE and SM in This study mainly discussed the variation in and SM in different climate zones of the terrestrial ecosystem. The terrestrial ecosystem in China has eight climatic zones, namely, the marginal tropical zone (A), south subtropics zone (B), middle subtropics zone (C), north subtropics zone (D), warm temperate zone (E), middle temperate zone (F), cold temperate zone (G), temperate zone in Qinghai–Tibet Plateau (H), and sub-frigid zone in Qinghai–Tibet Plateau (I) (Figure 1). As the terrestrial ecosystems of China span multiple climatic zones, its climatic growth seasons also have signiﬁcant spatial differences, ranging from the shortest two-month growth season to the longest twelve-month growth season (Figure 1).


Introduction
The water use efficiency (WUE) of an ecosystem is the ratio of carbon sequestration to water consumption, which couples the terrestrial carbon and water cycles and closely links photosynthesis and evapotranspiration (ET) processes in the ecosystem [1,2]. Although there are many ways to express the WUE of an ecosystem, the ratio of the net primary production (NPP) to ET is widely used [3]. The variation of WUE controlled by both biotic and abiotic factors [4,5] and recent research showed that a vapor pressure deficit (VPD) and canopy conductance are two dominant factors of WUE in response to drought [6]. On the global scale, WUE is positively related with drought in the majority of regions [7]. In China, most forest ecosystems exhibit strong drought resistance by improving their WUE [8]. Under the influence of the negative feedback relationship between ET and air temperature [9,10], the variation in WUE is usually positively correlated with air

Data
Monthly latent heat flux (LE), sensible heat flux (H)and ground heat flux (G), ET, surface net radiation (Rn), and SM (four layers: 0-10, 10-40, 40-100, 100-200 cm) datasets were obtained from the NASA Global Land Data Assimilation System (GLDAS) Noah model (GLDAS_NOAH025_M_2.0). We used the dataset from 1985 to 2014 with a spatial resolution of 0.25° × 0.25° (downloaded from https://disc.gsfc.nasa.gov/datasets, 10/10/2020). In this study, SM in the 0-100 cm soil layer was regarded as the SM of root zone based on previous research [35][36]. The monthly NPP dataset (spatial resolution of 1 km×1 km) of China's terrestrial ecosystems north of 18°N from 1985 to 2014 was obtained from the Global Change Research Data Publishing & Repository of China. The NPP data were calculated using the Carnegie-Ames-Stanford Approach (CASA) model based on the monthly meteorological data of China's land from 1985 to 2015, national soil texture data, and land cover and vegetation index data from the Moderate Resolution Imag-ingSpectroradiometer (MODIS) and the Advanced Very High Resolution Radiometer (AVHRR) remote sensing images. CASA model is a process-based remote sensing model for estimating NPP of terrestrial ecosystem, which is driven by grid datasets of climate, radiation, soil, and remote sensing vegetation index [37]. A comparative analysis with the measured and simulated data from previous studies [38][39][40][41][42] showed that the NPP dataset used in this study hadgood accuracy [43].

Data
Monthly latent heat flux (LE), sensible heat flux (H) and ground heat flux (G), ET, surface net radiation (R n ), and SM (four layers: 0-10, 10-40, 40-100, 100-200 cm) datasets were obtained from the NASA Global Land Data Assimilation System (GLDAS) Noah model (GLDAS_NOAH025_M_2.0). We used the dataset from 1985 to 2014 with a spatial resolution of 0.25 • × 0.25 • (downloaded from https://disc.gsfc.nasa.gov/datasets, accessed on 10 October 2020). In this study, SM in the 0-100 cm soil layer was regarded as the SM of root zone based on previous research [35,36]. The monthly NPP dataset (spatial resolution of 1 km × 1 km) of China's terrestrial ecosystems north of 18 • N from 1985 to 2014 was obtained from the Global Change Research Data Publishing & Repository of China. The NPP data were calculated using the Carnegie-Ames-Stanford Approach (CASA) model based on the monthly meteorological data of China's land from 1985 to 2015, national soil texture data, and land cover and vegetation index data from the Moderate Resolution Imaging Spectroradiometer (MODIS) and the Advanced Very High Resolution Radiometer (AVHRR) remote sensing images. CASA model is a process-based remote sensing model for estimating NPP of terrestrial ecosystem, which is driven by grid datasets of climate, radiation, soil, and remote sensing vegetation index [37]. A comparative analysis with the measured and simulated data from previous studies [38][39][40][41][42] showed that the NPP dataset used in this study hadgood accuracy [43].
To verify the NPP and ET (GLDAS) data, this study collected the monthly observation data from eight flux towers in China from 2003 to 2010, which were provided by ChinaFLUX (http://www.cnern.org.cn/, accessed on 15 September 2021). In addition, this study also collected the monthly observation flux data from the other four flux towers from the stations, and the collected data length of these four stations was 2004-2006, 2011-2014, 2008-2011, and 2013-2014, respectively. This study also used SM data based on microwave remote sensing in China from 2002 to 2011 (provided by the National Tibetan Plateau Data Center (TPDC) of China; http://data.tpdc.ac.cn, accessed on 2 August 2021) to verify the SM data from GLDAS. The spatial and temporal resolutions of the SM dataset (from the TPDC) were 0.25 • and monthly, respectively. The dataset was calculated based on the high spatial and temporal resolution surface meteorological dataset and an improved land surface assimilation system to drive the Simple Biosphere model (SiB2) and assimilate the brightness temperature observed by the Advanced Microwave Scanning Radiometer Earth Observing System sensor (AMSR-E) satellite [44,45]. The verifications showed that the RMSE of this SM dataset (from the TPDC) was approximately 5% volumetric soil water content [45,46].
The soil texture dataset was obtained from the TPDC of China (http://data.tpdc.ac.cn, accessed on 13 May 2020). The dataset was calculated based on the 1:1,000,000 scale soil map and 8595 soil profiles of China's second national soil survey and the regional land and climate simulation standard of the U.S. Department of Agriculture (soil texture classification standard is shown in Table 1). This study also used the data of soil hydraulic properties in China, which were provided by the National Cryosphere Desert Data Center of China (http://www.crensed.ac.cn/portal/, accessed on 21 August 2020). This study extracted the growing season by using the base temperature threshold of 10 • C (≥10 • C) [47] based on the multi-year monthly average air temperature. The required grid data of monthly average temperature in China (1985China ( -2014 were downloaded from the "National Tibetan Plateau Data Center" of China (http://data.tpdc.ac.cn, accessed on 13 May 2020). However, in the Qinghai Tibet Plateau, the maximum monthly average temperature in most areas is usually lower than 10 • C. Therefore, according to the local actual situation and previous studies [48,49], the growth season was defined as June to September.

Mathematical Expression between SM and EF
Previous studies have proposed an effective analysis method for the functional relationship between the evaporation fraction (EF) and SM [50,51]. This method expressed the dependence of EF on SM (θ) as follows: where θ is volumetric soil moisture content (m 3 /m 3 ), θ c is the soil moisture at the critical point, and θ r is the soil wilting coefficient. EF is the evaporative fraction, which is calculated as follows: where LE is the latent heat flux and R n is the surface net radiation. When θ < θ r , EF is not always zero owing to the presence of hygroscopic SM [52]. Therefore, EF min is defined as the EF value corresponding to a cumulative frequency percentage of 5%. In addition, when θ > θ c , energy rather than SM became the dominant limiting factor for ET. Thus, the study similarly defined EF max as the EF value corresponding to 95% of the cumulative frequency percentage. Because this study only focused on the functional relationship between SM and EF in the range of EF min ≤ EF ≤ EF max , Equation (1) could be expressed as follows: This study estimated the raster map of θ r using soil hydraulic parameters and soil texture datasets [53]. Although the fitting between EF and θ could obtain the parameter θ c , it could result in great uncertainty. Thus, this study used the soil field capacity (θ f ) to replace θ c and introduced a dimensionless parameter A into the calculation scheme for correction. Finally, Equation (3) was deformed as follows: The gridded θ f value was also calculated using soil hydraulic parameters and soil texture datasets [53]. First, we used the monthly SM and EF data from 1985 to 2000 to fit the multi-year average monthly (January to December) parameter A. Then, the monthly SM from 2001 to 2014 was estimated using the EF data during the same period and the fitted parameter A. The GLDAS SM data finally verified the estimated soil moisture data during 2001-2014.

The Response of WUE to SM Variation
Generally, WUE is calculated as follows: Based on Equation (2), ET could be expressed as follows: where L is the latent heat of vaporization, which is the ratio of LE to ET. Thus, the elasticity of WUE to the variation of SM could be divided into two contributions based on Equations (4) to (6), as follows: The first term on the right side of the equation describes the impact of SM variation on the ET process and the further influence on the WUE. This term was calculated based on Equations (3), (5) and (6). The second term on the right of the equation is the sensitivity of WUE to SM variation via the NPP pathway, which was estimated based on Equation (5) and the following statistical method: where ε is the sensitivity coefficient; NPP and θ are the long-term average values of the NPP and SM in growing season of a pixel; NPP i and θ i are the growing season values of NPP and SM in each year, respectively. Thus, the following formula could be used to evaluate the average impact of SM variation on WUE in the past 30 years: where ∆θ is the percentage change in SM in the growing season from 1985 to 2014. Thus, ∆WUE was the average influence (change percentage) of SM variation on WUE. Similarly, the impact of the maximum variation of SM on WUE was calculated as follows: where ∆θ max is the maximum rate of change in SM during the growing season from 1985 to 2014, which was calculated as the difference between the maximum and minimum values of SM during the period. ∆WUE max is the maximum possible impact of SM variation on WUE during the growing season in past 30 years.

Data Verification
The analysis of this study was based on the output datasets of the climate model and remote sensing model, including the ET and SM datasets from GLDAS, and the NPP dataset from the Global Change Research Data Publishing & Repository of China. The verification results showed that the monthly ET of the GLDAS dataset and the observed ET data in 12 flux towers in China had a good consistence (Figure 2b). The determination coefficient (R 2 ), root mean square error (RMSE), mean absolute percent error (MAPE), and Nash-Sutcliffe efficiency coefficient (NSE) between the simulated and measured ET values were 0.79, 14.78 mm, 31.32%, and 0.77, respectively. In addition, the verification results of NPP data showed that, although the estimated monthly NPP data used in this study underestimated the actual NPP, the R 2 , RMSE, MAPE, and NSE between the estimations and observations were 0.97, 69.02 g C·m −2 ·mo −1 , 35.09%, and 0.42, respectively, which indicated that the NPP dataset used in this study was also consistent with the observations (Figure 2a). This study compared the SM data based on microwave remote sensing data assimilation in China (2002-2011) with GLDAS SM data (Figure 2c), owing to the lack of measured SM data. The results showed a good consistency, and the R 2 , RMSE, MAPE, and NSE values were 0.96, 0.05 m 3 ·m −3 , 18.25%, and 0.59, respectively. Furthermore, The monthly SM from 2001 to 2014 was estimated using Equation (4) and then verified using SM data of GLDAS. The overall trend of change in estimated SM was close to the GLDAS SM data, and the determination coefficient (R 2 ) was 0.97 ( Figure 2d). Meanwhile, the RMSE, MAPE, and NSE of this simulation were 0.01 m 3 ·m −3 , 18%, and 0.92, respectively. The fitting error was primarily due to the overestimation of SM, especially when the SM was low. This study may have overestimated SM in some arid regions or during dry seasons. Although there were still some simulation errors, the datasets used in this study were generally reliable and could ensure the accuracy of the results. Remote Sens. 2021, 13, x FOR PEER REVIEW 7 of 18

Spatiotemporal Variation of WUE
In terms of the growing season, the range of WUE in all terrestrial ecosystems was 0.02-19.26 g C/Kg H2O, with an average value of 1.05 g C/Kg H2O and a standard deviation of 1.45 g C/Kg H2O. Furthermore, the high-value areas of WUE were mainly located in the arid regions of Northwest China (Figure 3a). The spatial variation in WUE in the different climatic zone was also distinct ( Figure 3c). The sub-frigid zone in the Qinghai-Tibet Plateau(PSF) and the temperate zone in the Qinghai-Tibet Plateau (PTMP) had the highest WUE, with maximum values of 1.80 g C/Kg H2O and 1.03 g C/Kg H2O, respectively. The north subtropics (NST) had the lowest WUE with a value of 0.48 g C/Kg H2O. The WUE in the warm temperate zone (WTMP), middle temperate zone (MTMP), cold temperate zone (CTMP), south subtropics zone (SST),marginal tropical(MT) zone, and middle subtropics zone (MST) was0.98, 0.87, 0.64, 0.60, 0.55, and 0.51 g C/Kg H2O, respectively. In general, the change in WUE was correlated with the SM (coefficient of determination: R =0.56), and a higher SM usually led to a lower WUE (Figure 3d).

Spatiotemporal Variation of WUE
In terms of the growing season, the range of WUE in all terrestrial ecosystems was 0.02-19.  During the growing season, WUE showed a downward trend in 43.99% of the total area. Among these areas, WUE showed extremely significant (Sig = 0.99), significant (Sig = 0.95), and insignificant (Sig = 0.90) downward trends in 35.03%, 6.11%, and 2.58% of the total terrestrial ecosystem area, respectively. In contrast, WUE also showed an increasing trend in 34.91% of the total area of the terrestrial ecosystem. Among these areas, the percentage with extremely significant, significant, and insignificant increases in WUE was 28.07%, 4.94% and 1.91%, respectively. Thus, in the past 30 years, the WUE of the terrestrial ecosystem has generally shown a downward trend, especially in the arid areas of Northwest China, the Mongolian Plateau, and the southern Qinghai-Tibet Plateau. WUE also showed a significant downward trend in some areas of the monsoon region of China (Figure 3b).

Elasticity of WUE to SM
The variation in SM can, simultaneously, directly affect the NPP and ET, which, ultimately, leads to the changes in WUE. On the ET pathway, the sensitivity coefficient of the WUE to SM ( ) at the growing season scale ranged from −18.49 to −0.04, with an average of −3.02 (Figure 4a). Generally, the negative impact of SM on WUE via the ET pathway in the arid region of northwestern China was the most significant, followed by that in the southeast monsoon region. The impact in the central China region was rela- During the growing season, WUE showed a downward trend in 43.99% of the total area. Among these areas, WUE showed extremely significant (Sig = 0.99), significant (Sig = 0.95), and insignificant (Sig = 0.90) downward trends in 35.03%, 6.11%, and 2.58% of the total terrestrial ecosystem area, respectively. In contrast, WUE also showed an increasing trend in 34.91% of the total area of the terrestrial ecosystem. Among these areas, the percentage with extremely significant, significant, and insignificant increases in WUE was 28.07%, 4.94% and 1.91%, respectively. Thus, in the past 30 years, the WUE of the terrestrial ecosystem has generally shown a downward trend, especially in the arid areas of Northwest China, the Mongolian Plateau, and the southern Qinghai-Tibet Plateau. WUE also showed a significant downward trend in some areas of the monsoon region of China (Figure 3b).

Elasticity of WUE to SM
The variation in SM can, simultaneously, directly affect the NPP and ET, which, ultimately, leads to the changes in WUE. On the ET pathway, the sensitivity coefficient of the WUE to SM ( ∂WUE ∂ET ∂ET ∂θ ) at the growing season scale ranged from −18.49 to −0.04, with an average of −3.02 (Figure 4a). Generally, the negative impact of SM on WUE via the ET pathway in the arid region of northwestern China was the most significant, followed by that in the southeast monsoon region. The impact in the central China region was relatively   In the NPP pathway, the sensitivity coefficient of WUE to SM ( ∂WUE ∂NPP ∂NPP ∂θ ) had a more apparent spatial variation, and the value ranged from −68.66 to 43.19 with an average of 0.13 (Figure 4c). Thus, SM via the NPP pathway positively affected WUE in 55.71% of the area, whereas it showed a negative influence in the remaining 44.29% of the area. In addition, SM via the NPP pathway often had a more significant impact (larger absolute value) on WUE in the southeast monsoon region and northwestern arid region of China. In different climatic zones (Figure 4d), the impact of SM on WUE via the NPP pathway showed a positive effect in MTMP, WTMP, and PTMP, with elasticity coefficients of 0.74, 0.50, and 0.33, respectively. In contrast, the impacts of SM on WUE via the NPP pathway were negative in the other climatic zones, and the elasticity coefficients ranged from −1.39 to −0.12. Generally, the area percentage of the region in which SM had a positive effect on WUE via the NPP pathway ranged from 28.10% to 63.44% in different climatic zones.
Under the dual effects of ET and NPP pathways, the final SM significantly impacted WUE. The elasticity coefficient of WUE to SM ranged from −74.01 to 39.81, with an average value of −2.87 (Figure 4e). Because the impact of SM on WUE via the ET pathway was completely negative, the comprehensive influence of SM on WUE also showed a more remarkable negative effect. In the terrestrial ecosystem, SM negatively and positively impacted WUE in 84.62% and 15.38% of the total area, respectively. For each climatic zone, the average value of the effect of SM on WUE was also negative, and the negative effect was the most significant in the NST, which had an average elasticity coefficient of −4.27. The average value was the lowest in PTMP, which had an average elasticity coefficient of −1.93. For the other climatic zones, the elasticity coefficient of WUE to SM ranged from −4.03 to −2.01. Generally, the NPP pathway mainly determined the trend of change in WUE in different climatic zones, which was because the correlation coefficient between ∂WUE ∂θ and ∂WUE ∂NPP ∂NPP ∂θ (R = 0.84) was more significant than that between ∂WUE ∂θ and ∂WUE ∂ET ∂ET ∂θ (R = 0.69). In contrast, the ET pathway mostly dominated the influence level (absolute value of the elasticity coefficient of WUE to SM) of SM on the WUE (Figure 4b,d,f).

WUE Variability Caused by SM
In this study, the absolute value of the percentage change in the WUE (∆WUE) was used as an index to analyze the WUE variability. Figure 5a shows that the variation of SM led to significant WUE variability, especially in northeast and northwest China. At the growing season scale, the average ∆WUE due to the variation in SM ranged from 0% to 190.86%, with an average value of 28.02% (Figure 5a). For different climatic zones (Figure 5c), SM in the CTMP had the most significant effect on WUE, causing an average variation in WUE of approximately 53.90%. In contrast, SM in the PTMP had the least impact on WUE, causing only, approximately, a 11.96% variability in WUE. For the other climatic zones, ∆WUE caused by SM variation in MTMP, WTMP, NST, MT, SST, MST, and PSF were 38.93%, 27.75%, 23.04%, 18.39%, 14.78%, 12.84%, and 12.30%, respectively. The spatial trend of the maximum WUE variability (∆WUE max ) was similar to the variability of the ∆WUE (Figure 5b). It also had the greatest change in northeast and northwest China. Furthermore, the ∆WUE max ranged from 0% to 758.78%, with an average value of 109.29%. These results indicated that the maximum change in SM over the past 30 years doubled the variability in WUE. In different climatic zones, the trend of change in ∆WUE max was the same as that for ∆WUE, but the value of ∆WUE max was larger than that of ∆WUE. ∆WUE max was generally 3.7-4.4 times greater than ∆WUE (Figure 5c). The feedback relationship between SM and WUE was affected by the climatic factors and vegetation cover. The ratio of vegetation cover between forest and non-forest was also an essential factor affecting the feedback relationship between WUE and SM (Figure 5d). The lower the forest cover, the more significant the change in WUE caused by changes in SM. As the forest cover increased, the variability of WUE caused by changes in SM rapidly decreased. When the cover ratio of forest to non-forest vegetation was higher than one, the variability of the WUE caused by changes in the SM tended to stabilize. The WUE of the forest ecosystem had strong resistance to variation in soil water. In contrast, the WUE of non-forest vegetation was more sensitive to the changes in SM.

The Reliabilityof Evaluation Results
In most cases, SM dominates ET, especially in water-limited regions [54][55]. To explore the coupling relationship between SM and ET, researchers proposed a simple parameterization schematic, in which the SM and ET of a forest presented a linear relationship [56]. Recently, studies have proposed a conceptual framework for the dependence of the evaporative fraction (EF) on SM. It is suitable to analyze the soil moisture control on energy partitioning [50][51]. Based on this coupling relationship between SM and EF, this study established an evaluation framework for the feedback relationship between the WUE and SM. The uncertainty of this evaluating framework mainly origi- The feedback relationship between SM and WUE was affected by the climatic factors and vegetation cover. The ratio of vegetation cover between forest and non-forest was also an essential factor affecting the feedback relationship between WUE and SM (Figure 5d). The lower the forest cover, the more significant the change in WUE caused by changes in SM. As the forest cover increased, the variability of WUE caused by changes in SM rapidly decreased. When the cover ratio of forest to non-forest vegetation was higher than one, the variability of the WUE caused by changes in the SM tended to stabilize. The WUE of the forest ecosystem had strong resistance to variation in soil water. In contrast, the WUE of non-forest vegetation was more sensitive to the changes in SM.

The Reliability of Evaluation Results
In most cases, SM dominates ET, especially in water-limited regions [54,55]. To explore the coupling relationship between SM and ET, researchers proposed a simple parameterization schematic, in which the SM and ET of a forest presented a linear relationship [56]. Recently, studies have proposed a conceptual framework for the dependence of the evaporative fraction (EF) on SM. It is suitable to analyze the soil moisture control on energy partitioning [50,51]. Based on this coupling relationship between SM and EF, this study established an evaluation framework for the feedback relationship between the WUE and SM. The uncertainty of this evaluating framework mainly originated from two aspects: (1) the reliability of the datasets used and (2) the uncertainty of the parameterized expression of EF and SM.
Previous studies have shown that evaporation, LE, H, and SM data in the GLDAS dataset have good applicability in China [57,58]. Our verification results (Figure 2) also showed that the datasets used in this study were, overall, reliable, despite some overestimation or underestimation. The analysis method for the coupling relationship between EF and SM has already been applied [51]; however, this framework needs to fit two parameters (namely, θ c and θ r ), which may cause uncertainty. Given this, the study directly used two determined parameters, namely, θ f , and θ r , provided by the dataset of soil hydraulic parameters in China (http://data.tpdc.ac.cn, accessed on 13 May 2020). This study also introduced a dimensionless parameter A. Thus, this study only fit one parameter (A), which could avoid the uncertainty caused by the need to fit two parameters in the original equation and eliminate the error caused by possible overfitting. This study used the monthly SM and EF data from 1985 to 2000 to fit the multi-year average monthly parameter A. The monthly SM from 2001 to 2014 was estimated using Equation (4) and then verified using SM data of GLDAS (Figure 2d). The overall trend of change in estimated SM was close to the GLDAS SM data, and the R 2 , RMSE, MAPE, and NSE of this simulation were 0.99, 0.01 m 3 ·m −3 , 18%, and 0.98, respectively. Although there were still some simulation errors, the estimation was generally reliable. Furthermore, the trend of change in the estimated value was consistent with the GLDAS SM data; thus, the functional relationship between SM and EF was credible overall.

The Coupling Relationship between Soil Moisture and WUE
Studies have shown that WUE is sensitive to changing environments [59]. Generally, many factors, such as the CO 2 concentration, nitrogen deposition, climatic factors, vegetation, human drivers, and SM, significantly influence WUE [6,[60][61][62]. SM is one of the critical factors affecting WUE [21,63]. However, most studies have not thoroughly explored the mechanisms and accurate quantification of the effects of SM on WUE. Our results showed that SM had a completely negative influence on WUE via the ET pathway. However, SM also directly influenced WUE via the NPP pathway, which included negative (44.29% of the total areas) and positive effects (55.71% of the total areas); the overall positive effect was dominant. The areas in which SM had a negative impact on WUE via NPP were mainly distributed in the extremely arid desert areas of Northwest China, the cold temperate areas in Northeast China, and humid areas in Southeast China, which might have been caused by the negative impact of SM on NPP [64].
Under the dual effects of SM, there was a clear negative feedback relationship between WUE and SM, which was generally consistent with existing research [7,31]. The increase in SM often led to a decrease in WUE in most areas. The strength of the feedback relation between SM and WUE may have partly depended on the hydrothermal conditions of the local environment. However, the most significant feature of the spatial variability of the WUE in this study was that WUE showed a downward trend in the arid area of northwestern China. In contrast, there was an upward trend in the humid eastern region of China. Therefore, China's terrestrial ecosystem WUE may be more dominated by ET than by NPP [21]. Because the NPP of tropical and subtropical vegetational zones is higher, it is generally approximately 8-10 times that of an arid desert and alpine vegetation areas in northwest China [65]. However, the WUE in the tropical and subtropical vegetation zones was smaller than that of the desert and alpine vegetation areas, which may have been due to excessive ET reducing the WUE. Furthermore, studies have shown that whether the effect of SM on WUE is positive or negative depends primarily on whether the ET or NPP process dominates the WUE process [21]. This study used the absolute value of ∂WUE ∂ET ∂ET ∂θ and ∂WUE ∂NPP ∂NPP ∂θ to determine which pathway of ET or NPP dominated the feedback relationship between SM and WUE. The results showed that the ET pathway dominated the feedback relationship between SM and WUE in 72% of the total area. The region in which the NPP pathway determined the feedback relationship only accounted for approximately 28% of the total area, especially in the Qinghai-Tibet Plateau. Therefore, terrestrial ecosystems can partially offset the influences of soil dryness by improving WUE.
In addition, this adjustment should be more remarkable in areas where ET controls the WUE process, such as in the arid regions of the northwest.

Plausible Changing Trend and Adjustment of WUE
Under the background of climate change, especially with the increase in temperature and drought, the trend of change in WUE is uncertain. Studies have shown that a temperature rise would increase the NPP of global terrestrial ecosystems [66,67] However, the continuous temperature rise may eventually lead to a significant decrease in NPP [51], mainly owing to the limitation of the optimum temperature for plant growth [68,69]. Thus, WUE may show two completely different changes in future scenarios of a continuous temperature rise. First, the temperature rise in the early stage could result in a general increase in NPP. Meanwhile, the temperature rise in arid areas (water-limited areas) might not contribute much to ET, but may cause a linear increase in humid areas [65,66]. Therefore, under the early warming scenario, WUE in the arid region of northwestern China could tend to increase significantly. In the humid area, WUE may only slightly change. In contrast, in the late warming period, WUE in the arid northwestern region may become smaller, owing to the decrease in NPP. Correspondingly, the WUE in the humid eastern area may decline because of the possible continued increase in ET. In addition, future climate changes may lead to an increase in the frequency, intensity, and duration of drought events [70,71]. Therefore, drought will cause a decrease in SM, which will affect the trend of change in WUE.
Although studies have shown that an ecosystem can respond to increased drought by increasing WUE [7,8], this issue may not be straightforward. Because the WUE of different vegetation types has various responses to drought, the WUE may be positively or negatively correlated with drought and show a quadratic function relationships [59]. In the arid region, the decrease in SM caused by drought could lead to a much greater decline in ET compared with the change in NPP (increasing or decreasing trend) ( Figure 4). Finally, a drought may lead to an increase in WUE in most areas of the arid region. In humid areas, the energy conditions mainly determine the evaporation process. Therefore, the decline in SM in the early drought stage might not cause a significant change in ET. At this time, the trend of change in WUE might mainly depend on the changing in NPP, which could show an upward or downward trend. With the aggravation of drought, SM could decrease significantly, which could cause a decline in ET. As a result, in the humid area, the original upward trend of WUE may be strengthened, and the original downward trend of WUE may be stabilized or even reversed in some cases.
To enhance the adaptability of ecosystems to future climate change or mitigate the adverse effects of climate change, it may be necessary to introduce the following suggestions:(1) In arid regions, the existing desert vegetation and grassland should be prevented from shrinking or degrading owing to intense human factors because these vegetation have a high ecosystem resilience [72]. (2) In arid areas with sparse vegetation, greening should be achieved by expanding grassland areas rather than artificial afforestation. C4 plants, which are usually dominated by herbaceous plants, have a higher internal WUE and can better adapt to the natural conditions of limited water and sufficient energy [73]. With the possible aggravation of drought in humid areas, expanding the C4 vegetation area can also be considered.

Conclusions
This study developed an analysis framework to quantitatively assess the influence of SM on WUE. Based on the analysis framework, this paper found that: (1) SM had an overall negative effect on the WUE in the terrestrial ecosystems during the last 30 years. In the NPP pathway, SM significantly affected WUE, and the sensitivity coefficient ranged from −68.66 to 43.19. However, in the ET pathway, SM completely negatively affected WUE. Under the dual effects of ET and NPP pathways, SM negatively affected WUE in 84.62% of the area, whereas it showed a positive influence on WUE in the remaining 15.38% of the area. Thus, although SM had a positive effect on the WUE of individual regions, it still had an overall negative effect on the WUE of each climatic zone. (2) The variation in SM could lead to a significant WUE variability, especially in Northeast and Northwest China. During the growing season in the past 30 years, the average and maximum variability of WUE caused by variation in SM were 28.02% and 109.29%, respectively. Regarding different climatic zones, SM in CTMP had the greatest effect on WUE, whereas it had the least impact on WUE in PTMP. In addition, the lower the forest covers, the higher the variability of WUE caused by changes in SM. The WUE of forest ecosystems is more resistant to changes in SM, whereas the WUE of non-forest vegetation is more sensitive to changes in SM.
These findings deepened our understanding of WUE changes and their underlying mechanisms, thereby offering important insight for predicting the response of ecosystems to climate change.

Conflicts of Interest:
The authors declare no conflict of interest.