Evaluation of different evapotranspiration products in the middle Yellow River Basin, China

Actual evapotranspiration (ET_a) is one of the important hydrological processes (Brutsaert 2005; Jung et al. 2010). ET_a determines the partitioning of available energy on the land surface into latent and sensible heat flux and can influence regional and global climates (Bonan et al. 1992). Global land ET_a returns approximately 60% of annual land precipitation to the atmosphere (Hetherington & Woodward 2003; Koster et al. 2004; Oki & Kanae 2006), and the ratio is even larger (>95%) in arid and semi-arid regions (Moiwo & Tao 2015). Acceleration or intensification of the hydrological cycle has occurred with global warming, and ET_a has also been affected (Brutsaert & Parlange 1998; Jung et al. 2010). However, the existing observations of ET_a primarily focus on point observations, such as Bowen ratio–energy balances (Shen et al. 2002), flux towers (Fisher et al. 2008), lysimeters (Xu & Chen 2005; Trajkovic 2010), and large aperture scintillometers (Hemakumara et al. 2003), whereas direct observations of ET_a at regional or global scales are still lacking.

Different methods have been used to estimate regional ET_a, such as complementary relationship (Brutsaert 2005; Gao et al. 2012; Ma et al. 2015), the model tree ensemble (MTE) method (Shutov et al. 2006; Jung et al. 2010), the remote sensing-based method (Mu et al. 2011), land surface models (LSMs) (Lakshmi et al. 2011; Yin et al. 2013; Cai et al. 2014), and the reanalysis model (Nakaegwa 2008; Kobayashi et al. 2015). The MTE method, based on eddy covariance, employs a machine-learning algorithm, which is trained by evaporation measurements from the global Flux Network (FLUXNET) database, gridded global meteorological data, and the remotely sensed fraction of absorbed photosynthetically active radiation (FPAR). Largely independent of theoretical model assumptions, the MTE method is considered to perform best in capturing ET_a, and has been widely used for validating and parameterizing other models (Cleugh et al. 2007; Fisher et al. 2008). The remote sensing method provides an unprecedented opportunity to monitor spatiotemporal variability in ET_a in two basic ways: (1) via vegetation information (such as the leaf area index (LAI) or normalized difference vegetation index (NDVI)) and (2) the land surface temperature (LST). The LAI or NDVI, derived from remote sensing data, is used for calculating surface resistance in the Penman–Monteith (PM) algorithm to estimate ET_a (Zhang et al. 2010; Mu et al. 2011). Remote sensed LST is used to calculate latent heat (LE) as the residual of surface net radiation (R_n), soil heat flux (G), and sensible heat flux (H) (i.e., LE = R_n−G−H) (Bastiaanssen 2000; EI Haj EI Tahir et al. 2012; Gokmen et al. 2012). LSMs are excellent in capturing ET_a variation (Xue et al. 2013; Cai et al. 2014) and are widely recommended for estimating the hydrological components for large areas, such as surface energy, water fluxes and their response to near-surface atmospheric forcing (Rodell et al. 2004; Sheffield et al. 2006). Atmospheric reanalysis data, with past observations and state-of-the-art numerical weather prediction systems, can provide various surface hydrological variables, such as soil moisture, precipitation, and ET_a (Hogue et al. 2006; Ebita et al. 2011). Although these methods provide a potential means to estimate ET_a, it is still a challenge to accurately obtain ET_a due to the differences of models and forcing data (Cai et al. 2014; Long et al. 2014). In addition, the impact of different model types and forcing data on ET_a products performance needs to be detected. The traditional water balance method provides a useful tool for addressing these questions. Evaluating ET_a products in a closed basin (Rodell et al. 2004; Li et al. 2014) and analyzing the forcing data and models will benefit the improvement of each model in sensitive regions of climate change.

The middle Yellow River Basin (MYRB), located in the Loess Plateau, endures the most severe soil erosion and water loss in the world (approximately nine to 21 times greater than most major rivers in the world) (Huang et al. 2003). To address the environmental crises and promote ecological restoration, the Chinese government has implemented a number of national strategies since 1999, such as the Natural Forest Conversation Program and the Grain to Green Project (GTGP) (Zhang et al. 2000; FAO 2010), which have significantly impacted hydrological processes (Zhang et al. 2007, 2013). In addition, the MYRB is located in the transition zone of semi-arid and semi-humid regions of China, and is particularly sensitive to climate change (Yang et al. 2010; Wang et al. 2012; Zhu et al. 2015). Various studies have documented significant changes in evapotranspiration in the region (Zhang et al. 2011b; Liu et al. 2012a). Although the spatiotemporal variation of ET_a can be detected by various ET_a datasets, it is unclear which one is optimal for semi-arid and semi-humid regions. Therefore, the evaluation of ET_a products and investigation into the variation of ET_a in the MYRB are essential for understanding hydrological processes, and are of great value for regional sustainable water resource management.

