A Hybrid Model Coupling Physical Constraints and Machine Learning to Estimate Daily Evapotranspiration in the Heihe River Basin

Abstract

Accurate estimation of surface evapotranspiration (ET) in the Heihe River Basin using remote sensing data is crucial for understanding water dynamics in arid regions. In this paper, by coupling physical constraints and machine learning for hybrid modeling, we develop a hybrid model based on surface conductance optimization. A hybrid modeling algorithm, two physical process-based ET algorithms (Penman–Monteith-based and Priestley–Taylor-based ET algorithms), and three pure machine learning algorithms (Random Forest, Extreme Gradient Boosting, and K Nearest Neighbors) are comparatively analyzed for estimating the ET. The results showed that, in general, the machine learning model optimized by parameters was able to better predict the surface conductance of the hybrid model. Driver analyses showed that radiation, normalized difference vegetation index (NDVI), and air temperature had high correlations with ET. The hybrid model had a better prediction performance for ET than the other five models, and it improved the R₂ of the two physical process-based algorithms to 0.9, reduced the root mean square error (RMSE) to 0.5 mm/day, reduced the BIAS to 0.2 mm/day, and improved the Kling–Gupta efficiency (KGE) to 0.9. The hybrid model outperformed the others across different time scales, displaying lower BIAS, RMSE, and higher KGE. Spatially, its ET patterns aligned with regional vegetation changes, with superior accuracy in annual ET estimation compared to the other models. Comparison with other ET products shows that the estimation results based on the hybrid model have better performance. This approach not only improves the accuracy of ET estimation but also improves the understanding of the physical mechanism of ET estimation by pure machine learning models. This study can provide important support for understanding ET and hydrological processes under different climatic and biotic vegetation in other arid and semi-arid regions.

1. Introduction

Land surface evapotranspiration (ET) is an important component of the energy balance and hydrological cycle [1,2] and plays an important role in the climate system, energy balance processes, and the carbon cycle [3]. Accurate monitoring and estimation of ET are not only crucial for water resource management but also for regional and global climate and hydrological cycle modeling [4]. Remote sensing technology is an effective means to monitor ET [5], and several remote sensing ET products have emerged, such as the medium-resolution MODIS global ET product [6], the Global Land Evaporation Amsterdam Model (GLEAM) [7], the Global LAnd Surface Satellite (GLASS) [8], ETMonitor [9], and the Breathing Earth System Simulator (BESS) [10]. Currently, fresh water is becoming increasingly scarce in many parts of the world, and it is important to improve the effectiveness measure of agricultural irrigation water [11]. The Heihe River Basin (HRB) belongs to the arid and semi-arid region of Northwest China, where irrigation is often required for agricultural production, and the portion of consumed water used for irrigated agriculture in these areas ranges from about 80% to 95% [12] and is particularly evident in Northwest China. Irrigation water use efficiency (IWUE) is often used by scientific researchers as a measure to characterize the irrigation effectiveness of farmland [13]. The quantification of ET can provide an important basis for improving irrigation efficiency and water management.

Quantifying IWUE at the farm scale first requires more accurate and high spatial and temporal resolution estimates of surface ET [14,15], and these measurements can provide a basis for decision-making in irrigation management [16,17]. In the ET study in the HRB, satellite sensor data were mainly used to estimate ET, such as the Moderate Resolution Imaging Spectroradiometer (MODIS), Advanced Very-High-Resolution Radiometer (AVHRR), and Visible Infrared Imaging Radiometer (VIIRS) to provide daily thermal infrared (TIR) data, but their resolution is mainly 1 km, which is too coarse for estimating ET at the field scale [18]. On the other hand, the Landsat series of thematic mappers and augmented line thematic mappers, as well as thermal infrared sensors and the Terra Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER), provide higher resolution TIR data (100 m), but they have a long round trip period, and clouds can interfere with the acquisition of imagery [19]. Therefore, estimating ET using data that combine the fine spatial resolution of Landsat-type data with the high temporal resolution of MODIS-type data will provide significant benefits for water resource management at the field scale [20,21,22,23].

In situ stations can provide accurate point-scale measurements of ET, but they are typically sparsely distributed throughout the study area, including in the middle and lower reaches of the HRB. There are significant limitations to accurately obtaining regional ET for the HRB. Remote sensing observations combined with process-based flux equations to quantify regional ET in the region is an effective strategy [24]. Conventional process-based algorithms mainly including the Peman–Monteith (P-M)-based [25] equations, Priestley–Taylor (P-T)-based equations [26], Surface Energy Balance System (SEBS) [27], Single Source models [28], dual-source models [29], and empirical/semi-empirical formulation methods [30,31] have been widely used to produce global and regional satellite ET products [6,32,33,34,35]. However, regional high-resolution ET quantification remains a challenging task due to the complex physical and biological controls on ET and the highly variable geophysical conditions in the HRB [36,37]. The multi-year mean (1982–2016) based on Penman–Monteith–Leuning version 2 (PMLV2) is 353 ± 24 mm/year. Process-based ET estimation models vary considerably when applied to the HRB [38,39,40,41]. One possible reason is that along with the substantial variation in ET over time, the static parameters of process-based physical models limit the dynamic capture of ET in different plant functional types, especially when heterogeneous surfaces are involved. With the increasing availability of satellite and ground-based observations, machine learning methods are increasingly used in hydrological studies, especially for regional and global scale ET prediction [42,43,44,45,46]. Data-driven machine-learning-based ET models can collect patterns and trends from satellite and meteorological data streams and automatically extract spatio-temporal features to provide better predictions [47]. However, these data-oriented estimation models for ET may not be suitable for derivation for regions with a different subsurface and sparse data because they are empirically based [48,49]. Therefore, machine learning models can be further modified and used as a basis to develop and optimize ET models to better fit the space of hydrological science applications.

The coupling of process-based models and machine learning provides alternative ways of modeling complex phenomena derived from natural physical processes [50]. Physical models and machine learning methods can be seen as complementary to each other rather than separate scientific paradigms [47]. Similar studies have already been conducted. For example, Zhao et al. [49] proposed a physically constrained hybrid model to preserve the surface energy budget of the ET process, which uses a modified PM equation as a loss function for an artificial neural network (ANN) to estimate ET at the site scale. However, an important question arose as to whether the hybrid model could be hybridized with a machine learning model to replace process-specific model formulations at the regional scale. Koppa et al. [51] trained a deep learning model to embed ET coercion into a process-based ET algorithm to produce a global hybrid ET model algorithm. However, the hybrid model focused on machine learning modeling of vegetation ET stress and was not applicable to sparsely vegetated, high-elevation regions. Also, these hybrid models need further testing and validation in the HRB, especially site-scale validation under high-resolution imagery.

For process-based physical modeling ET algorithms, Tan et al. [52] argued that global ET can be simulated more efficiently if reliable surface conductance (Gs) is obtained and quantified as a control on the environment, allowing for combination with meteorological variables [52]. The parameterization of Gs faces many challenges, as it is not only regulated by the physical environment but also varies between species [53,54]. Significant variations in leaf conductance, leaf-size-dependent boundary layer conductance, and canopy-structure- and roughness-dependent aerodynamic conductance can occur between biomes and species [55,56]. Wang et al. [30] used linearly fitted regressions to compute surface conductance and the complexity of the surface conductance mechanisms empirically. However, there is no consensus on a clear paradigm for describing complex feedback mechanisms, which leads to substantial uncertainty in estimating feedback [41]. This provides an opportunity for the coupling of machine learning and process-based modeling. In addition, P-M-based ET algorithms have been shown to have substantial errors in many studies, with the main challenge being to determine surface conductivity for different climatic, vegetation, and moisture conditions [6]. This is because surface conductance varies greatly with actual environmental and biological conditions, making it difficult to calibrate [57]. Despite studies conducted to estimate surface conductance [58,59,60], ET estimation shows great variability among different models using surface conductance due to insufficient theoretical foundations and large uncertainties. Therefore, machine learning methods provide an opportunity of integration to replace the traditional expression for surface conductance. The traditional expression for surface conductance not only requires more variables, but the calculation does not better fulfill the real geobiological conditions at the surface. However, the different hybrid models have not been substantially compared, and the link between the constrained physical mechanisms and the coupling performance is not yet clear.

The HRB is one of the world’s more hydrologically vulnerable regions, typical of many arid and semi-arid areas [61]. In recent years, water scarcity and conflicts between consumers in the HRB have become more intense due to increased irrigation water use in the middle reaches and environmental flow requirements in the lower reaches [61]. As a first step towards resolving these conflicts in the HRB, accurate monitoring of ET with high spatial and temporal resolution is necessary in the region [62]. In addition, analyzing experimental data in the HRB will make the region a test for additional examination of hydrologically related issues in other arid and semi-arid regions.

Based on the above introduction, whether it is the P-M equation based on large-leaf theory or the empirical/semi-empirical P-T equation and surface energy balance models, they lack strong physical constraints on vegetation ET or more accurate aerodynamic conductance. Meanwhile, current pure data-oriented machine learning models generally select the factors affecting ET for fitting and lack the study of aerodynamic conductivity factors. This study explores a more promising hybrid modeling strategy for ET and emanation, applying machine learning to complement the physical model. In this paper, we couple physical constraints and machine learning methods to develop a hybrid ET estimation model based on surface conductance optimization and compare and analyze three kinds of models, one process-based physical algorithm, one data-oriented pure machine learning algorithm, and one model that couples physical constraints and machine learning algorithms, based on EC observations at flux tower sites in the upper and middle reaches of the HRB; firstly, we explore how the machine learning algorithm fits the performance of surface conductance in the physical framework and perform parameter optimization of the machine learning model to obtain the best fitting capability, analyze the correlation between input variables such as atmospheric drive data and optical sensor data and site-observed ET, select variables with high correlation with ET as inputs to the machine learning model, and then input them into the machine learning model to predict ET. The ET accuracies predicted by the six models are compared based on the station observation data, and the trend comparison analysis and deviation analysis of the predicted ET are conducted; the performance of all the models is also analyzed on three-time scales: daily, monthly, and annual. Finally, the surface ET in the study area is mapped at a high spatial resolution to analyze their spatial and temporal performances at a relatively high resolution, and the performance of different models is analyzed by the comparison of Taylor diagrams; the ET predicted by the hybrid models are analyzed in comparison with other remote sensing ET products. This will provide support for the optimal allocation and management of water resources in the HRB and similar semi-arid and arid regions.

2. Model Description

This study presents a hybrid model of physical constraints and machine learning based on surface conductance optimization. This model has an explicit physical framework where a parameter or variable component is modeled using a machine learning model. The input variables for the pure machine learning model are incident solar radiation (Rs), air temperature (Ta), relative humidity (RH), NDVI, and minimum air temperature ($T_{min}$). These variables were chosen because they are components of the key parameter species of the moisture ET mechanism and were also filtered by the machine learning prediction of ET impact factor analysis in Section 4.2 below. The model involves a wide range of data sources, including meteorological data, optical satellite data, and others. The logic of the ET algorithm is depicted in Figure 1:

2.1. Surface-Conductance-Based Machine Learning

To be able to utilize meteorological data, Wang et al. [30] provided a semi-empirical formulation based on the P-M framework by using both meteorological information and remotely sensed observations to monitor ET globally on an interannual scale:

(1)$\lambda \mathrm{E}=\lambda {\mathrm{E}}_E+\lambda {\mathrm{E}}_A$

(2)$\lambda {\mathrm{E}}_E=\frac{\mathsf{\Delta }}{\mathsf{\Delta }+\gamma }\cdot \left(R_n-G\right)\cdot g_s$

(3)$\lambda {\mathrm{E}}_A=\frac{\mathsf{\Delta }}{\mathsf{\Delta }+\gamma }\cdot \mathrm{V}\mathrm{P}\mathrm{D}\cdot g_a\cdot g_s$

In Equation (1), $\lambda \mathrm{E}$ is the latent heat flux ($\mathrm{W}/{\mathrm{m}}^2$), $\lambda {\mathrm{E}}_E$ is the energy control on the ET ($\mathrm{W}/{\mathrm{m}}^2$), $\lambda {\mathrm{E}}_A$ is the atmospheric control on the ET ($\mathrm{W}/{\mathrm{m}}^2$), $\lambda$ is the latent heat of evaporation (J/kg), and $\mathsf{\Delta }$ is the slope of the saturated vapor pressure versus temperature curve ($\mathrm{K}\mathrm{p}\mathrm{a}/\mathrm{℃}$). In Equations (2) and (3), $R_n$ is the net radiation, $G$ is the soil heat flux ($\mathrm{W}/{\mathrm{m}}^2$), where $g_s$ and $g_a$ are the surface and aerodynamic conductivity (m/s), respectively, and VPD is the vapor pressure difference (Pa) of the air. The semi-empirical formulation of this model is as follows:

(4)$\lambda \mathrm{E}=\left[\frac{\mathsf{\Delta }}{\mathsf{\Delta }+\gamma }\cdot \left(R_n-G\right)+\frac{\gamma }{\mathsf{\Delta }+\gamma }\cdot \mathrm{V}\mathrm{P}\mathrm{D}\cdot g_a\right]\cdot g_s$

where $g_a$ in Equation (2) can be parameterized by wind speed ($\mathrm{W}\mathrm{s}$) [63,64,65]. Wang et al. [40] used relative humidity deficit (RHD) and vegetation index (VI) for surface conductance parameterization. The algorithm worked well when $R_n$ and stomatal conductance data were available, and the vegetation was not under water stress [8,31]:

(5)$g_a=0.26\cdot \left(1+0.54\cdot \mathrm{W}\mathrm{s}\right)$

To improve the above empirical equations, we used a machine learning model fitted with surface conductance as the input variable into the empirical expressions. To train the machine learning model, the input variables included Rs, Ta, RH, NDVI, and $T_{min}$. We had no available ground observations for $g_s$, and we used the latent heat flux (LE) obtained from EC observations to infer $g_s$ based on Equations (1)–(4); the inferred $g_s$ served as the target variable for training the machine learning model, and the daily ET estimates for this algorithm were generated at a spatial resolution of 100 m. We also used the latent heat flux from the EC observation towers to infer $g_s$ from the energy balance equations, and the inferred $g_s$ served as the target variable for training the machine learning model.

2.2. Data-Oriented Machine Learning Model

We chose the Random Forest (RF) model as a regression model for fitting the surface conductance and a pure machine learning model for predicting ET. The RF model is a machine learning method for integrated learning. Predictions are made by constructing multiple decision trees and combining them. Also, the RF model is a collection of many classification and regression models whose main function is to produce accurate predictions and generally does not overfit the data. It also belongs to the combination of tree predictors, which depend on the values of independently sampled random vectors and have a similar distribution for all trees in the forest [66,67]. That is, the model combines the ease of interpretation and robustness of decision trees with the high performance and resistance to overfitting of integrated learning. RF-based machine learning models have been shown to predict ET overall better than generalized regression neural networks (GRNNs) in southwest China [68].

Meanwhile, we used Extreme Gradient Boosting (XGBoost 1.6.1) and K Nearest Neighbors (KNNs) algorithms as references for pure machine learning models. These three pure machine learning models were trained and fitted using the same set of input combinations and the same set of data. Fan et al. [69] used XGBoost to study the daily reference ET across China in small-scale numerical weather prediction, and they showed that the Extreme Gradient Boosting method can reduce the BIAS of predicting daily ET, which is the same as that in this paper, which also investigates the ET at the daily scale. Similarly, Sevim Seda Yamac et al. [70] used the KNN model to study crop ET and compared artificial neural network and adaptive enhancement models, and they found that the KNN model performs the best, which is why KNN was used in this study.

In this study, we selected data from three sites, Arou, Dashalong, and Daman, in 2013 and three sites, Yingke, Linze, and Zhangye, in 2015 to be used for training and validation. Among them, Arou and Dashalong are two independent validation sites, and the other four sites are training sites. During the training process, we selected 30% of the data as validation data and 70% of the data as training data. At the same time, we selected six sites in other available years of data to predict and analyze model performance.

2.3. Comparison with Process-Based ET Algorithm

In this paper, a performance study was conducted using two process-based ET models driven by satellite and meteorological variables, modeled as the P-M algorithm and P-T algorithm. Each is described in detail below.

• (1). P-M algorithm

We used the algorithm recommended by the Food and Agriculture Organization of the United Nations (FAO) for calculating ET adoption, which is a basic P-M equation [25].

(6)${\mathrm{E}\mathrm{T}}_0=\frac{0.408\mathsf{\Delta }\left(R_n-G\right)+\gamma \frac{900}{T_{mean}+273}{U}_2\left(e_s-e_a\right)}{\mathsf{\Delta }+\gamma \left(1+0.34U_2\right)}$

In Equation (6), ${\mathrm{E}\mathrm{T}}_0$ is the reference crop ET (mm/day), $R_n$ is the net radiation ($\mathrm{M}\mathrm{J}/{\mathrm{m}}^2\mathrm{d}\mathrm{a}\mathrm{y}$), G is the soil heat flux density ($\mathrm{M}\mathrm{J}/{\mathrm{m}}^2\mathrm{d}\mathrm{a}\mathrm{y}$), $T_{mean}$ is the mean atmospheric temperature (°C), $e_s$ is the saturated vapor pressure (kPa), $e_a$ is the actual vapor pressure (kPa), Δ is the slope of the saturated vapor pressure curve (kPa/°C), $\gamma$ is the hygrometer constant (kPa/°C), and $U_2$ is the wind speed (m/s) at a height of 2 m. All parameters needed to calculate the ET follow FAO56 Chapter III. $\gamma$ is the hygrometer constant (kPa/°C), and $U_2$ is the wind speed (m/s) at a height of 2 m. The parameters required for the calculation of ET follow the methodology and procedures of FAO56 Chapter III. We mainly used $R_n$, Ta, $T_{min}$, RH, and normalized vegetation index (NDVI) as forcing data to produce ET_PM products. Among them, the Rn was calculated by us using a similar method to that of Wu et al. [71]:

(7)$R_n=\left(1-albedo\right)\cdot \mathrm{R}\mathrm{s}+\left({\epsilon }_a-1\right)\cdot \sigma \cdot \frac{T_{max}^4+T_{min}^4}{2}$

In Equation (7), albedo is the surface albedo, Rs is the downward shortwave radiation ($\mathrm{W}/{\mathrm{m}}^2$), ${\epsilon }_a$ is the surface specific emissivity, $\sigma$ is the Stefan-Boltzmann constant (5.67 × 10⁻⁸ $\mathrm{W}/\left({\mathrm{m}}^2{\mathrm{K}}^4\right)$), and $T_{max}$ and $T_{min}$ are the maximum and minimum air temperatures, respectively (°C). Wu et al. [71] provided the equation, ${\epsilon }_a$, calculated as follows:

(8)${\epsilon }_a=0.95+0.01LAI$

For the soil heat flux G, Su and Lima et al. [27,72] considered the best formula for the soil heat flux as follows:

(9)$\frac{G}{R_n}=\frac{LST}{albedo}\left(0.0038·albedo+0.0074·{albedo}^2\right)\left(1-0.98{NDVI}^4\right)$

In Equation (9), G is the soil heat flux, albedo is the surface albedo, LST is land surface temperature (°C), NDVI is the normalized vegetation index, and $R_n$ is calculated in the above equation.

• (2). P-T algorithm

(10) $\lambda {\mathrm{E}}_{\mathrm{P}\mathrm{T}}=\alpha \cdot \frac{\mathsf{\Delta }}{\mathsf{\Delta }+\gamma }$

In Equation (10), $\alpha$ is the correction factor, which is generally between 1.0 and 1.3, and Yao et al. [73] considered that $\alpha$ was selected to be 1.26. The forcing data for the P-T equation were mainly $R_n$, Ta, RH, Press, and NDVI. Here, Rn was also calculated using Equation (7) above.

3. Data and Model Validation

3.1. Study Domain and In Situ Observations

The HRB is located in northwestern China, geographically between 97°E and 102°E in longitude and 38°N and 43°N in latitude, with a total area of about 142,900 square kilometers. The climate of the region is influenced by the circulation of the westerly wind belt in the middle and high latitudes and the polar cold air masses, showing complex and variable climatic characteristics. The geomorphology in the basin is diverse, presenting diversified natural landscapes such as ice/tundra, forests, meadows, artificial/natural oases, deserts, and lakes from the upper to the middle and lower reaches in turn [74]. This region is characterized by both cold and arid zones, with significant surface heterogeneity, especially the mountainous cryosphere and the extreme aridity of the river tailrace areas forming a clear contrast. The region has become an important commercial grain producer in the northwest. Irrigation water in the irrigation area mainly comes from surface water of the main stream of the Heihe River and is supplemented by groundwater. The HRB has become one of the most important areas for climate research in northwest China because of its representativeness and rich research base.

The Heihe integrated observatory network was established in the HRB in 2007, during the WATER (Watershed Allied Telemetry Experimental Research) experiment (2007–2011), and was completed in 2013 during the HiWATER (Heihe Watershed Allied Telemetry Experimental Research) experiment (2012–2015) [75]. At present, many scholars in the HRB have carried out long-term eco-hydrological experiments and an integrated observation network established in the HRB, accumulating a large amount of hydrometeorological observation data [76,77], and there are currently 15 observation stations in operational operation in the HRB, including 3 super-stations and 12 ordinary stations, which cover a variety of subsurface in the upper, middle, and lower reaches of the HRB, providing an important support and foundation for an in-depth understanding of the regional water cycle and ground–air energy exchange [78]. The study area was selected as the range of Landsat series imagery, and in the data preprocessing stage, we evaluated the availability and temporal continuity of the measured data at the stations in the middle and upper reaches of the Heihe River, taking into account that this study requires the use of data of different temporal periods for the training and validation of the model, as well as for the study of evapotranspiration in the middle and upper reaches of the Heihe River under different surface covers. In this study, the site-scale validation and analysis studies were mainly carried out at six sites in the middle and upper reaches (Figure 2), and the site details are shown in Table 1.

3.2. Satellite and Meteorological Datasets