The main objectives of this paper were to evaluate the performance of different ET_a products and to investigate the causes of different performances in the MYRB, which was not only meaningful for the selection of appropriate ET_a products for regional water resource management, but also significant for improving the ET_a models in semi-arid and semi-humid regions. The study objectives are achieved through the following steps: (1) to estimate the reference ET_a based on the traditional water balance method in the MYRB; (2) to evaluate the performance of different ET_a products used in the MYRB against the reference ET_a; (3) to analyze the reasons for differences in ET_a model performance; and (4) to investigate the performance of each model in capturing the temporal variations of ET_a.

Four different categories of ET_a products (Table 1) were evaluated in our study: (1) LSMs, (2) reanalysis, (3) MTE method, and (4) diagnostic models (DMs). Five ET_a products, coming from LSMs in the Global Land Data Assimilation System (GLDAS) with the spatial resolution of 1° and time range from 1979, except for Noah2 from 1948, were employed. AET_Noah1 and AET_Noah2 came from the Noah model, AET_CLM from the community land model (CLM), AET_VIC from the variable infiltration capacity (VIC) and AET_MOS from the mosaic model (MOS) (Rodell et al. 2004; Ferreira et al. 2013). The AET_Noah2 was Land Surface Version 2 (for 1948–2010), forced by the Princeton meteorological forcing data (Sheffield et al. 2006), while the others were Version 1 forced by constantly updated meteorological data (http://disc.sci.gsfc.nasa.gov/hydrogy/data-holdings). The reanalysis of 55 years' ET_a production, with the spatial resolution of 1.25° and time range from 1958 to 2014, was the second Japanese global atmospheric reanalysis project, and improved on many of the deficiencies of the first 25-year reanalysis (namely, AET_JRA55 and AET_JRA25, respectively) (Nakaegwa 2008). The main objective of AET_JRA55 was to produce a comprehensive atmospheric dataset suitable for studying multi-decadal variability and climate change (http://jra.kishou.go.jp/JRA-55/index_en.html). The Jung global land ET_a product, with the spatial resolution of 0.5° and time range from 1982 to 2011, was data-driven based and is compiled using data from 198 global flux monitoring towers (FLUXNET), remote sensing data, meteorological observations, and the MTE (a machine-learning algorithm, data download from http://www.bgc-jena.mpg.de/geodb/projects/Data.php; Jung et al. 2010). ET_a from ZhangKe (http://www.ntsg.umt.edu/project/et; Zhang et al. 2011a) and MODIS (http://www.ntsg.umt.edu/project), with the spatial resolution of 8 km and 1 km, respectively, and time range 1983–2006 and 2000–2013, respectively, were DMs and used remote sensed vegetation information as an important input using the PM method (Zhang et al. 2010; Mu et al. 2011). Additional details on the nine ET_a products employed in our study are shown in Table 1.

Monthly streamflow data from 1960 to 2013 at two gauges, the inlet at Hekou Town and the outlets at Sanmen Gorge, were obtained from the Yellow River Hydrological Bureau. The average annual streamflow was aggregated from monthly data. Monthly meteorological data of 73 national meteorological stations (Figure 1), from the National Climate Center of the China Meteorological Administration (CMA) were used in this study. These data included precipitation, air temperature, wind speed, vapor pressure, and sunshine duration.

The monthly grid precipitation data were averaged from three interpolation results by the inverse distance weighted (IDW) (Qian et al. 2006), spline (New et al. 2000), and kriging (Xu et al. 2006) methods to eliminate the uncertainty from different interpolation methods. Precipitation with a resolution of 0.05° from 1960 to 2014 was used to calculate the reference ET_a with the water balance methods.

Soil moisture at 10 cm, 20 cm, and 50 cm depths was obtained from the agrometeorological stations of the CMA for 1992 to 2013, and the volumetric soil moisture was calculated from the mass percentage using the gravimetric technique (Li et al. 2005). These data were used to evaluate the performance of soil moisture products from four models including MOS, Noah2, CLM, and VIC produced by the Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences (VIC-IGSNRR, http://hydro.igsnrr.ac.cn/resources.html). The four soil moisture products were summed to the 50 cm depth, and compared with the observed data.

The rank-based non-parametric Mann–Kendall statistic test (Kendall 1948) has been commonly used to detect trends due to its robustness for non-normally distributed data, which have been frequently encountered in hydro-climatic time-series (Liu et al. 2011a, 2011b, 2012b, 2013). Assuming a normal distribution at the significance level of P = 0.05, a positive Mann–Kendal statistic Z larger than 1.96 indicates a significant increasing trend, while a Z lower than −1.96 indicates a significant decreasing trend. Critical Z values of ±1.64 and ±2.58 are used for the significance of P = 0.1 and 0.01, respectively. The difference for factors between 1960–1990 and 1991–2013 was analyzed using one-way analysis of variance (ANOVA), the least significant difference test, and the F test (P < 0.05) in SPSS 19.0.

TSMC can be assumed to be negligible at multiyear or a longer time scale (Figure 3(a)). However, this assumption cannot be applied to evaluate ET_a at finer temporal scales (e.g., annual or monthly), as neglecting the interannual TWSC in the water balance method would significantly undermine the robustness for validation of interannual ET_a variation (Han et al. 2015). The variations of TWSC in our study region ranged from −31.1 mm yr⁻¹ to 32.9 mm yr⁻¹, accounting for 6% to 7% of the multiyear average precipitation. Thus, ΔS cannot be neglected when applying the water balance method to estimate the reference ET_a in the MYRB. After that, the reference ET_a was calculated based on the water balance method following Equation (1). The estimated annual reference ET_a showed a large variation for the last 50 years (Figure 3(b)). The multiyear average reference ET_a was approximately 466.9 mm yr⁻¹, accounting for more than 90% of the precipitation.

The MTE method, based on eddy covariance and highly accurate for measuring exchanges of carbon dioxide, water vapor, and energy between the biosphere and atmosphere as well as validating the other models (Fisher et al. 2008; Vinukollu et al. 2011), was employed in the AET_JUNG. Therefore, AET_JUNG was expected to have a higher capacity for capturing the reference ET_a variation than other models; however, the performance of AET_JUNG was moderate with underestimation of the reference ET_a by 11.8% and had a larger BIAS/RMSE (−54.7/76.9 mm yr⁻¹) and a relatively low CORR (0.35). The two DMs (AET_MODIS and AET_ZhangKe) both performed adequately, underestimating (19.5% and 18.9%, respectively) the ET_a variation with relatively large BIAS/RMSEs (−94.0/103.6 mm yr⁻¹ for AET_MODIS and −86.1/106.1 mm yr⁻¹ for AET_ZhangKe). Although AET_MODIS performed poorly with regard to BIAS/RMSE, it accurately reproduced the interannual variation of the reference ET_a, with the highest CORR (0.66). Therefore, AET_MODIS is recommended when higher spatial resolution is required. In contrast, interannual variation of AET_ZhangKe was poorly associated with the lowest CORR. However, AET_JRA55 from the reanalysis model had poor performance, strongly overestimating the ET_a by 25.7% with the largest BIAS/RMSE (120.6/135.0 mm yr⁻¹).

Therefore, compared to the other three categories, the LSMs performed the best in capturing the reference ET_a variation, especially for the Noah model. AET_JUNG had a good performance in BIAS and RMSE. The DMs (i.e., AET_MODIS and AET_ZhangKe) did not perform well compared to the LSM and MTE models. However, AET_JRA55 from the reanalysis model performed poorly in capturing the variation and amplitude of the reference ET_a.

There are several reasons for the annual differences in ET_a, including model types (e.g., LSM for GLDAS, machine-learning algorithm for the MTE method, the PM method for the DMs), forcing data (listed in Table 3) and parameter optimization. Compared to other ET_a categories, the LSMs, especially the Noah and Mosaic models, captured the variability of the reference ET_a fairly well. In the LSMs, the impact of soil moisture stress on bare soil evaporation and plant transpiration was considered. Soil evaporation happens only when the top soil can draw water from lower layers (Chen et al. 1996). Plant transpiration was also constrained by canopy interception and conductance (Koster & Suarez 1994). The Noah and Mosaic models could successfully simulate the energy and water flux, as reported by several studies (Chen et al. 1996; Koster & Suarez 1996). The VIC model was a macro-scale, semi-distributed, hydrological model (Liang et al. 1994), which calculated the infiltration, evaporation, and runoff variation within each grid cell (Franchini & Pacciani 1991). Therefore, it depicted the hydrological process with high accuracy and has been widely used in large river basins around the world (Sheffield et al. 2006; Livneh et al. 2013). This model also exhibited good performance in the MYRB.

AET_JRA55 was from the Japanese 55-year reanalysis (JRA-55) dataset of the Japan Meteorological Agency (Kobayashi et al. 2015). Although AET_JRA55 had adopted a series of strategies to improve the capacity for simulating ET_a, it still used the same forcing field as JRA-25 (i.e., atmospheric parameters such as pressure, temperature, humidity and wind, precipitation, downward solar and long wave radiation fluxes, and total cloud cover) (Nakaegwa 2008; Kobayashi et al. 2015), which can cause uncertainties in the reanalysis dataset. Further analysis of the precipitation (Pre_JRA55 in Figure 5(a)) and radiation (Ra_JRA55 in Figure 5(b)) showed that both were overestimated (19.8% and 19.7%, respectively) and had relatively low CORR (0.67 and 0.66, respectively). This can explain the large overestimation of AET_JRA55 and relatively low performance in capturing the variation of the reference ET_a. Therefore, great caution should be paid to hydrological variables before using them in reanalysis.

The AET_JUNG was trained with in situ flux tower measurement with an average duration of only 2 years (Jung et al. 2010). These measurements were strongly dependent on spatial gradients when reproducing the trends and variation in ET_a. In addition, the flux towers were sparse throughout our study region (Jung et al. 2010, Supplementary Figure 1), which may introduce some uncertainties to AET_JUNG in capturing the reference ET_a variation. Both of the ET_a from DMs were underestimated. The AET_ZhangKe did not consider the evaporation from canopy interception (Zhang et al. 2010), which may be a reason for the difference between the reference ET_a and AET_ZhangKe. The Chinese government adopted a series of strategies to restore the ecosystem since 1999 (Zhang et al. 2000; FAO 2010) which had significantly increased vegetation coverage (Chen et al. 2015; Wang et al. 2015) with a higher LAI and water storage capacity of the canopy (Zhang et al. 2007). Therefore, without considering the canopy interception, ET_a could be underestimated. In addition, insufficient accuracy of remotely sensed radiation employed in the AET_ZhangKe models may contribute to the uncertainty of ET_a (Zhang et al. 2012). Although AET_MODIS had a limited time span, it showed better performance in capturing the reference ET_a variability, which may result from the improvement of the ET_a algorithm (Mu et al. 2007, 2011). In addition, the forcing data in the AET_MODIS model maintained the spatial coherence. For example, the FPAR, LAI, albedo, enhanced vegetation index (EVI), and land cover were all from a single MODIS product, which could also explain the improvement in the estimation accuracy of AET_MODIS.

Evaporation from canopy interception and the uppermost soil layer and transpiration by the plant were all highly related to precipitation. Precipitation was essential for deriving ET_a and considered as the main input parameter in models. Especially in semi-arid and semi-humid regions where water supply limited ET_a, the degree of ET_a model dependence on precipitation could reflect its performance in estimating ET_a (Zhang et al. 2012). The CORR of precipitation and the nine ET_a products are listed in Table 4.

All of the ET_a from LSMs showed significant correlation (p < 0.01, Table 4) with precipitation and larger CORR (>0.52), which indicates that LSMs were more dependent on precipitation than the other methods. Precipitation, one of the key forcing fields in LSMs (Table 3), was improved significantly by using the near real-time data, which were derived from the geostationary satellite infrared cloud-top temperature measurements and various merged microwave observation datasets (Turk et al. 2000). In addition, the water balance method and improved techniques (i.e., adjusting the soil moisture and precipitation) were applied to constrain the LSMs (Koster & Suarez 1996; Rodell et al. 2004). Direct evaporation from the bare soil in the LSMs was closely related to the soil water content via the soil water diffusivity equation for moisture transport (Zeng & Decker 2009), which provided an effective and direct coupling between groundwater and surface water. Therefore, the variation of LSMs' ET_a followed the variation of precipitation and responded to changes in soil moisture in the MYRB.

The EC-based AET_JUNG and reanalysis AET_JRA55 exhibited low CORR with precipitation (0.30 and 0.28, respectively). The precipitation variable in the JRA55 model was derived from the reanalysis dataset, which was produced by combining the observation and forecast model results and had numerous uncertainties (Bosilovich et al. 2011). AET_JUNG was a data-driven estimation of ET_a and had low dependence on precipitation. There was a large difference in CORR with precipitation for the two DMs, high correlation for AET_MODIS (CORR = 0.62, P < 0.05) and low correlation for AET_ZhangKe (CORR = 0.12, P > 0.05). Both of the DMs used remotely sensed vegetation information (i.e., LAI or NDVI) and the PM method forced by meteorological data without precipitation (Mu et al. 2009; Zhang et al. 2010). Even without precipitation in AET_MODIS, there was a general agreement, which may be a result of a series of improvements in the terrestrial evapotranspiration algorithm (Mu et al. 2007, 2011), for instance, separating dry canopy surfaces from wet ones and dividing soil surface into saturated wet surface and moist surface, which was explicitly linked to precipitation. Therefore, we could conclude that the more dependent on precipitation of ET_a models, the better the performance in capturing the reference ET_a in semi-arid and semi-humid transition regions. The correlation between precipitation and ET_a from models can be considered as a standard to evaluate models' performance to some extent.

The reference ET_a exhibited large variations from 1960 to 2013 in the MYRB (Figure 3(b)). To obtain the temporal trend of the reference ET_a and detect which products can capture the tendencies, the study time span was defined from 1983 to 2013, considering the time consistency of each ET_a product (Table 1). In addition, there was a point (around 1999, Figure 3(b)) where the reference ET_a trend changed (i.e., decreasing from 1983 to 1998 and then increasing from 1999 to 2013). Some previous studies showed global land ET_a decreasing from the 1980s to the late 1990s (Jung et al. 2010; Mueller et al. 2013; Zeng et al. 2014) and then increasing significantly during the last two decades of the 20th century (Jung et al. 2010; Zhang et al. 2012). Therefore, 1999 was chosen as the change year to detect ET_a temporal variation across the MYRB. The annual mean ET_a and trend of the reference ET_a associated with ET_a products are listed in Table 5.

The annual average reference ET_a changed abruptly around 1999 with a mean value of 464.7 mm yr⁻¹ from 1983 to 2013, and the mean values for 1983–1998 (pre-1999) and 1999–2013 (post-1999) were 455.5 mm yr⁻¹ and 474.6 mm yr⁻¹, respectively (Table 5). Although the ET_a from products showed significant interannual variability (Figure 4), there was no significant difference between pre-1999 and post-1999 by ANOVA (p > 0.1) test. The annual trend in the reference ET_a slope was found to change from negative (−4.54 mm yr⁻²) for the pre-1999 period to positive (6.50 mm yr⁻²) for the post-1999 period (Table 5). Although there was a large increase in amplitude during the post-1999 period, it was offset by the negative trend of pre-1999, resulting in no significant trend for the MYRB during the entire period from 1983 to 2013 (1.04 mm yr⁻², p > 0.1).

Regarding the ET_a products for the pre-1999 period, except for AET_CLM, all the ET_a from LSMs, particularly AET_VIC, captured the reference ET_a negative trend, whereas the other products showed a positive trend. For the post-1999 period, all the ET_a products reproduced the reference ET_a's positive trend, especially for AET_JUNG and AET_MOS (5.14 and 4.17 mm yr⁻², respectively). For the entire period (1983 to 2013), the ET_a from Noah, JUNG, and MODIS showed positive trends, while the other methods showed negative trends. AET_JUNG performed best in capturing the trend amplitude. The Noah model, in general, exhibited good performance in capturing the ET_a trend for all the periods, while AET_VIC was best for pre-1999, and AET_JUNG was best for capturing the ET_a amplitude in both the post-1999 and 1983–2013 periods.

Further analysis revealed that the annual average ET_a made up 92.7% of the total precipitation for the pre-1999 period and 96.4% for the post-1999 period. The abrupt change in ET_a may have resulted from the series of strategies mandated by the Chinese Government since 1999 that aimed to restore the ecosystem, prevent severe soil erosion, and ensure sustainable development. These strategies led to increasing the land surface roughness (Liu & Zhang 2013) and then changed the streamflow and ET_a in the MYRB (Huang et al. 2003; Zhang et al. 2007).

Potential uncertainties from various sources may undermine the confidence of this study. First, the reference ET_a, based on the water balance method and with the assumption that the deep groundwater storage was not changed at annual scale, can cause uncertainties. Soil moisture was influenced by field capacity, soil bulk density, location (longitude and latitude), meteorological elements, natural factors, and human activities, all of which could cause various uncertainties (Li et al. 2005; Liu et al. 2008). The validation of soil moisture for all layers needs to be tested in this region; however, we only used six moisture stations at the 50 cm depth from 1992 to 2012. In addition, human activities, such as reservoir construction, irrigation and extraction of groundwater, may impact deep groundwater. Second, the time scale difference used for evaluating ET_a performance may impact conclusions. For instance, the AET_Noah2 started in 1948 was only evaluated from 1960. In contrast, the ET_a products from other LSMs were evaluated from 1979, and an even shorter time scale for AET_MODIS and AET_ZhangKe.

Despite the above uncertainties, the performance of the nine ET_a products was evaluated based on the water balance method in the MYRB. We also recommend that more attention be paid to the ET_a and runoff response to human activities and climate change (Zhang et al. 2013; Yang et al. 2015a), including quantitatively analyzing the relative contributions of different factors (e.g., irrigation, water diversion, global warming, land use change, and CO₂ enrichment) to ET_a. The mechanisms influencing the spatiotemporal variation in ET_a and its components (e.g., interception, soil evaporation, and transpiration) at different time scales (Koster et al. 2004; Good et al. 2015) should be further investigated; this will be beneficial for water resources management and irrigation scheduling in the MYRB. In addition, projections of the ET_a spatiotemporal trend in the future under different scenarios (Meehl et al. 2007; IPCC 2013) were useful for informing policy, such as whether to restart the GTGP in the Loess Plateau (Chen et al. 2015). More reliable knowledge about ET_a will be helpful to better understand hydrological processes and improve regional sustainable development.

The multiyear actual evapotranspiration (reference ET_a) was evaluated with the water balance method, which considered the soil moisture change derived from the VIC-IGSNRR model. The performance of the four categories of global ET_a products was investigated in the MYRB. The results indicated that the total water storage change (ΔS), estimated from the soil moisture of the VIC-IGSNRR model, varied significantly and the TWSC cannot be neglected when evaluating the ET_a by the water balance method at annual scale. The ET_a from LSMs performed excellently in capturing the variation of the reference ET_a among the four categories of products. The ET_a from Noah (especially for Noah1) and MOS demonstrated the best performance. AET_Noah2 is recommended for application in large regions and at a half-century scale, whereas AET_Noah1 and AET_MOS are recommended for the past three decades and AET_MODIS for the past 10-year scale. In addition, a good correlation was found between LSMs and precipitation, indicating that in arid and semi-arid regions, considering the precipitation in ET_a models can significantly improve the performance in capturing the variation of the reference ET_a. The AET_JUNG model, which was based on eddy covariance, performed relatively well, with a small negative BIAS and a little weak relationship to the reference ET_a, which demonstrated its dependence on observation data and independence from precipitation. The two DMs demonstrated fair performance with small underestimation. The AET_MODIS model, although a short time span, showed a stronger capacity for capturing the reference ET_a than the other models. However, the AET_JRA55 model performed poorly, with the largest positive BIAS and RMSE; this was caused by higher precipitation and radiation input.

The variation of the reference ET_a shifted around 1999. Further analysis found that the reference ET_a decreased for the pre-1999 period (−4.54 mm yr⁻²) and increased for the post-1999 period (6.50 mm yr⁻²). However, the increase during the post-1999 period was offset by the decrease during the pre-1999 period, resulting in no significant change during the entire period from 1983 to 2013. Generally, the ET_a from the Noah model, especially Noah1, better captured the ET_a trend and amplitude for the three periods. The AET_VIC model best captured the ET_a amplitude for the pre-1999 period, while the AET_JUNG model was best for the post-1999 and 1983–2013 periods. The underlying reasons for ET_a changes may be the implementation of several forest strategies. Comprehensive evaluation and comparison of various ET_a products was valuable for the application of ET_a products for regional water resource management and guiding future research on ET_a estimation. Limited by the coarse temporal and spatial resolution of ET_a products in this study, more accurate regional ET_a models with higher temporal and spatial resolution forcing data would be necessary to improve the estimation performance of ET_a in further work.

This research was financially supported by the Natural Science Foundation of China (No. 41330529, 41501032, 41571024, and 41301496). We thank two anonymous reviewers for their comments and suggestions to improve the quality of the manuscript and Xiaomang Liu, Guotao Dong, Kaiwen Wang, and Xiaocong Liu for their great help regarding the clarity of this paper.