After years of development, Ma et al. [19] investigated the ET by temporally and spatially fusing ETM+ data and MODIS data and then studying the ET in the HRB. We used the weather-driven data and temperature data therein, as well as optical data from optical sensors. These data have the advantage of high spatial and temporal resolution, which is beneficial for minimizing the effects of scale. The data series provides meteorological forcing data, mainly including maximum air temperature ($T_{max}$, minimum air temperature ($T_{min}$), atmospheric pressure (Press), RH, wind speed ($\mathrm{W}\mathrm{s}$), and downward shortwave radiation (Rs), and the fused products also include NDVI and albedo data. This dataset has great advantages for high-resolution quantitative ET estimation in the HRB, provides multi-source data characterization for the study of the HRB using different data sources, can effectively predict continuous daily ET at Landsat-like scales (100 m resolution) [19], and can better test the hybrid models and process-based physical models’ performance. To be able to validate this study at the site scale, we visited the National Tibetan Plateau Data Center (https://www.tpdc.ac.cn/) on 23 March 2024 to obtain meteorological and flux data for each site and to calculate the energy balance equation to derive site ET data for model validation and comparison.

To investigate the optimization effect in this study, we selected five widely used and popular remote sensing ET products (Table 2), including the GLEAM and DTD based on the surface energy balance equation, the MODIS16 and ETMonitor based on the physiological and ecological characteristics of vegetation, and the integrated product GLASS. In this study, the time scales of all the products were standardized to the monthly scale, and the spatial resolution was resampled to the 1 km resolution to facilitate the direct comparison of the different products.

3.3. Model Validation and Evaluation Metrics

• (1). Validation metrics

We chose to use multiple metrics to evaluate all models. First, we used the correlation coefficient ($R^2$) to assess the correlation between the models and actual observations to determine the degree of model fit. Second, we used RMSE and BIAS (BIAS) to measure the model’s predictive accuracy and variance from observations. However, to assess model performance more comprehensively, we introduced the KGE metric [81]. This is because KGE combines correlation, relative variability ratio ($\alpha$), and mean ratio ($\beta$) to facilitate a better understanding of the mismatch between estimates and observations. Where ${\sigma }_p$ and ${\mu }_p$ are the standard deviation and mean of the prediction, respectively, and ${\sigma }_o$ and ${\mu }_o$ are the observed values, the model performs better when the KGE is close to one. In the study of ET on the Tibetan Plateau, Shang et al. [82] evaluated their model well using metrics such as KGE.

The $R^2$ is given by the following:

(11)$R^2=1-\frac{{\sum }_i{\left(\widehat{y}-y_i\right)}^2}{{\sum }_i\left({\overline{y}}_i-y_i\right)}$

The formula for the RMSE is as follows:

(12)$RMSE=\sqrt{\frac{1}{m}\sum _{i=1}^m{\left(y_i-\widehat{y}\right)}^2}$

The BIAS formula is as follows:

(13)$BIAS=\frac{1}{m}\sum _{i=1}^m\left(y_i-\widehat{y}\right)$

The KGE equation is as follows:

(14)$KGE=1-\sqrt{{\left(r-1\right)}^2+{\left(\alpha -1\right)}^2+{\left(\beta -1\right)}^2}$

(15)$\alpha =\frac{{\sigma }_p}{{\sigma }_o}$

(16)$\beta =\frac{{\mu }_p}{{\mu }_o}$

• (2). Machine learning model interpretability analysis

To gain insight into the performance of machine learning models and the interpretability of their predictions, we employed the Shapley value as an analytical tool. The Shapley value is a method based on the theory of cooperative games for evaluating the extent to which each feature contributes to the predictions of a model. K. Aas et al. [83] even aggregated various shapely calculations to evaluate models, and by calculating the Shapley value of each feature for each prediction, we can quantify the degree of influence of different features on the model predictions. This approach not only helps us understand how the model makes decisions but also identifies key features that can guide subsequent feature engineering and model optimization. In this study, we used Shapley values to analyze the performance of the machine learning model we built on different datasets and explore the impact of individual features on model performance. Through this analytical approach, we can understand the prediction process of the model more comprehensively, which provides valuable reference and guidance for further model improvement.

For sample $x_i$, the value of feature j is $x_{ij}$. The model’s predicted value for this sample is $y_i$, and the model’s baseline (usually the mean of all sample target variables) is $y_{base}$. This can be expressed as follows:

(17)$y_i=y_{base}+f\left(x_{i1}\right)+f\left(x_{i2}\right)+\cdots +f\left(x_{ij}\right)$

Here, f($x_{ij}$) is the SHAP value of $x_{ij}$.

4. Result

4.1. Model Optimization for Gs

To be able to evaluate the performance of machine-learning-simulated surface conductance in a hybrid model, we tested the model’s generalization performance out-of-sample by selecting one EC site for independent validation under each of the three different subsurface mat coverages in this study area. We used grid search validation to tune the parameters of the model to obtain the best model parameters. As shown in Figure 3, the validation results show that the machine learning model (RF) with tuned parameters (Figure 3a) predicts better Gs, and its fitted regression line is closer to the 1:1 line than the prediction in Figure 3b, with a KGE of 0.87, an R2 of 0.92 (p < 0.01), a BIAS of 0.00 m/s, and an RMSE of 0.05 m/s. In contrast to the machine learning model without the tuned-parameters scheme machine learning model (Figure 3b), $R^2$ is 0.86, BIAS is 0.03, RMSE is 0.08, and KGE is 0.74. The results in Figure 3 show that the machine learning model with tuning parameters performs well at the out-of-sample independent validation site with better generalization ability.

4.2. Machine Learning Variable Screening

To explore the performance of the pure machine learning model, we plotted a heat map of the feature variables of the pure machine learning model, in which we used the field site observation of ET data as the target variable and the other meteorologically driven data, and optical data as the feature variables. Because pure machine learning models are strongly influenced by the quality of sample data, we wished to explore the main drivers of the prediction performance of pure machine learning models and retain the main influencing factors as input variables to data-oriented pure machine learning algorithms for predicting ET. As shown in Figure 4, among the seven main input variables, the calculated Rs, $T_{min}$, and NDVI had high correlations with the observed ET at the flux site, with correlations of 0.79, 0.90, and 0.87, respectively, while the correlation of other meteorological driving data, such as Press, RH, and Ws, was weak. Among them, the correlations of Ws and Press with measured ET were lower than 0.1, while RH was negatively correlated with measured ET, which was generally consistent with the study of Han et al. [84]. For the input variables for predicting the surface conductance in the Random Forest model in the hybrid model, Rs, NDVI, $T_{min}$, Ta, and RH were selected as input variables to predict the surface conductance with reference to the results of shang et al. [82], since the surface conductance data were not available as measured data but were derived by derivation.

4.3. Model Regression Analysis of Verification Site

To explore the predictive performance of the hybrid models at each site in more depth, we analyzed the fitting performance of each model and observation at the site scale and plotted scatter density regressions. Each of the six EC tower sites in this study area has more than half a year of available ET observations from 2011 to 2016, which are used to compare the numerical scatter density regressions between three kinds of different models. Figure 5a shows that on the Daman site where the underlayment is farmland, the hybrid model has an $R^2$ of 0.92, a KGE as high as 0.89, a BIAS of 0.17, and an RMSE of 0.78, with the KGE being the best-performing model among the six, which is very informative for fitting hydrological models for ET. For the Arou site with forested underlayment as in Figure 6b, the correlation coefficient predicted by the hybrid model is higher at 0.92, which is a little lower than the correlation of ET_PM and KNN, but the hybrid model has a greater advantage for the KGE, which is an indicator of the superiority of hydrological models. For the Dashalong site with a grassy sub-surface as in Figure 6c, the ET observation values are low relative to the farmland area, and the machine learning model predictions, although with high correlation coefficients, deviate from the accurate values, which are not as accurate as the hybrid model predictions with higher KGE. Meanwhile, the KGE of the data-oriented pure machine learning model is less effective. For the area where the subsurface is wetland and there is a desert around, the ET predicted by the hybrid model at the Zhangye site as shown in Figure 6a is fitted better, with a KGE of 0.77, and the deviation from the observed value of the on-site EC tower is 0.55. Although this absolute value is relatively large, it is still better than the two process-based physical models and the pure data-oriented machine learning model. The RMSE of the predictive performance of all six models at this site is larger, and the hybrid model is more competitive. For the Yingke site, which is closer to the Daman site as shown in Figure 5c, the correlation coefficient between the predicted values of the hybrid model and the in situ observations is as high as 0.94, with a KGE of 0.81; the machine learning model performs well at this site, with a better-fitting correlation coefficient and an overall performance similar to that of the hybrid model. At the Linze site, which also belongs to the farmland area as shown in Figure 5b, the correlation coefficient between the predicted values of the hybrid model and the ET field observations is as high as 0.93, which is a significant advantage over the pure machine learning model of 0.76. In comparison, the process-based physical algorithms are also more adapted to the local vegetation and environment. These site-based scatter density regression plots show that while the current pure machine learning models have a low BIAS in some areas, the pure machine learning models still have high RMSE in some areas. The prediction performance of pure physical process-based models is variable and unstable across sites and surface vegetation types, and their static parameters limit the variation in ET across different vegetation types and heterogeneous surfaces [44]. In contrast, the hybrid model predicts more stable performance at different sites with better generalization performance, which enhances the ability to predict ET under different vegetation types with stronger advantages.

4.4. Temporal Assessment of Different Models

We wanted to further investigate the ability of different models to simulate daily ET, and the annual variation in ground-based measurements and the ET estimates from the six models, which are annually correlated, are plotted in Figure 7. As shown in Figure 7, for the Arou site where the surface vegetation type is forest, the hybrid model (RF_Gs) can fit the trend and values of measured ET well during the summer months, most notably during the time from May to August. The hybrid model was able to accurately extract not only the low values corresponding to the time period, but also accurately predict the peaks and trend changes for that time period. Similarly, for the Dashalong site, where the surface vegetation cover is grass, the ET at the site was relatively small in general, and the hybrid model and the process-based physics algorithm ET_PM underestimated the ET in the field from April to October. The pure data-oriented machine learning approach, on the other hand, overestimated ET at this site. For the Daman site, which has an area of agricultural land as surface cover, none of the models accurately predicted the single peak trend with relatively low values from March to May. The process-based physical algorithm ET_PM showed more overestimation from February to July for this site. For the Zhangye site, all models underestimated the fit to the February to July data, and the hybrid model predicted a fluctuating trend from March to May in the face of a small peak at this site, but did not fit the field-observed ET well. For the Yingke site, the mixed model showed some overestimation from June to August. The results of these site-based trend evaluations indicate that although current ET estimation models deviated somewhat from calculated ET from in situ flux tower observations, they generally shared the same trends and showed good agreement with observations. Among the six validation sites, the RF_Gs hybrid model enhanced ET estimation in areas with plant functional types in the middle and upper HRB compared to ET_PM and ET_PT. The pure machine learning models showed similar fitting ability to the process-based pure physical models at the site scale.

4.5. Time-Scale Performance Analysis of Different Models

Figure 8 and Figure 9 summarize the comparison of the performance evaluation metrics ($R^2$, BIAS, RMSE, KGE) of the hybrid model (RF_Gs), the pure machine learning models (RF, XGBoost, KNN), and the process-based physical models (ET_PM, ET_PT) at the daily, monthly, and annual scales and across all sites. The hybrid model was overall better than the other five models on daily, monthly, and annual scales. Furthermore, it can be concluded from the results that the three pure machine learning models have a high degree of similarity, as do the process-based physical models ET_PM and ET_PT. All the data in the figure were computed by comparing them with the observed data, whether on daily, monthly, or annual scales. For the deviations, the process-based physical models all exhibited negative values, with anomalies as high as −1.25 mm/day, whereas the pure machine learning models had negative deviations, but the absolute values of the deviations for all three were smaller compared to the process-based physical models. For the daily ET estimates, the median values of the deviations of the pure machine learning models ranged from −0.5 mm/day to −0.25 mm/day, which is a relatively good performance, whereas the median values of the deviations of the process-based physical models were mainly centered on −0.75 mm/day, due to the variations in the crop coefficients with the surface vegetation and the biochemical conditions, and −0.50 mm/day, due to the variations in the crop coefficient with the surface vegetation and biochemical conditions, to −0.50 mm/day, which is roughly twice that of the pure machine learning model. On the contrary, for the RMSE, the pure physical process-based models had larger values, with the ET_PT model mainly concentrating the RMSE above 1.5 mm/day due to the coarser semi-empirical formulas and coefficients, and the ET_PM model, although affected by the constancy of crop coefficients, mainly concentrating the RMSE at 1.5 mm/day, which leads us to conclude that in the middle and upper reaches of the HRB, the ET_PM needs to be adapted to accommodate surfaces with different vegetation cover types. While the pure machine learning model is in the middle, the main median values are all concentrated at 1.0 mm/day, which has less impact for the summer and fall seasons, but the estimation accuracy needs to be strengthened for the winter and spring seasons, when ET is small, and the hybrid model shows a relatively good performance. The error values of the process-based physical and machine learning models are relatively large due to the relatively large error values on the daily time scale and are even larger due to the significant cumulative effect on the monthly and annual time scales.

Since the plotted images are a mixture of data from all sites, for the $R^2$, which assesses the model fitting performance, and the KGE, which assesses the model performance, a single model with high $R^2$ at a single site does not show a strong advantage in general, but rather a hybrid model shows a much better performance in general because of the overall good prediction of ET at individual sites for different vegetation covers and shows strong superiority. For daily ET estimation, the pure machine learning model XGBoost performs very consistently and can maintain an $R^2$ of around 0.9 across all sites, while RF performs second best, but with a slightly higher model $R^2$ than the XGBoost model. As for the two process-based physical models, we believe that because the estimation is affected by the non-clear-sky weather conditions of the remotely sensed weather-driven data, and in the case of the satellite sensors being affected by severe and extreme weather, etc., some outliers are still unavoidable, despite the data fusion process.

4.6. Comparative Analysis of Model Residuals

To analyze the performance of the models from other different perspectives and to compare the differences, we plotted Taylor plots for six sites with long time series to quantify the prediction accuracies of the six models based on the computed model predictions versus the observed values of the six sites. As shown in Figure 10, for the Daman site, where the subsurface is covered with agricultural land, all three machine learning models have small standard deviation and excellent RMSE, while the process-based physical model performs weaker than the three pure machine learning models and the hybrid model. We believe that the reason for this result is that the large-leaf theory of the process-based physical model cannot accurately detect the ET of crops in this region, while the pure machine learning model can better extract the characteristics and change trend of ET of farmland crops in this region during the training process. At the same time, the hybrid model still has stronger competitiveness and greater advantages. For the Arou site with a relatively high elevation and grassland cover type, the RMSE of all six models is superior, but the hybrid model has a stronger prediction ability, due to the combination of the mechanism of the physical model and the learning ability of the machine learning model, and performs the best in terms of RMSE and correlation coefficient. For the Dashalong site, where the surface vegetation cover is grass, the ET_PM performs well and is very close to the performance of the hybrid model. On the contrary, the ET_PT performs relatively poorly, which we attribute to the fact that the semi-empirical coefficients used in this study solidify the predictive power of the model and its ability to generalize decreases accordingly. The three pure machine learning models are similar across all sites in Figure 10, with an overall good agreement in terms of deviation from observations and correlation for the prediction of ET in this study area. While the process-based physical and hybrid models have similar predictive ability at some sites, they are even worse than the data-oriented machine learning models at other sites, which we attribute to the strong automatic feature extraction and robust nonlinear simulation performance of the pure data-oriented machine learning models for sample data at the site scale.

4.7. Impact of the EC Flux Tower Density and Sample Representativeness

In previous studies, ET is usually validated at the site scale. However, the density and absolute number of site distributions at different ground covers may affect the predictive performance of the model. For pure machine learning models, they are very dependent on the quality of the input sample data and are also strongly influenced by the site distribution. For process-based physical models, their performance depends on the accuracy and spatial and temporal resolution of the input data, while the remote sensing data generate errors during the scale conversion process that are transferred to the physical model. Among the validation sites used in this study, the site density of farmland is greater than that of other surface vegetation type areas. The KGE of the hybrid model at Zhangye, Yingke, and Daman sites is as high as 0.8, and the KGE of the pure data-oriented machine learning model is more variable. And it is closer to the hybrid model at Yingke and Daman sites compared to the process-based physical model, while at the Zhangye site, the KGE is lower than 0.6, while the process-based physical model has a smaller KGE. The performance of the KGE of the six models at different time scales for all surface covers in this study area is shown in Figure 9. From the six sites in Figure 11, it can be reasoned that the process-based physical model ET_PT has an unstable KGE for hydrological model assessment metrics under different surface covers due to its semi-empirical coefficients. In contrast, the ET_PM model with a stronger physical mechanism has a stronger relative predictive performance and can better respond to the nonlinear characteristics and physical mechanisms of ET. The pure data-oriented machine learning model is greatly affected by the absolute number of sample data and site density under different subsurface layers and performs poorly in the Zhangye site with wetland and the Dashalong site with grassland. We attribute this to the small amount of sample data at the Dashalong site in this study and the unavailability of data at the beginning of the year at the Zhangye site that affected the site predictions of the machine learning model.

We also compared the performance of the hybrid model with some common ET products (DTD, GLEAM, MODIS, GLASS, ETMonitor). The results show that the hybrid model performed relatively well at the six flux tower sites. As shown in Figure 12, overall, the hybrid model had a KGE of 0.79, ETMonitor had a KGE of 0.54, DTD had a KGE of 0.50, GLASS had a KGE of 0.46, GLEAM had a KGE of 0.41, and MODIS had a KGE of 0.33. For the product species compared, ETMonitor performed the best in the forest region, DTD in the farmland region, and in the grassland region, the KGEs of all the compared products were relatively low. Secondly, GLASS products outperformed GLEAM products in the farmland and grassland areas, while GLEAM products outperformed GLASS products in the forest area, but both of them performed relatively well. MODIS products performed relatively poorly in all comparisons.

4.8. Map of the HRB ET from the Hybrid Model

Figure 13 shows the spatial distribution of ET in the middle and upper reaches of the HRB for the Landsat series and Moderate Resolution Imaging Spectroradiometer (MODIS) fusion data under all models. The ET in the upstream region is relatively high with increasing elevation. It can be seen that the three machine learning models exhibit consistent spatial distribution patterns, but with different predictions for desert regions. In general, ET is relatively high in the upstream region and farmland area and relatively low in the surrounding desert and bare land. These findings are consistent with those of Ma et al. [19]. The maximum value of ET was around 1000 mm/year under the process-based physical model, but the area with the maximum value of ET exceeding 1000 mm/year was relatively large under the pure machine learning model. Both the process-based physical model and the hybrid model showed that the minimum value of ET was less than 200 mm/year, while the pure machine learning model predicted a minimum value of ET slightly greater than 200 mm/year.

The spatial variability of all ET models is consistent with changes in vegetation as well as climate. For the upper of HRB, as well as the midstream farmland area, water resources are more abundant, and ET is higher. Despite the great consistency in spatial patterns, the six ET estimates show regional differences and uncertainties. The annual ET of the six models can be categorized into the pure machine learning model annual ET belonging to the same relatively large range, while the process-based physical model and the hybrid model have another range of relatively small annual ET. The spatial prediction of ET by the pure machine learning model is significantly weaker in sparsely vegetated areas, with larger errors. The process-based physical and hybrid models, on the other hand, were in good agreement with the spatial distribution pattern of ET in the HRB studied by Ma et al. [19].

5. Discussion

5.1. Performance of the Hybrid ET Model

By integrating physical mechanism constraints with machine learning methods, the hybrid model retains the logical processes of physical mechanisms and produces more accurate ET estimates. The ET estimates from the hybrid model are in good agreement with ground measurements and are consistent with annual variations in ground measurements. As shown in Figure 7 at the Dashalong site, the hybrid model underestimated daily ET during the growing season, and ET_PM similarly underestimated ET during the growing season, mainly because during the growing season, ET peaked at the site due to intensive irrigation and strong human intervention in water resources, while the pure machine learning model overestimated annual ET. In the growing seasons of the Daman and Zhangye sites, the hybrid models were able to fit the trend and peak of ET better, i.e., the hybrid models were able to over-represent the irrigation-induced changes in temperature and moisture. Therefore, the hybrid model can better fit the ET in agricultural areas under intensive irrigation conditions.

For grassland, ET_PT overestimated the daily ET during the plant growing season from May 2015 to September 2015 in the one-year ET estimation, while the ET_PM model and the hybrid model showed a trend of less underestimation. Previous studies have shown that for complex topography and heterogeneous landscapes, physical model estimates show some variability and uncertainty when compared to measured calculated values from in situ EC flux towers [38]. The hybrid model, on the other hand, was closer to the data observed in the field, which reflects that the coupled machine learning method hybrid model had better spatio-temporal feature extraction and simulation capabilities [47].

For the area with forested ground cover, the RF_Gs model was not only able to predict its peak trend in the growing season but also able to fit the fluctuating trend in the growing season under some special weather conditions, which reflects that the hybrid model can simulate the ET changes under spatio-temporal features very accurately in the upstream vegetation area. The ET_PT model, on the other hand, suffered from a serious overestimation of daily ET in the forest region, which suggests that a simplified empirical formula such as PT is not applicable in regions with high vegetation indices. The pure machine learning model also had some overestimation, but the fit was competitive. This indicates that the pure data-oriented machine learning model has some automatic feature extraction capability and some nonlinear simulation performance.

Overall, traditional process-based physical models tend to perform poorly in arid regions where moisture conditions are not well represented by the surrounding arid areas [5,85]. This result is basically the same as that of Shang et al. [82]. In addition, the underestimation of ET_PM and ET_PT mainly occurs during the canopy dieback season or in areas with sparse vegetation cover, while the hybrid model performs better in all cases. The validation results suggest that the hybrid model improves the simulation capability of the physical model and can be used for a wider range of ecosystems in the HRB. For areas outside of the HRB where there are available flux tower sites, the parameters can be re-tuned to train the model to fit the values of ET observed in the field.

5.2. Interpretability of Hybrid Machine Learning Models

While improvements to the performance of a model are important, it is not enough to improve just one parameter of the model. In hydrometeorological applications, the interpretability and logic of machine learning models are as important as the prediction accuracy. In this study, we use Shapley values to analyze the interpretability of machine learning coupling in hybrid models. We quantify the proportion of Shapley values for each variable as a whole, as well as its contribution to the model. As shown in Figure 14, the results indicate that the atmospherically driven data and temperature data, as well as the optical sensor data, exhibit different contributions to the hybrid model. For the hybrid model, the contribution of downward shortwave radiation (Rs) to the hybrid model is 38%. This is followed by the normalized vegetation index (NDVI) contribution of 24% for the hybrid model. From the perspective of process-based biophysical parameters, RH accounts for 20%, and these three comprise the main influential components of the hybrid model predictions in this study. The sum of air temperature data (Temperature) also accounts for 18%, indicating that temperature plays an important role in surface conductance in this study. The study shows that surface conductance is mainly affected by moisture (RH), Ta, and VPD [86,87]. Mu et al. [6] attempted to solve for canopy conductance using water cover and leaf area index. Wang et al. [30], on the other hand, parameterized surface conductance as a linear function of RH and vegetation index (NDVI) to estimate ET in a global region. It can be concluded from the interpretable display of the hybrid model of Shapley values (RF_Gs) that the present hybrid algorithm has a common correlation with the above algorithms.

5.3. Uncertainties and Future Research Directions

Although the hybrid model has a better performance than the other models, at the same time, it can have a better predictive performance on the spatio-temporal patterns and changes in the region under study. The hybrid model still has some limitations and uncertainties in our approach. First, the representation of ground-based in situ measurements is not optimal at the pixel scale. Also, due to the spatial scale of 100 m resolution, the conversion to the site scale still does not accurately judge the performance of the model due to the influence of surrounding pixels. Because the HRB is characterized by a variety of geographic and climatic features, the flux sites used in this study were few and unevenly distributed, although the hybrid model performed well in the middle and upper reaches of the study area. The sparse distribution of flux sites weakens the spatial representativeness of the observation samples [74]. In addition, the scale mismatch between satellite and in situ observations may lead to 5% to 25% uncertainty in the spatio-temporal estimates [5]. Ground-based observations are assumed to be real targets in the hybrid modeling process; however, there are no real ET observations for satellite image element scales with resolutions greater than 100 m, especially over heterogeneous surfaces [88]. In addition, flux observation methods are unable to measure large eddies, leading to energy imbalances, which can lead to errors of ~5–20% [89,90]. We collected data from the main flux observation sites in the upper-middle region and selected sites that generally fulfilled at least nine months of continuous observation values and were corrected to still have at least six months of usable data. We collected data from 2011 to 2016 for our experimental study. All selected ground-based observation sites were corrected for quality control and optimization related to energy imbalance. The training and validation sites took into account the participation of different surface vegetation cover types in limited areas and conditions to ensure that the ground observations were representative.

Second, there is some uncertainty about the synergy between physical laws and machine learning. Although many algorithms for ET have been developed since the 1960s; there is still no consensus on the best way to parameterize ET at the regional scale [41]. Process-based physical modeling algorithms show great variability and uncertainty in different terrestrial biome ecosystems [91]. In this study, we did not have ground-based observations of surface conductance to correct our parameter variables. The target variables used for model training were inferred from energy balance equations and process-based physical models. Uncertainties in the inferred results that may be caused by the physical model can directly affect the hybrid model’s prediction accuracy. Similarly, models based on pure machine learning are extremely likely to have better performance locally and fail to scale to other regions, thus falling into local optimization, which leads to overfitting [92,93]. In our study, the hybrid model performs a grid search and cross-validation step in the modeling of machine learning to avoid transition fitting.

For data sources, a wide range of data are involved due to the large number of variables required to compute ET. Complex data sources contain large uncertainties. Mixedness models incorporate a variety of variables that are not identical, including meteorological and satellite data. According to research, comparisons between regional reanalysis data and observations at the micrometeorological scale show large BIASes [94]. The analysis showed that the largest differences between the input variables for ET came from the radiative dataset [95].

Finally, we hope that future research will explore some of the potential advantages and benefits beyond accuracy. Insufficient evaluation of satellite remote sensing models against ground-based observations can affect researchers’ understanding of regional ET. Due to concurrent scale effects, studies lack real values at the satellite image element scale to test remote-sensing-based ET estimates. In the face of this problem, improving the accuracy of ground-based observations may not be the only goal in the development of hybrid models. In addition to the validation of in situ-observation-based values, future studies should conduct more comprehensive assessments, such as the performance of energy balance and water balance, the capture of irrigated agriculture and intensive precipitation, as well as feedback to extreme weather events, that can provide future perspectives and avenues for exploring the structure and potential benefits of hybrid models.

6. Conclusions

This thesis provides a hybrid model for estimating ET coupled with physical constraints and machine learning methods for estimating surface ET in the middle and upper reaches of the HRB. The model has the advantages of being based on physical mechanisms (theoretical foundation, interpretability) and combining machine learning methods (data adaptability). The results of this study are summarized below:

Evaluation metrics used in this study indicate that the hybrid model enhances the performance of the process-based ET algorithm. In areas where the ground cover is agricultural land, the hybrid model improves the $R^2$ between daily ET predictions and ground observations from the process-based ET algorithm to 0.9 and reduces the BIAS to 0.2, while increasing the KGE to 0.8 and above. Under other surface vegetation types such as grassland and forest, the KGE of the hybrid model stabilizes at around 0.8, greatly improving the generalization ability of the model.

Second, although the accuracy of data-oriented pure machine learning models is higher in some areas, they cannot replace physical models, especially for areas with heterogeneous subsurface and sparse and unstudied vegetation. Pure machine learning is strongly influenced by sample data. Hybrid models are more advantageous compared to traditional models in sparsely vegetated areas. The hybrid model outperformed the traditional process-based ET algorithm for the grassland site in this study. For sparsely vegetated areas, the more stable the coupled physical mechanisms are, the better the generalization ability of the hybrid model is.

Overall, ET estimates from all models showed consistent spatial patterns and distributions across the upper and middle HRB region with changes in elevation and local microclimate, as well as surface plant functional types. All ET models showed spatial variations in estimates consistent with variations in the distribution of surface plant functional types. The multi-year average annual ET of the six models ranged from 423.78 to 536.83 mm/year. Compared with the other models, the hybrid model (RF_Gs) was better accordingly for the distribution of surface vegetation types. Mainly, the hybrid model can better capture the transition and demarcation parts of the southern grassland and forest regions. The hybrid model demonstrates a notably reduced tendency for overestimating ET in wetland and farmland regions compared to the standalone machine learning models.

Comparison with the other products showed that the ET products based on the DTD and ETMonitor were overall better than the other three ET products in terms of accuracy; among the remaining products, GLASS and GLEAM were the next most accurate, and the MODIS product was the worst. The hybrid model performed best in this study.

In this study, a hybrid model is constructed by coupling physical constraints and machine learning models, which provides a machine learning-based approach to the accurate quantification of surface conductance and, at the same time, can provide a parameter correction scheme for the prediction of ET by the dual-source models of Peman-Monteith, based on the large-leaf theory, and SEABL, based on the surface energy balance model. The hybrid model retains the physical framework while combining the ability of the machine learning model to extract features and construct nonlinear fits. This is a promising path. Although the hybrid model performs better than other models in terms of spatial and temporal patterns and ET variability, there are still some limitations and uncertainties in the methodology of this study. First, the values measured in the field are not the best comparative data at the image element scale. In addition, due to the complex geographic and climatic characteristics of the HRB, the EC flux sites used in this study are limited and unevenly distributed. The ET process is influenced by the local microclimate and the feedback mechanism is complex, and not all relevant parameters can be easily derived from a certain principle or a certain formula. In future studies, for the multidimensional data sources involved in the model (e.g., radiation, vegetation, temperature, aerodynamics, water vapor pressure, or spatio-temporal autocorrelation information), further research is needed to better incorporate the machine learning approach so that better dynamic dependencies can be constructed between the model parameters and the fusion framework.

Author Contributions

Conceptualization, X.L. and F.X.; methodology, X.L. and F.X.; validation, F.X.; formal analysis, X.L. and F.X.; investigation, X.L., F.X., and J.D.; resources, X.L., T.X., L.S., Y.M. and Z.X.; data curation, L.S., Y.M. and Z.X.; writing—original draft preparation, X.L.; writing—review and editing, X.L., F.X., J.D., J.W. (Jinjie Wang), T.X., Z.P., Z.X., X.H., Y.Z. and Z.X.; visualization, X.L., Z.L., D.W., J.W. (Jiaxing Wei), X.H. and Y.Z.; supervision, X.L.; project administration, X.L.; funding acquisition, X.L. and J.D. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data presented in this study were provided by the National Tibetan Plateau Data Center (https://data.tpdc.ac.cn) on 23 March 2024.

Acknowledgments

The authors would like to thank all the scientists, engineers, and students who participated in WATER and HiWATER field campaigns. We appreciate all reviewers and editors for their comments on this paper.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. The framework chart of predicting the ET based on the hybrid model coupled with physical constraints and machine learning. (SE_PM is the abbreviation of the semi-empirical model proposed by Wang et al. [30]).

Figure 2. The spatial location of the HRB and the location of stations used in this study.

Figure 3. Comparison of inferred surface conductance values with those predicted at independent sites using methods with (a) and without (b) optimization parameters.

Figure 4. Heat map of pure machine learning model features (where correlations less than 0.05 are drawn ×).

Figure 5. The comparison between the ET values estimated by six models and the ET values observed on the ground at Daman Station (a), Linze Station (b), and Yingke Station (c).

Figure 6. The comparison between the ET values estimated by six models and the ET values observed on the ground at Zhangye Station (a), Arou Station (b), and Dashalong Station (c).

Figure 7. Comparison of time series changes in ground observation ET, hybrid model estimation ET, pure machine learning estimation ET, and physical model estimation ET at six ground observation sites.

Figure 8. Boxplots of all models at daily (upper), monthly (middle), and annual (under) scales for BIAS (left) and RMSE (right).

Figure 9. Boxplots of all models at daily (upper), monthly (middle), and annual (under) scales for [Forumla omitted. See PDF.] (left) and KGE (right).

Figure 10. Taylor analysis error plots of all models at all sites.

Figure 11. Comparison of hybrid model, pure physical model, and pure machine learning model KGE under different plant functional types (where Daman represents farmland, Arou represents forest, and Dashalong represents grassland).

Figure 12. Comparison of hybrid ET models and commonly known ET products for different plant functional types. All: all plant functional types.

Figure 13. Average ET predicted by six models in the middle and upper reaches of the HRB (2011–2016).

Figure 14. Analytical plot of the main factor contributions of machine learning models for predicting surface conductance in hybrid models.

Table 1

Specifications of the stations in the HRB used in this study.

Region	No.	Station	Longitude	Latitude	Elevation	Landscape	Time of Period
Upstream	1	Arou	100.4643	38.0473	3033	grassland	January 2013–December 2016
	2	Dashalong	98.9406	38.8399	3739	grassland	August 2013–December 2016
Midstream	3	Daman	100.3722	38.8555	1556	cropland	June 2012–December 2016
	4	Zhangye	100.4464	38.9751	1460	wetland	June 2012–December 2016
	5	Yingke	100.4103	38.8571	1519	cropland	January 2011–December 2011
	6	Linze	100.1408	39.3281	1252	cropland	January 2013–December 2014

Table 2

Specific characteristics of the coarse remote sensing ET products used in this study.

Product Category	Product Name	Algorithm	Time Period	Time Resolution	Spatial Resolution	Unit	References
Surface evaporation products based on energy balance	GLEAM	P-T	1980–2017	Daily	0.25°	$\mathrm{m}\mathrm{m}/\mathrm{d}\mathrm{a}\mathrm{y}$	[79]
	DTD	TSEB	2012–2016	Daily	1 km	$\mathrm{m}\mathrm{m}/\mathrm{d}\mathrm{a}\mathrm{y}$	[80]
Surface evapotranspiration products based on vegetation physiological and ecological characteristics	MODIS	P-M	2000–2017	8 day	500 m	$\mathrm{J}/{\mathrm{m}}^2\mathrm{d}\mathrm{a}\mathrm{y}$	[6]
	ETMonitor	S-W	2009–2016	Daily	1 km	$\mathrm{m}\mathrm{m}/\mathrm{d}\mathrm{a}\mathrm{y}$	[9]
Integrated products	GLASS	BMA	2012–2016	8 day	1 km	$\mathrm{W}/{\mathrm{m}}^2$	[8]

References

1. Kool, D.; Agam, N.; Lazarovitch, N.; Heitman, J.L.; Sauer, T.J.; Ben-Gal, A. A review of approaches for evapotranspiration partitioning. Agric. For. Meteorol.; 2014; 184, pp. 56-70. [DOI: https://dx.doi.org/10.1016/j.agrformet.2013.09.003]

2. Trenberth, K.E.; Fasullo, J.T.; Kiehl, J. Earth’s Global Energy Budget. Bull. Am. Meteorol. Soc.; 2009; 90, pp. 311-324. [DOI: https://dx.doi.org/10.1175/2008BAMS2634.1]

3. Jung, M.; Reichstein, M.; Ciais, P.; Seneviratne, S.I.; Sheffield, J.; Goulden, M.L.; Bonan, G.; Cescatti, A.; Chen, J.; de Jeu, R. et al. Recent decline in the global land evapotranspiration trend due to limited moisture supply. Nature; 2010; 467, pp. 951-954. [DOI: https://dx.doi.org/10.1038/nature09396] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/20935626]

4. Fisher, J.B.; Melton, F.; Middleton, E.; Hain, C.; Anderson, M.; Allen, R.; McCabe, M.F.; Hook, S.; Baldocchi, D.; Townsend, P.A. et al. The future of evapotranspiration: Global requirements for ecosystem functioning, carbon and climate feedbacks, agricultural management, and water resources. Water Resour. Res.; 2017; 53, pp. 2618-2626. [DOI: https://dx.doi.org/10.1002/2016WR020175]

5. Li, Z.-L.; Tang, R.; Wan, Z.; Bi, Y.; Zhou, C.; Tang, B.; Yan, G.; Zhang, X. A Review of Current Methodologies for Regional Evapotranspiration Estimation from Remotely Sensed Data. Sensors; 2009; 9, pp. 3801-3853. [DOI: https://dx.doi.org/10.3390/s90503801] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/22412339]

6. Mu, Q.; Zhao, M.; Running, S.W. Improvements to a MODIS global terrestrial evapotranspiration algorithm. Remote Sens. Environ.; 2011; 115, pp. 1781-1800. [DOI: https://dx.doi.org/10.1016/j.rse.2011.02.019]

7. Miralles, D.G.; Holmes, T.R.H.; De Jeu, R.A.M.; Gash, J.H.; Meesters, A.G.C.A.; Dolman, A.J. Global land-surface evaporation estimated from satellite-based observations. Hydrol. Earth Syst. Sci.; 2011; 15, pp. 453-469. [DOI: https://dx.doi.org/10.5194/hess-15-453-2011]

8. Yao, Y.; Liang, S.; Li, X.; Hong, Y.; Fisher, J.B.; Zhang, N.; Chen, J.; Cheng, J.; Zhao, S.; Zhang, X. et al. Bayesian multimodel estimation of global terrestrial latent heat flux from eddy covariance, meteorological, and satellite observations. J. Geophys. Res. Atmos.; 2014; 119, pp. 4521-4545. [DOI: https://dx.doi.org/10.1002/2013JD020864]

9. Hu, G.; Jia, L. Monitoring of Evapotranspiration in a Semi-Arid Inland River Basin by Combining Microwave and Optical Remote Sensing Observations. Remote Sens.; 2015; 7, pp. 3056-3087. [DOI: https://dx.doi.org/10.3390/rs70303056]

10. Jiang, C.; Ryu, Y. Multi-scale evaluation of global gross primary productivity and evapotranspiration products derived from Breathing Earth System Simulator (BESS). Remote Sens. Environ.; 2016; 186, pp. 528-547. [DOI: https://dx.doi.org/10.1016/j.rse.2016.08.030]

11. van der Kooij, S.; Zwarteveen, M.; Boesveld, H.; Kuper, M. The efficiency of drip irrigation unpacked. Agric. Water Manag.; 2013; 123, pp. 103-110. [DOI: https://dx.doi.org/10.1016/j.agwat.2013.03.014]

12. Battude, M.; Al Bitar, A.; Brut, A.; Tallec, T.; Huc, M.; Cros, J.; Weber, J.-J.; Lhuissier, L.; Simonneaux, V.; Demarez, V. Modeling water needs and total irrigation depths of maize crop in the south west of France using high spatial and temporal resolution satellite imagery. Agric. Water Manag.; 2017; 189, pp. 123-136. [DOI: https://dx.doi.org/10.1016/j.agwat.2017.04.018]

13. Ahadi, R.; Samani, Z.; Skaggs, R. Evaluating on-farm irrigation efficiency across the watershed: A case study of New Mexico’s Lower Rio Grande Basin. Agric. Water Manag.; 2013; 124, pp. 52-57. [DOI: https://dx.doi.org/10.1016/j.agwat.2013.03.010]

14. Bastiaanssen, W.G.M.; Noordman, E.J.M.; Pelgrum, H.; Davids, G.; Thoreson, B.P.; Allen, R.G. SEBAL Model with Remotely Sensed Data to Improve Water-Resources Management under Actual Field Conditions. J. Irrig. Drain. Eng.; 2005; 131, pp. 85-93. [DOI: https://dx.doi.org/10.1061/(ASCE)0733-9437(2005)131:1(85)]

15. Yang, Y.; Shang, S.; Jiang, L. Remote sensing temporal and spatial patterns of evapotranspiration and the responses to water management in a large irrigation district of North China. Agric. For. Meteorol.; 2012; 164, pp. 112-122. [DOI: https://dx.doi.org/10.1016/j.agrformet.2012.05.011]

16. de Fraiture, C.; Wichelns, D. Satisfying future water demands for agriculture. Agric. Water Manag.; 2010; 97, pp. 502-511. [DOI: https://dx.doi.org/10.1016/j.agwat.2009.08.008]

17. Fereres, E.; Soriano, M.A. Deficit irrigation for reducing agricultural water use. J. Exp. Bot.; 2007; 58, pp. 147-159. [DOI: https://dx.doi.org/10.1093/jxb/erl165] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/17088360]

18. Allen, R.G.; Tasumi, M.; Morse, A.; Trezza, R. A Landsat-based energy balance and evapotranspiration model in Western US water rights regulation and planning. Irrig. Drain. Syst.; 2005; 19, pp. 251-268. [DOI: https://dx.doi.org/10.1007/s10795-005-5187-z]

19. Ma, Y.; Liu, S.; Song, L.; Xu, Z.; Liu, Y.; Xu, T.; Zhu, Z. Estimation of daily evapotranspiration and irrigation water efficiency at a Landsat-like scale for an arid irrigation area using multi-source remote sensing data. Remote Sens. Environ.; 2018; 216, pp. 715-734. [DOI: https://dx.doi.org/10.1016/j.rse.2018.07.019]

20. Cammalleri, C.; Anderson, M.C.; Gao, F.; Hain, C.R.; Kustas, W.P. A data fusion approach for mapping daily evapotranspiration at field scale. Water Resour. Res.; 2013; 49, pp. 4672-4686. [DOI: https://dx.doi.org/10.1002/wrcr.20349]

21. Cammalleri, C.; Anderson, M.C.; Kustas, W.P. Upscaling of evapotranspiration fluxes from instantaneous to daytime scales for thermal remote sensing applications. Hydrol. Earth Syst. Sci.; 2014; 18, pp. 1885-1894. [DOI: https://dx.doi.org/10.5194/hess-18-1885-2014]

22. Cammalleri, C.; Anderson, M.C.; Gao, F.; Hain, C.R.; Kustas, W.P. Mapping daily evapotranspiration at field scales over rainfed and irrigated agricultural areas using remote sensing data fusion. Agric. For. Meteorol.; 2014; 186, pp. 1-11. [DOI: https://dx.doi.org/10.1016/j.agrformet.2013.11.001]

23. Eswar, R.; Sekhar, M.; Bhattacharya, B.K.; Bandyopadhyay, S. Spatial Disaggregation of Latent Heat Flux Using Contextual Models over India. Remote Sens.; 2017; 9, 949. [DOI: https://dx.doi.org/10.3390/rs9090949]

24. Chen, J.M.; Liu, J. Evolution of evapotranspiration models using thermal and shortwave remote sensing data. Remote Sens. Environ.; 2020; 237, 111594. [DOI: https://dx.doi.org/10.1016/j.rse.2019.111594]

25. Monteith, J.L. Evaporation and environment. Symp. Soc. Exp. Biol.; 1965; 19, pp. 205-234.

26. Priestley, C.H.B.; Taylor, R.J. On the Assessment of Surface Heat Flux and Evaporation Using Large-Scale Parameters. Mon. Weather. Rev.; 1972; 100, pp. 81-92. [DOI: https://dx.doi.org/10.1175/1520-0493(1972)100<0081:OTAOSH>2.3.CO;2]

27. Su, Z. The Surface Energy Balance System (SEBS) for estimation of turbulent heat fluxes. Hydrol. Earth Syst. Sci.; 2002; 6, pp. 85-100. [DOI: https://dx.doi.org/10.5194/hess-6-85-2002]

28. Bastiaanssen, W.G.M.; Menenti, M.; Feddes, R.A.; Holtslag, A.A.M. A remote sensing surface energy balance algorithm for land (SEBAL). 1. Formulation. J. Hydrol.; 1998; 212–213, pp. 198-212. [DOI: https://dx.doi.org/10.1016/S0022-1694(98)00253-4]

29. Norman, J.M.; Kustas, W.P.; Humes, K.S. Source approach for estimating soil and vegetation energy fluxes in observations of directional radiometric surface temperature. Agric. For. Meteorol.; 1995; 77, pp. 263-293. [DOI: https://dx.doi.org/10.1016/0168-1923(95)02265-Y]

30. Wang, K.; Dickinson, R.E.; Wild, M.; Liang, S. Evidence for decadal variation in global terrestrial evapotranspiration between 1982 and 2002: 1. Model development. J. Geophys. Res. Atmos; 2010; 115, [DOI: https://dx.doi.org/10.1029/2009JD013671]

31. Wang, K.; Dickinson, R.E.; Wild, M.; Liang, S. Evidence for decadal variation in global terrestrial evapotranspiration between 1982 and 2002: 2. Results. J. Geophys. Res. Atmos; 2010; 115, [DOI: https://dx.doi.org/10.1029/2010JD013847]

32. Chen, X.; Massman, W.J.; Su, Z. A Column Canopy-Air Turbulent Diffusion Method for Different Canopy Structures. J. Geophys. Res. Atmos.; 2019; 124, pp. 488-506. [DOI: https://dx.doi.org/10.1029/2018JD028883]

33. Fisher, J.B.; Lee, B.; Purdy, A.J.; Halverson, G.H.; Dohlen, M.B.; Cawse-Nicholson, K.; Wang, A.; Anderson, R.G.; Aragon, B.; Arain, M.A. et al. ECOSTRESS: NASA’s Next Generation Mission to Measure Evapotranspiration From the International Space Station. Water Resour. Res.; 2020; 56, e2019WR026058. [DOI: https://dx.doi.org/10.1029/2019WR026058]

34. Shengzhe, L. Study on Cognitive Optical Network Structure and Self-optimization with the Application of Artificial Intelligence Technology. Mod. Electron. Technol.; 2020; 3, [DOI: https://dx.doi.org/10.26549/met.v4i1.3100]

35. Martens, B.; Miralles, D.G.; Lievens, H.; van der Schalie, R.; de Jeu, R.A.M.; Fernández-Prieto, D.; Beck, H.E.; Dorigo, W.A.; Verhoest, N.E.C. GLEAM v3: Satellite-based land evaporation and root-zone soil moisture. Geosci. Model. Dev.; 2017; 10, pp. 1903-1925. [DOI: https://dx.doi.org/10.5194/gmd-10-1903-2017]

36. Allen, R.G.; Pereira, L.S.; Howell, T.A.; Jensen, M.E. Evapotranspiration information reporting: I. Factors governing measurement accuracy. Agric. Water Manag.; 2011; 98, pp. 899-920. [DOI: https://dx.doi.org/10.1016/j.agwat.2010.12.015]

37. Smith, W.K.; Dannenberg, M.P.; Yan, D.; Herrmann, S.; Barnes, M.L.; Barron-Gafford, G.A.; Biederman, J.A.; Ferrenberg, S.; Fox, A.M.; Hudson, A. et al. Remote sensing of dryland ecosystem structure and function: Progress, challenges, and opportunities. Remote Sens. Environ.; 2019; 233, 111401. [DOI: https://dx.doi.org/10.1016/j.rse.2019.111401]

38. Han, C.; Ma, Y.; Wang, B.; Zhong, L.; Ma, W.; Chen, X.; Su, Z. Long-term variations in actual evapotranspiration over the Tibetan Plateau. Earth Syst. Sci. Data; 2021; 13, pp. 3513-3524. [DOI: https://dx.doi.org/10.5194/essd-13-3513-2021]

39. Ma, N.; Zhang, Y. Increasing Tibetan Plateau terrestrial evapotranspiration primarily driven by precipitation. Agric. For. Meteorol.; 2022; 317, 108887. [DOI: https://dx.doi.org/10.1016/j.agrformet.2022.108887]

40. Wang, W.; Li, J.; Yu, Z.; Ding, Y.; Xing, W.; Lu, W. Satellite retrieval of actual evapotranspiration in the Tibetan Plateau: Components partitioning, multidecadal trends and dominated factors identifying. J. Hydrol.; 2018; 559, pp. 471-485. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2018.02.065]

41. Yuan, L.; Ma, Y.; Chen, X.; Wang, Y.; Li, Z. An Enhanced MOD16 Evapotranspiration Model for the Tibetan Plateau During the Unfrozen Season. J. Geophys. Res. Atmos.; 2021; 126, e2020JD032787. [DOI: https://dx.doi.org/10.1029/2020JD032787]

42. Bodesheim, P.; Jung, M.; Gans, F.; Mahecha, M.D.; Reichstein, M. Upscaled diurnal cycles of land–atmosphere fluxes: A new global half-hourly data product. Earth Syst. Sci. Data; 2018; 10, pp. 1327-1365. [DOI: https://dx.doi.org/10.5194/essd-10-1327-2018]

43. Jung, M.; Koirala, S.; Weber, U.; Ichii, K.; Gans, F.; Camps-Valls, G.; Papale, D.; Schwalm, C.; Tramontana, G.; Reichstein, M. The FLUXCOM ensemble of global land-atmosphere energy fluxes. Sci. Data; 2019; 6, pp. 1-14. [DOI: https://dx.doi.org/10.1038/s41597-019-0076-8] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/31133670]

44. Li, X.; Zhang, S.; Liu, X.; Wang, X.; Zhou, A.; Liu, P. Important role of MAMs in bifurcation and coherence resonance of calcium oscillations. Chaos Solitons Fractals; 2018; 106, pp. 131-140. [DOI: https://dx.doi.org/10.1016/j.chaos.2017.11.018]

45. Lu, X.; Zhuang, Q. Evaluating evapotranspiration and water-use efficiency of terrestrial ecosystems in the conterminous United States using MODIS and AmeriFlux data. Remote Sens. Environ.; 2010; 114, pp. 1924-1939. [DOI: https://dx.doi.org/10.1016/j.rse.2010.04.001]

46. Trajkovic, S.; Todorovic, B.; Stankovic, M. Forecasting of reference evapotranspiration by artificial neural networks. J. Irrig. Drain. Eng.; 2003; 129, pp. 454-457. [DOI: https://dx.doi.org/10.1061/(ASCE)0733-9437(2003)129:6(454)]

47. Reichstein, M.; Camps-Valls, G.; Stevens, B.; Jung, M.; Denzler, J.; Carvalhais, N.; Prabhat,. Deep learning and process understanding for data-driven Earth system science. Nature; 2019; 566, pp. 195-204. [DOI: https://dx.doi.org/10.1038/s41586-019-0912-1] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30760912]

48. Shang, K.; Yao, Y.; Liang, S.; Zhang, Y.; Fisher, J.B.; Chen, J.; Liu, S.; Xu, Z.; Zhang, Y.; Jia, K. et al. DNN-MET: A deep neural networks method to integrate satellite-derived evapotranspiration products, eddy covariance observations and ancillary information. Agric. For. Meteorol.; 2021; 308–309, 108582. [DOI: https://dx.doi.org/10.1016/j.agrformet.2021.108582]

49. Zhao, W.L.; Gentine, P.; Reichstein, M.; Zhang, Y.; Zhou, S.; Wen, Y.; Lin, C.; Li, X.; Qiu, G.Y. Physics-Constrained Machine Learning of Evapotranspiration. Geophys. Res. Lett.; 2019; 46, pp. 14496-14507. [DOI: https://dx.doi.org/10.1029/2019GL085291]

50. Karpatne, A.; Atluri, G.; Faghmous, J.H.; Steinbach, M.; Banerjee, A.; Ganguly, A.; Shekhar, S.; Samatova, N.; Kumar, V. Theory-guided data science: A new paradigm for scientific discovery from data. IEEE Trans. Knowl. Data Eng.; 2017; 29, pp. 2318-2331. [DOI: https://dx.doi.org/10.1109/TKDE.2017.2720168]

51. Koppa, A.; Rains, D.; Hulsman, P.; Poyatos, R.; Miralles, D.G. A deep learning-based hybrid model of global terrestrial evaporation. Nat. Commun.; 2022; 13, pp. 1-11. [DOI: https://dx.doi.org/10.1038/s41467-022-29543-7] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35395845]

52. Tan, Z.-H.; Zhao, J.-F.; Wang, G.-Z.; Chen, M.-P.; Yang, L.-Y.; He, C.-S.; Restrepo-Coupe, N.; Peng, S.-S.; Liu, X.-Y.; da Rocha, H.R. et al. Surface conductance for evapotranspiration of tropical forests: Calculations, variations, and controls. Agric. For. Meteorol.; 2019; 275, pp. 317-328. [DOI: https://dx.doi.org/10.1016/j.agrformet.2019.06.006]

53. Meinzer, F.C.; Goldstein, G.; Holbrook, N.M.; Jackson, P.; Cavelier, J. Stomatal and environmental control of transpiration in a lowland tropical forest tree. Plant Cell Environ.; 1993; 16, pp. 429-436. [DOI: https://dx.doi.org/10.1111/j.1365-3040.1993.tb00889.x]

54. Condit, R.; Pitman, N.; Leigh, E.G.; Chave, J.; Terborgh, J.; Foster, R.B.; Núñez, P.; Aguilar, S.; Valencia, R.; Villa, G. et al. Beta-Diversity in Tropical Forest Trees. Science; 2002; 295, pp. 666-669. [DOI: https://dx.doi.org/10.1126/science.1066854] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/11809969]

55. Roberts, J.; Osvaldo, M.R.C.; De Aguiar, L.F. Stomatal and Boundary-Layer Conductances in an Amazonian terra Firme Rain Forest. J. Appl. Ecol.; 1990; 27, pp. 336-353. [DOI: https://dx.doi.org/10.2307/2403590]

56. Roberts, J.; Cabral, O.M.R.; Fisch, G.; Molion, L.C.B.; Moore, C.J.; Shuttleworth, W.J. Transpiration from an Amazonian rainforest calculated from stomatal conductance measurements. Agric. For. Meteorol.; 1993; 65, pp. 175-196. [DOI: https://dx.doi.org/10.1016/0168-1923(93)90003-Z]

57. Ma, N.; Zhang, Y.; Xu, C.Y.; Szilagyi, J. Modeling actual evapotranspiration with routine meteorological variables in the data-scarce region of the Tibetan Plateau: Comparisons and implications. J. Geophys. Res. Biogeosci; 2015; 120, pp. 1638-1657. [DOI: https://dx.doi.org/10.1002/2015JG003006]

58. Kochendorfer, J.; Castillo, E.G.; Haas, E.; Oechel, W.C.; Paw U, K.T. Net ecosystem exchange, evapotranspiration and canopy conductance in a riparian forest. Agric. For. Meteorol.; 2011; 151, pp. 544-553. [DOI: https://dx.doi.org/10.1016/j.agrformet.2010.12.012]

59. Nishida, K.; Nemani, R.R.; Running, S.W.; Glassy, J.M. An operational remote sensing algorithm of land surface evaporation. J. Geophys. Res. Atmos.; 2003; 108, D9. [DOI: https://dx.doi.org/10.1029/2002JD002062]

60. Zhao, W.; Liu, B.; Chang, X.; Yang, Q.; Yang, Y.; Liu, Z.; Cleverly, J.; Eamus, D. Evapotranspiration partitioning, stomatal conductance, and components of the water balance: A special case of a desert ecosystem in China. J. Hydrol.; 2016; 538, pp. 374-386. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2016.04.042]

61. Cheng, G.; Li, X.; Zhao, W.; Xu, Z.; Feng, Q.; Xiao, S.; Xiao, H. Integrated study of the water–ecosystem–economy in the Heihe River Basin. Natl. Sci. Rev.; 2014; 1, pp. 413-428. [DOI: https://dx.doi.org/10.1093/nsr/nwu017]

62. Li, X.; Cheng, G.; Ge, Y.; Li, H.; Han, F.; Hu, X.; Tian, W.; Tian, Y.; Pan, X.; Nian, Y. et al. Hydrological Cycle in the Heihe River Basin and Its Implication for Water Resource Management in Endorheic Basins. J. Geophys. Res. Atmos.; 2018; 123, pp. 890-914. [DOI: https://dx.doi.org/10.1002/2017JD027889]

63. Parlange, M.B.; Eichinger, W.E.; Albertson, J.D. Regional scale evaporation and the atmospheric boundary layer. Rev. Geophys.; 1995; 33, pp. 99-124. [DOI: https://dx.doi.org/10.1029/94RG03112]

64. Shuttleworth, W.J. Putting the “vap” into evaporation. Hydrol. Earth Syst. Sci.; 2007; 11, pp. 210-244. [DOI: https://dx.doi.org/10.5194/hess-11-210-2007]

65. Yang, H.; Yang, D.; Lei, Z. Seasonal variability of the complementary relationship in the Asian monsoon region. Hydrol. Process.; 2013; 27, pp. 2736-2741. [DOI: https://dx.doi.org/10.1002/hyp.9400]

66. Breiman, L. Random Forests. Mach. Learn.; 2001; 45, pp. 5-32. [DOI: https://dx.doi.org/10.1023/A:1010933404324]

67. Matin, S.S.; Chelgani, S.C. Estimation of coal gross calorific value based on various analyses by random forest method. Fuel; 2016; 177, pp. 274-278. [DOI: https://dx.doi.org/10.1016/j.fuel.2016.03.031]

68. Feng, Y.; Cui, N.; Gong, D.; Zhang, Q.; Zhao, L. Evaluation of random forests and generalized regression neural networks for daily reference evapotranspiration modelling. Agric. Water Manag.; 2017; 193, pp. 163-173. [DOI: https://dx.doi.org/10.1016/j.agwat.2017.08.003]

69. Fan, J.; Wu, L.; Zheng, J.; Zhang, F. Medium-range forecasting of daily reference evapotranspiration across China using numerical weather prediction outputs downscaled by extreme gradient boosting. J. Hydrol.; 2021; 601, 126664. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2021.126664]

70. Yamaç, S.S.; Todorovic, M. Estimation of daily potato crop evapotranspiration using three different machine learning algorithms and four scenarios of available meteorological data. Agric. Water Manag.; 2020; 228, 105875. [DOI: https://dx.doi.org/10.1016/j.agwat.2019.105875]

71. Wu, B.; Liu, S.; Zhu, W.; Yan, N.; Xing, Q.; Tan, S. An Improved Approach for Estimating Daily Net Radiation over the Heihe River Basin. Sensors; 2017; 17, 86. [DOI: https://dx.doi.org/10.3390/s17010086] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/28054976]

72. Lima, C.E.S.d.; Costa, V.S.d.O.; Galvíncio, J.D.; Silva, R.M.d.; Santos, C.A.G. Assessment of automated evapotranspiration estimates obtained using the GP-SEBAL algorithm for dry forest vegetation (Caatinga) and agricultural areas in the Brazilian semiarid region. Agric. Water Manag.; 2021; 250, 106863. [DOI: https://dx.doi.org/10.1016/j.agwat.2021.106863]

73. Yao, Y.; Liang, S.; Cheng, J.; Liu, S.; Fisher, J.B.; Zhang, X.; Jia, K.; Zhao, X.; Qin, Q.; Zhao, B. et al. MODIS-driven estimation of terrestrial latent heat flux in China based on a modified Priestley–Taylor algorithm. Agric. For. Meteorol.; 2013; 171–172, pp. 187-202. [DOI: https://dx.doi.org/10.1016/j.agrformet.2012.11.016]

74. Liu, S.; Li, X.; Xu, Z.; Che, T.; Xiao, Q.; Ma, M.; Liu, Q.; Jin, R.; Guo, J.; Wang, L. et al. The Heihe Integrated Observatory Network: A Basin-Scale Land Surface Processes Observatory in China. Vadose Zone J.; 2018; 17, 180072. [DOI: https://dx.doi.org/10.2136/vzj2018.04.0072]

75. Li, X.; Cheng, G.; Liu, S.; Xiao, Q.; Ma, M.; Jin, R.; Che, T.; Liu, Q.; Wang, W.; Qi, Y. et al. Heihe Watershed Allied Telemetry Experimental Research (HiWATER): Scientific Objectives and Experimental Design. Bull. Am. Meteorol. Soc.; 2013; 94, pp. 1145-1160. [DOI: https://dx.doi.org/10.1175/BAMS-D-12-00154.1]

76. Xu, Z.; Liu, S.; Li, X.; Shi, S.; Wang, J.; Zhu, Z.; Xu, T.; Wang, W.; Ma, M. Intercomparison of surface energy flux measurement systems used during the HiWATER-MUSOEXE. J. Geophys. Res. Atmos.; 2013; 118, 23. [DOI: https://dx.doi.org/10.1002/2013JD020260]

77. Liu, S.; Xu, Z.; Song, L.; Zhao, Q.; Ge, Y.; Xu, T.; Ma, Y.; Zhu, Z.; Jia, Z.; Zhang, F. Upscaling evapotranspiration measurements from multi-site to the satellite pixel scale over heterogeneous land surfaces. Agric. For. Meteorol.; 2016; 230–231, pp. 97-113. [DOI: https://dx.doi.org/10.1016/j.agrformet.2016.04.008]

78. Liu, S.; Xu, Z.; Che, T.; Li, X.; Xu, T.; Ren, Z.; Zhang, Y.; Tan, J.; Song, L.; Zhou, J. et al. A dataset of energy, water vapor, and carbon exchange observations in oasis–desert areas from 2012 to 2021 in a typical endorheic basin. Earth Syst. Sci. Data; 2023; 15, pp. 4959-4981. [DOI: https://dx.doi.org/10.5194/essd-15-4959-2023]

79. Miralles, D.G.; De Jeu, R.A.M.; Gash, J.H.; Holmes, T.R.H.; Dolman, A.J. Magnitude and variability of land evaporation and its components at the global scale. Hydrol. Earth Syst. Sci.; 2011; 15, pp. 967-981. [DOI: https://dx.doi.org/10.5194/hess-15-967-2011]

80. Song, L.; Liu, S.; Kustas, W.P.; Nieto, H.; Sun, L.; Xu, Z.; Skaggs, T.H.; Yang, Y.; Ma, M.; Xu, T. et al. Monitoring and validating spatially and temporally continuous daily evaporation and transpiration at river basin scale. Remote Sens. Environ.; 2018; 219, pp. 72-88. [DOI: https://dx.doi.org/10.1016/j.rse.2018.10.002]

81. Gupta, H.V.; Kling, H.; Yilmaz, K.K.; Martinez, G.F. Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling. J. Hydrol.; 2009; 377, pp. 80-91. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2009.08.003]

82. Shang, K.; Yao, Y.; Di, Z.; Jia, K.; Zhang, X.; Fisher, J.B.; Chen, J.; Guo, X.; Yang, J.; Yu, R. et al. Coupling physical constraints with machine learning for satellite-derived evapotranspiration of the Tibetan Plateau. Remote Sens. Environ.; 2023; 289, 113519. [DOI: https://dx.doi.org/10.1016/j.rse.2023.113519]

83. Aas, K.; Jullum, M.; Løland, A. Explaining individual predictions when features are dependent: More accurate approximations to Shapley values. Artif. Intell.; 2021; 298, 103502. [DOI: https://dx.doi.org/10.1016/j.artint.2021.103502]

84. Han, J.; Wang, J.; Zhao, Y.; Wang, Q.; Zhang, B.; Li, H.; Zhai, J. Spatio-temporal variation of potential evapotranspiration and climatic drivers in the Jing-Jin-Ji region, North China. Agric. For. Meteorol.; 2018; 256–257, pp. 75-83. [DOI: https://dx.doi.org/10.1016/j.agrformet.2018.03.002]

85. Marshall, M.; Tu, K.; Andreo, V. On Parameterizing Soil Evaporation in a Direct Remote Sensing Model of ET: PT-JPL. Water Resour. Res.; 2020; 56, e2019WR026290. [DOI: https://dx.doi.org/10.1029/2019WR026290]

86. Gan, R.; Zhang, Y.; Shi, H.; Yang, Y.; Eamus, D.; Cheng, L.; Chiew, F.H.S.; Yu, Q. Use of satellite leaf area index estimating evapotranspiration and gross assimilation for Australian ecosystems. Ecohydrology; 2018; 11, e1974. [DOI: https://dx.doi.org/10.1002/eco.1974]

87. Leuning, R.; Zhang, Y.Q.; Rajaud, A.; Cleugh, H.; Tu, K. A simple surface conductance model to estimate regional evaporation using MODIS leaf area index and the Penman-Monteith equation. Water Resour. Res.; 2008; 44, 10. [DOI: https://dx.doi.org/10.1029/2007WR006562]

88. Kalma, J.D.; McVicar, T.R.; McCabe, M.F. Estimating Land Surface Evaporation: A Review of Methods Using Remotely Sensed Surface Temperature Data. Surv. Geophys.; 2008; 29, pp. 421-469. [DOI: https://dx.doi.org/10.1007/s10712-008-9037-z]

89. Foken, T. The energy balance closure problem: An overview. Ecol. Appl.; 2008; 18, pp. 1351-1367. [DOI: https://dx.doi.org/10.1890/06-0922.1]

90. Mahrt, L. Computing turbulent fluxes near the surface: Needed improvements. Agric. For. Meteorol.; 2010; 150, pp. 501-509. [DOI: https://dx.doi.org/10.1016/j.agrformet.2010.01.015]

91. Jiménez, C.; Martens, B.; Miralles, D.M.; Fisher, J.B.; Beck, H.E.; Fernández-Prieto, D. Exploring the merging of the global land evaporation WACMOS-ET products based on local tower measurements. Hydrol. Earth Syst. Sci.; 2018; 22, pp. 4513-4533. [DOI: https://dx.doi.org/10.5194/hess-22-4513-2018]

92. Ma, L.; Liu, Y.; Zhang, X.; Ye, Y.; Yin, G.; Johnson, B.A. Deep learning in remote sensing applications: A meta-analysis and review. ISPRS J. Photogramm. Remote Sens.; 2019; 152, pp. 166-177. [DOI: https://dx.doi.org/10.1016/j.isprsjprs.2019.04.015]

93. Ball, J.E.; Anderson, D.T.; Chan, C.S. Comprehensive survey of deep learning in remote sensing: Theories, tools, and challenges for the community. J. Appl. Remote Sens.; 2017; 11, 42609. [DOI: https://dx.doi.org/10.1117/1.JRS.11.042609]

94. Rienecker, M.M.; Suarez, M.J.; Gelaro, R.; Todling, R.; Bacmeister, J.; Liu, E.; Bosilovich, M.G.; Schubert, S.D.; Takacs, L.; Kim, G.K. et al. MERRA: NASA’s modern-era retrospective analysis for research and applications. J. Clim.; 2011; 24, pp. 3624-3648. [DOI: https://dx.doi.org/10.1175/JCLI-D-11-00015.1]

95. Badgley, G.; Fisher, J.B.; Jiménez, C.; Tu, K.P.; Vinukollu, R. On uncertainty in global terrestrial evapotranspiration estimates from choice of input forcing datasets. J. Hydrometeorol.; 2015; 16, pp. 1449-1455. [DOI: https://dx.doi.org/10.1175/JHM-D-14-0040.1]