Abstract

Estimation of ET0 in high-elevation (above sea level, ASL) areas of Tibet presents unique challenges: scarcity of monitoring stations, short-time coverage of meteorological data, low-oxygen and low-pressure environment, and strong solar radiation. In this study, altitude factors and modified temperature constants are utilized to improve the Hargreaves (HS) model based on 30-year daily meteorological data from 20 typical sites in Tibet. The improved model, Hargreaves-Elevation (HS-E) improved model, could provide better results at different time scales. Negative ET0 values were unavoidable in the HS model when applied to high-elevation areas. However, the HS-E model solved this problem and improved the accuracy of estimated ET0. In particular, the HS-E model performs better when the time scale becomes larger. Therefore, the HS-E model is highly recommended to estimate ET0 in high-elevation areas where the meteorological data is scarce, for example, in Tibet above 2,000 m.

INTRODUCTION

Crop evapotranspiration (ET), the only term that appears both in water balance equations and surface energy balance equations (Xu & Singh 2005), is an essential component of the hydrologic cycle. Therefore, crop ET plays a key role in applying rainfall-runoff and ecological models, calculating regional water balance, and determining irrigation water demands. The first step of crop evapotranspiration estimation is estimating the reference crop evapotranspiration (ET0). Inaccurate ET0 estimates could lead to error propagation in the model and finally unreliable groundwater recharge estimates (Allen et al. 1998). Although there are many methods to estimate ET0, each of them is generally limited by geographical and meteorological conditions (Liu & Pereira 2001; Irmak et al. 2003). The applicability of the Penman-Monteith model (PM) has been shown in many regions of different climate types (Jabloum & Sahli 2008; Tabari et al. 2013). Thus, the Food and Agriculture Organization (FAO) has adopted the PM equation as the standard method to estimate ET0. The main shortcoming of the PM model is its requirement of some hard-to-obtain meteorological factors (Valiantzas 2015). It is very difficult to obtain continuous and comprehensive data in many cases. Furthermore, installation and maintenance of meteorological stations in high-elevation terrain is very expensive and complicated. Therefore, it is crucial to find a simple and accurate method to estimate ET0 for high-elevation areas.

The existing empirical equations to estimate ET0 with limited data were categorized as temperature-based, radiation-based, evaporation-based, material conversion-based and mixed-type (Tabari et al. 2013). The simplified temperature-based models, which only require easy-to-access data, have been widely used (Almorox et al. 2015). The temperature-based Hargreaves-Samani model (HS) (Hargreaves & Samani 1985) calculated solar radiation using the highest temperature, lowest temperature and external radiation. Almorox et al. (2015) found that the HS model was the simplest and most accurate temperature-based model to use in arid and semi-arid areas after studying data from 4,362 meteorological stations to estimate ET0 using 11 temperature-based formulas. Er-Raki et al. (2010) also showed that the HS model was the most accurate temperature-based model in semi-arid regions of the Americas. Some scholars analyzed the applicability of commonly-used ET0 empirical models in various climate zones throughout China. Hu et al. (2011) evaluated the applicability of the HS model at monthly scale in different climate zones of China. In addition, Wang et al. (2008) validated the applicability of the HS model in some arid areas and semi-arid areas of China.

Although the HS model has been popular for estimating ET0 due to its simplicity, there is evidence that ET0 is underestimated in arid regions and overestimated in humid areas by the HS model. Martinez & Thepadia (2010) found that ET0 was overestimated by the HS model after an analysis of data from 72 stations in Florida. Yoder et al. (2005) estimated the daily and monthly ET0 in the Cumberland Plateau in the southeastern humid region of the United States and found that the HS model was more suitable for estimating the monthly ET0 compared with the PM equation.

Many ET0 empirical models were only applied in certain areas (Jabloum & Sahli 2008), which means that the parameters in the empirical models had to be recalibrated for new areas. For example, Allen et al. (1998) pointed out that the value of the temperature coefficient in the HS model was not reasonable at altitudes over 1,500 m. Annandale et al. (2002) stated that the temperature coefficient and temperature index in the HS model should be modified to account for atmospheric pressures in different areas, and the atmospheric pressure decreased with altitude (Yang & Zhang 2009).

However, elevation factors, which were overlooked before, can be easily obtained and are closely related to the hygrometer constant (γ), the zenith radiation (Ra, solar radiation absorbed by the top of the earth's atmosphere) and the net solar radiation (Rn, the energy absorbed by the Earth's surface) (Allen et al. 1998). Existing studies on simplified models mostly focus on temperature, relative humidity, wind speed, sunshine hours, rainfall and other meteorological factors. The main objective of this study is to provide a simple and accurate model to estimate ET0 based on the HS model, using elevation information.

STUDY AREA AND DATA

Study area

Tibet, located in the southwest of China and known as the ‘roof of the world’, has an average altitude of 4,000 m covering an area of 1.23 million km2. The high altitude leads to unique climate characteristics: low oxygen, low pressure (less than 2/3 atmospheric pressure than at sea level), long sunshine duration (more than 3,000 hours per year) and strong radiation (annual solar radiation ranging from 6,000 to 8,000 MJ m−2). The humidity and temperature vary widely, and there are distinctive dry and wet seasons in the high-altitude areas of Tibet (above 2,000 m) due to the frequent cold-heat exchanges in the near surface layer. According to the characteristics of Tibetan geomorphology and climate, this study fully considered the altitude variations of Tibet (the overall trend is higher in the west and lower in the east). Twenty representative Tibetan sites with altitude ranging from 2,000 m to 5,000 m were selected (Figure 1) for the study. Their geographical data are shown in Table 1.

Figure 1

Location of the 20 representative sites in Tibet.

Figure 1

Location of the 20 representative sites in Tibet.

Table 1

Geographic information and annual average temperature of the 20 representative sites in Tibet, China

Station Elevation (ASL, m) Northern latitude (°) Eastern longitude (°) Tmean (°C) 
Chayu 2,327.6 28.39 97.28 11.91 
Bomi 2,736.0 29.52 95.45 9.02 
Milin 2,950 29.13 94.13 8.84 
Linzhi 2,991.8 29.41 94.20 8.76 
Changdu 3,315.0 31.09 97.10 7.38 
Gongga 3,555.3 29.18 90.59 8.78 
Zedang 3,560 29.16 91.46 8.70 
Luolong 3,640 30.45 95.50 5.72 
Lasa 3,648.9 29.38 91.08 7.38 
Shigatse 3,836 29.17 88.52 6.21 
Pulan 3,900 30.17 81.15 2.66 
Lazi 4,000 29.05 87.35 6.54 
Jiangzi 4,040 28.55 89.36 5.51 
Damxung 4,200 30.29 91.06 1.81 
Shiquanhe 4,278.6 32.30 80.05 1.21 
Tingri 4,300 28.38 87.05 3.12 
Gerze 4,414.9 32.09 84.25 1.23 
Nagqu 4,507 31.29 92.04 −1.07 
Shenzha 4,672 30.93 88.70 −0.2 
Ando 4,800 32.04 91.87 −2.70 
Station Elevation (ASL, m) Northern latitude (°) Eastern longitude (°) Tmean (°C) 
Chayu 2,327.6 28.39 97.28 11.91 
Bomi 2,736.0 29.52 95.45 9.02 
Milin 2,950 29.13 94.13 8.84 
Linzhi 2,991.8 29.41 94.20 8.76 
Changdu 3,315.0 31.09 97.10 7.38 
Gongga 3,555.3 29.18 90.59 8.78 
Zedang 3,560 29.16 91.46 8.70 
Luolong 3,640 30.45 95.50 5.72 
Lasa 3,648.9 29.38 91.08 7.38 
Shigatse 3,836 29.17 88.52 6.21 
Pulan 3,900 30.17 81.15 2.66 
Lazi 4,000 29.05 87.35 6.54 
Jiangzi 4,040 28.55 89.36 5.51 
Damxung 4,200 30.29 91.06 1.81 
Shiquanhe 4,278.6 32.30 80.05 1.21 
Tingri 4,300 28.38 87.05 3.12 
Gerze 4,414.9 32.09 84.25 1.23 
Nagqu 4,507 31.29 92.04 −1.07 
Shenzha 4,672 30.93 88.70 −0.2 
Ando 4,800 32.04 91.87 −2.70 

Data collection

The daily meteorological data including the minimum and maximum temperatures, humidity, sunshine duration, wind speed, rainfall and elevation from the 20 sites over 30 years were obtained from the National Meteorological Information Center (NMIC, http://data.cma.cn). The daily and monthly ET0 were calculated using the Penman-Monteith equation. Then, the 30-year meteorological data were divided into two sub sets: the training dataset – the daily weather data of the 20 sites from 1981 to 2000 (n = 135144, no interpolation); and the testing dataset – the daily weather data of the 20 sites from 2001 to 2010 (n = 73,040, no interpolation). An elevation-adjusted HS model was developed using the training dataset, and validated using the testing dataset.

METHODS

FAO56 Penman-Monteith equation

The FAO stated that all calculation procedures for the PM equation can be standardized using available meteorological data and time scales, and the results can be used by project managers, project consultants, irrigation engineers, hydrologists, agronomists and meteorologists as reference standards. Taking into account the aerodynamic term (ET0aero) and radiation term (ET0grad), the PM equation is often used as a standard equation suitable for most climatic areas (Allen et al. 1998). The Penman-Monteith equation is 
formula
(1)
where ET0 PM is the ET0 (mm day−1) calculated by Penman-Monteith equation and will be used as the true value of ET0 in this paper; Rn is the net radiation at the crop surface (MJ m−2day−1); G is the soil heat flux density (MJ m−2 day−1, ); T is the air temperature at a 2 m height (°C); u2 is the wind speed at a 2 m height (m s−1); es and ea are the saturated and actual vapour pressure (kPa), respectively; Δ is the slope vapour pressure curve (kPa, °C) and γ is a psychrometric constant (kPa, °C).

Hargreaves-Samani equation

The HS equation estimates reference crop evapotranspiration (ET0) only using temperature data (Hargreaves & Samani 1985). 
formula
(2)
where ET0 HS is the ET0 calculated by HS equation (mm day−1); Tmean, Tmax, Tmin are the mean, maximum and minimum temperature (°C), respectively; Ra is the extraterrestrial radiation, which is calculated as a function of the day of the year and latitude and longitude (Allen et al. 1998).

The HS equation was tested based on lysimeter measurements over an eight-year period (Hargreaves & Samani 1985). The HS equation can utilize atmospheric temperature differences (Tmax–Tmin, indicating how much of the zenith radiation can reach the earth's surface) to estimate solar radiation. Many studies showed that the HS equation gave a reasonable ET0 throughout the world in the condition of some missing meteorological data. However, the empirical coefficients, such as the temperature coefficient (0.0023), temperature constant (−17.8) and temperature index (0.5), were obtained through regression analysis of many observations collected from around the world over many years. Thus, the parameters in the HS equation should be modified before it is used in a new area.

Methods to determine the influence factors on ET0 in high-elevation areas

The influence of elevation and temperature on the ET0

The analysis of the variation of ET0 with elevation and temperature showed that the annual ET0 increased gradually with increasing elevation (decreasing temperature) in areas below 4,000 m, but decreased with increasing elevation (decreasing temperature) in areas above 4,000 m (Figure 2 and Table 1). The annual highest cumulative value of ET0 (about 1,280 mm/year on average) was found in the Lazi area (about 4,000 m ASL), and the lowest cumulative value of ET0 (about 950 mm/year on average) was found in Chayu (the lowest altitude site selected in this study, about 2,000 m ASL) and the Ando area (the highest altitude site selected in this study, about 5,000 m ASL). In summary, elevation and temperature were very important in ET0 estimation, and this was why we modified the temperature constant and the elevation function to improve the HS model.

Figure 2

The influence of elevation on annual ET0 in the high elevation areas. The horizontal axis indicates the elevation in 20 representative sites (range: 2,000–5,000 m); the vertical axis represents the annual cumulative ET0.

Figure 2

The influence of elevation on annual ET0 in the high elevation areas. The horizontal axis indicates the elevation in 20 representative sites (range: 2,000–5,000 m); the vertical axis represents the annual cumulative ET0.

The ET0 and the annual average temperature were abnormal at the Damxung area (4,200 m above sea level), which has a lower annual ET0 compared to those of Tingri (4,300 m above the sea level) and Lazi (4,000 m above sea level) (Figure 2). This was due to the unique terrain and landscape of Damxung, which led to lower annual average temperature (1.81 °C) in Damxung than those of Tingri (3.12 °C) and Lazi (6.54 °C), and further caused the lower annual cumulative evapotranspiration in the Damxung area.

Principal component analysis

Principal component analysis (PCA) is a statistical method that is used to extract the most important information (i.e. variance) in a dataset (Wang 2015; Wang et al. 2018). The new variables obtained from the PCA are called principal components (PCs), which are linear combinations of the original variables, and each PC is orthogonal to other PCs.

When calculating ET0 using the PM equation, seven meteorological parameters were used, including the daily maximum temperature, daily minimum temperature, daily average temperature, daily average wind speed, daily average relative humidity and daily sunshine duration, all of which were obtained from the NMIC. The non-continuous data series contained geospatial data in three dimensions: longitude, latitude and altitude. Due to the large amount of continuous meteorological parameters needed for this complex calculation equation, a simpler alternative is often necessary when data is sparse. Statistical Product and Service Solutions (SPSS) software (IBM, USA) was used to analyze the PCs based on the geospatial data and 30-year meteorological data from the 20 study sites in Tibet. The original data was automatically standardized to eliminate the influence of the different variance scales of the variables.

The first three PCs (Table 2) had eigenvalues larger than 1, and the total variance explained by them was nearly 80%. Therefore, the first three PCs could be used to illustrate the main information contained in the original nine variables. The correlations between the first three PCs and the original variables were used to find the relationship between the PCs and the original meteorological indexes or geographic information (Table 3). There were significant correlations between the first principal component and four of the initial variables: the daily minimum temperature, daily mean temperature, daily maximum temperature and elevation. The second and third PCs only had significant correlation with sunshine duration and latitude, respectively.

Table 2

Eigenvalues and explained variance

Principal component Eigenvalue Explained variance (%) Cumulative explained variance (%) 
4.174 46.375 46.375 
1.612 17.906 64.281 
1.171 13.009 77.290 
0.957 10.633 87.923 
0.576 6.398 94.321 
0.318 3.529 97.850 
0.142 1.582 99.432 
0.046 0.515 99.947 
0.005 0.053 100 
Principal component Eigenvalue Explained variance (%) Cumulative explained variance (%) 
4.174 46.375 46.375 
1.612 17.906 64.281 
1.171 13.009 77.290 
0.957 10.633 87.923 
0.576 6.398 94.321 
0.318 3.529 97.850 
0.142 1.582 99.432 
0.046 0.515 99.947 
0.005 0.053 100 
Table 3

PC coefficients

Original index First principal component Second principal component Third principal component 
Minimum daily temperature 0.913 0.268 0.240 
Mean temperature 0.885 0.412 0.162 
Maximum daily temperature 0.846 0.488 0.070 
Height above sea level (elevation) −0.783 0.249 0.416 
Longitude 0.622 −0.547 −0.091 
Latitude −0.589 −0.085 0.677 
Mean of air relative humidity 0.587 −0.417 0.533 
Sunshine duration −0.378 0.703 −0.204 
Average wind speed −0.088 0.309 0.340 
Original index First principal component Second principal component Third principal component 
Minimum daily temperature 0.913 0.268 0.240 
Mean temperature 0.885 0.412 0.162 
Maximum daily temperature 0.846 0.488 0.070 
Height above sea level (elevation) −0.783 0.249 0.416 
Longitude 0.622 −0.547 −0.091 
Latitude −0.589 −0.085 0.677 
Mean of air relative humidity 0.587 −0.417 0.533 
Sunshine duration −0.378 0.703 −0.204 
Average wind speed −0.088 0.309 0.340 

The PCA analysis further confirmed that temperature and elevation were important indicators of ET0.

Hargreaves-Elevation model

Including an altitude function in the HS model

Through analysis of the PCs, elevation data was shown to be important due to its easy access and lack of need for continuous observations. Therefore, the elevation factor can be included in the HS equation. Elevations of the study sites were easily accessible from digital elevation model (DEM) or handheld GPS and remained constant in a very long period; therefore, it was included to improve the HS equation: 
formula
(3)
where f(H) is the altitude function to be solved; a is the temperature constant, which is often the minimum value of Tmean in the last 50 years in areas higher than 2,000 m to prevent the calculation of ET0 HS-E from generating negative values. This modified HS-E model represents an improvement on the HS model. Other variables are the same as in Equation (2).

An improved method to determine the temperature constant

FAO recommended that the temperature constant should be −17.8, when the unit of Ra was measured in mm d−1. However, Tmean is often lower than −17.8 °C in high-altitude areas, especially in areas above 3,500 m in January and December. This led to negative ET0 values calculated by the HS model. This result is contrary to common sense (actual evapotranspiration cannot be negative). Therefore, the method to determine the temperature constant needs to be improved.

The lowest daily Tmean of the 20 representative sites (Table 1, elevation ranging from 2,000 to 5,000 m) from 1960 to 2015 was −36.6 °C, which was observed at the Gerze region (4,414.9 m altitude) on December 26th, 1987. Therefore, the temperature constant was set to −36.6 to avoid negative values of ET0 calculated using the new model.

Derivation of the elevation function

Comparing Equations (1) and (3), it was obtained 
formula
(4)
Plugging in the parameters obtained from monthly summary of the training dataset (the meteorological data from the 20 typical sites from 1981 to 2000), the value of f(H) could be determined (Figure 3). Assuming f(H) is a linear function of H, linear regression analysis suggested that the elevation function is as follows: 
formula
(5)
Figure 3

The linear regression analysis of elevation (ASL) function (n = 12 × 20 × 20) (*** represents significance at the 0.1% level; the standard error = 1.7 × 10−4).

Figure 3

The linear regression analysis of elevation (ASL) function (n = 12 × 20 × 20) (*** represents significance at the 0.1% level; the standard error = 1.7 × 10−4).

Then, the HS-E model is: 
formula
(6)

Model calibration and validation

The new HS-E model was validated using the testing dataset (the meteorological data from the 20 typical sites from 2001 to 2010) at different time scales.

The mean relative error (MRE), the root mean square error (RMSE) and Nash–Sutcliffe efficiency (NSE) were used as statistical indices to show the applicability of the HS-E model (Cai et al. 2007). The calculation process of each indicator is as follows: 
formula
(7)
 
formula
(8)
 
formula
(9)
where N is the total number of error comparison points; Pi and Oi represent the predicted and measured values, respectively (i=1, 2, … , N).

The MRE is a percentage and is used to characterize the accuracy and error scope of the prediction; the RMSE is used to measure the deviation between the observed values and true values. The closer the MRE and RMSE are to zero, the higher reliability of the model. The closer the NES value is to 1, the higher the reliability of the model.

RESULTS

Performance on estimating daily ET0

Daily ET0 (n = 3,653 × 20) was estimated by the HS model and the HS-E model using the testing dataset. Then RMSE, MRE and NSE were calculated, using the results from PM equation as the true values, to evaluate the performance of the HS model and the HS-E model (Table 4).

Table 4

Performance of the HS model and the HS-E model on daily ET0 estimation in the 20 sites of high-elevation region

Station District Elevation (m) Method NSE RMSE (mm day−1MRE (%) ETHS/ET0 PM; ETHS-E/ET0 PM 
Chayu Linzhi 2,327.6 HS 0.71 0.60 18.36 1.16 
HS-E 0.85 0.44 −1.51 0.98 
Bomi Linzhi 2,736 HS 0.35 0.79 27.43 1.27 
HS-E 0.75 0.49 18.72 1.19 
Milin Linzhi 2,950 HS 0.69 0.66 18.41 1.18 
HS-E 0.74 0.6 1.13 1.24 
Linzhi Linzhi 2,991.8 HS 0.66 0.55 11.22 1.11 
HS-E 0.83 0.40 7.04 1.07 
Changdu Changdu 3,315 HS 0.69 0.66 18.41 1.18 
HS-E 0.74 0.60 1.13 1.24 
Gongga Shannan 3,555.3 HS 0.84 0.51 1.63 1.02 
HS-E 0.8 0.57 9.08 1.09 
Zedang Shannan 3,560 HS 0.81 0.54 −2.09 0.98 
HS-E 0.82 0.54 3.07 1.03 
Luolong Changdu 3,640 HS 0.8 0.55 −3.52 0.96 
HS-E 0.78 0.58 9.33 1.09 
Lasa Lasa 3,648.9 HS 0.82 0.53 −6.85 0.93 
HS-E 0.60 0.80 −13.49 0.87 
shigatse Shigatse 3,836 HS 0.84 0.51 −1.05 0.99 
HS-E 0.80 0.59 14.58 1.15 
Pulan Ali 3,900 HS 0.73 0.70 −22.17 0.78 
HS-E 0.88 0.46 −2.78 0.97 
Lazi Shigatse 4,000 HS 0.65 0.76 −9.70 0.90 
HS-E 0.73 0.66 4.56 1.05 
Jiangzi Shigatse 4,040 HS 0.8 0.52 −0.61 0.94 
HS-E 0.79 0.94 13.76 1.14 
Damxung Lasa 4,200 HS 0.79 0.53 −13.88 0.86 
HS-E 0.79 0.53 17.07 1.17 
Shiquanhe Ali 4,278.6 HS 0.74 0.85 −28.43 0.72 
HS-E 0.89 0.54 5.87 1.06 
Tingri Shigatse 4,300 HS 0.75 0.61 −13.87 0.86 
HS-E 0.83 0.50 13.98 1.14 
Gerze Ali 4,414.9 HS 0.64 0.87 −29.80 0.70 
HS-E 0.86 0.54 6.22 1.06 
Nagqu Nagqu 4,507 HS 0.74 0.56 −17.69 0.82 
HS-E 0.73 0.58 21.42 1.21 
Shenzha Nagqu 4,672 HS 0.41 0.89 −29.25 0.71 
HS-E 0.76 0.56 6.43 1.06 
Ando Nagqu 4,800 HS 0.51 0.75 −31.64 0.68 
HS-E 0.77 0.51 15.21 1.15 
Average HS 0.69 0.65 15.13(ABS) 0.94 
   HS-E 0.79 0.56 9.80 1.09 
Station District Elevation (m) Method NSE RMSE (mm day−1MRE (%) ETHS/ET0 PM; ETHS-E/ET0 PM 
Chayu Linzhi 2,327.6 HS 0.71 0.60 18.36 1.16 
HS-E 0.85 0.44 −1.51 0.98 
Bomi Linzhi 2,736 HS 0.35 0.79 27.43 1.27 
HS-E 0.75 0.49 18.72 1.19 
Milin Linzhi 2,950 HS 0.69 0.66 18.41 1.18 
HS-E 0.74 0.6 1.13 1.24 
Linzhi Linzhi 2,991.8 HS 0.66 0.55 11.22 1.11 
HS-E 0.83 0.40 7.04 1.07 
Changdu Changdu 3,315 HS 0.69 0.66 18.41 1.18 
HS-E 0.74 0.60 1.13 1.24 
Gongga Shannan 3,555.3 HS 0.84 0.51 1.63 1.02 
HS-E 0.8 0.57 9.08 1.09 
Zedang Shannan 3,560 HS 0.81 0.54 −2.09 0.98 
HS-E 0.82 0.54 3.07 1.03 
Luolong Changdu 3,640 HS 0.8 0.55 −3.52 0.96 
HS-E 0.78 0.58 9.33 1.09 
Lasa Lasa 3,648.9 HS 0.82 0.53 −6.85 0.93 
HS-E 0.60 0.80 −13.49 0.87 
shigatse Shigatse 3,836 HS 0.84 0.51 −1.05 0.99 
HS-E 0.80 0.59 14.58 1.15 
Pulan Ali 3,900 HS 0.73 0.70 −22.17 0.78 
HS-E 0.88 0.46 −2.78 0.97 
Lazi Shigatse 4,000 HS 0.65 0.76 −9.70 0.90 
HS-E 0.73 0.66 4.56 1.05 
Jiangzi Shigatse 4,040 HS 0.8 0.52 −0.61 0.94 
HS-E 0.79 0.94 13.76 1.14 
Damxung Lasa 4,200 HS 0.79 0.53 −13.88 0.86 
HS-E 0.79 0.53 17.07 1.17 
Shiquanhe Ali 4,278.6 HS 0.74 0.85 −28.43 0.72 
HS-E 0.89 0.54 5.87 1.06 
Tingri Shigatse 4,300 HS 0.75 0.61 −13.87 0.86 
HS-E 0.83 0.50 13.98 1.14 
Gerze Ali 4,414.9 HS 0.64 0.87 −29.80 0.70 
HS-E 0.86 0.54 6.22 1.06 
Nagqu Nagqu 4,507 HS 0.74 0.56 −17.69 0.82 
HS-E 0.73 0.58 21.42 1.21 
Shenzha Nagqu 4,672 HS 0.41 0.89 −29.25 0.71 
HS-E 0.76 0.56 6.43 1.06 
Ando Nagqu 4,800 HS 0.51 0.75 −31.64 0.68 
HS-E 0.77 0.51 15.21 1.15 
Average HS 0.69 0.65 15.13(ABS) 0.94 
   HS-E 0.79 0.56 9.80 1.09 

Note: ABS = absolute value.

First, the NSE value of the HS-E model was between 0.60 and 0.88 with an average of 0.79; while the NSE value of the HS model was between 0.35 and 0.84 with an average of 0.69 (Table 4). Thus, with NSE closer to 1, the simulation results of the HS-E model were more accurate than those of the HS model. Second, the RMSE of the HS-E model was between 0.40 and 0.80 mm day−1, with an average of 0.56 mm day−1; while the RMSE of the HS model was between 0.51 and 0.89 mm day−1, with an average of 0.65 mm day−1 (Table 4). Third, the MRE of the HS-E model was between 1.13% and 18.72% with an average (the average MRE absolute value) of 9.80%, while the MRE of the HS-E model was between −31.64% and 18.41% (the average absolute value was 15.13%) (Table 4). Therefore, with the RMSE and MRE closer to 0, the values of ET0 calculated by the HS-E model had a smaller deviation compared with the true values. Fourth, when plotting against the ET0 PM, the slope of ET0 HS-E ranged from 1.03 to 1.24 with an average of 1.09; while the ET0 HS ranged from 0.63 to 1.21. Therefore, the HS-E model was more stable and its estimation was closer to the true value. In summary, when estimating daily ET0, the HS-E model avoided the occurrence of negative values and provided estimations closer to true values calculated by the PM equation, compared with the HS model.

The performance on estimating monthly ET0

Monthly ET0 (n = 120 × 20) was also estimated by the HS model and the HS-E model using the testing dataset (Figure 4), and then the RMSE, MRE and NSE were calculated, using the results from PM equation as the true values, to evaluate the performance of the HS model and the HS-E model (Table 5). First, the NSE value of HS-E model for monthly ET0 was between 0.77 and 0.98 with an average of 0.89; while the NSE of HS model for monthly ET0 ranged from 0.06 to 0.93, with an average of 0.69 (Table 5). Therefore, the estimated results of the HS-E model were of higher quality than those of the HS model. Second, the average RMSE of 9.52 mm month−1 for ET0 obtained by the HS-E model was less than the average RMSE of 15.51 mm month−1 for ET0 obtained by the HS model. Third, the MRE of ET0 according to HS-E model varied from −4.26% (Pulan) to 22.10% (Changdu), with an average (the average of MRE absolute values) of 7.21%; while the average MRE of the HS model was 15.51%. Therefore, with the RMSE and MRE closer to 0, the values of ET0 calculated by the HS-E model had a smaller deviation compared with the true values.

Figure 4

Performance on estimating monthly ET0 by HS and HS-E models using data from 2000 to 2010.

Figure 4

Performance on estimating monthly ET0 by HS and HS-E models using data from 2000 to 2010.

Table 5

Performance of Monthly ET0 calculated by HS model and HS-E model in the 20 sites of high-elevation region

Station District Elevation (m) Method NSE RMSE (mm month−1MRE (%) ETHS/ET0 PM; ETHS-E/ET0 PM 
Chayu Linzhi 2,327.6 HS 0.65 15.66 16.85 1.17 
HS-E 0.90 8.21 −3.21 0.97 
Bomi Linzhi 2,736 HS 0.06 22.38 25.01 1.25 
HS-E 0.77 11.19 1.14 1.16 
Milin Linzhi 2,950 HS 0.33 16.45 13.50 1.14 
HS-E 0.83 8.72 9.06 1.09 
Linzhi Linzhi 2,991.8 HS 0.64 13.54 9.13 1.09 
HS-E 0.94 5.32 4.74 1.05 
Changdu Changdu 3,315 HS 0.71 17.11 16.86 1.17 
HS-E 0.78 14.84 22.10 1.22 
Gongga Shannan 3,555.3 HS 0.92 8.76 −0.23 1.00 
HS-E 0.94 8.66 5.98 1.09 
Zedang Shannan 3,560 HS 0.91 10.02 −4.41 0.96 
HS-E 0.94 8.03 0.02 1.00 
Luolong Changdu 3,640 HS 0.88 9.98 −5.92 0.94 
HS-E 0.95 7.24 6.18 1.06 
Lasa Lasa 3,648.9 HS 0.92 9.94 −8.29 0.92 
HS-E 0.94 8.23 −1.86 0.98 
shigatse Shigatse 3,836 HS 0.93 9.06 −2.95 0.97 
HS-E 0.91 10.43 11.25 1.11 
Pulan Ali 3,900 HS 0.76 18.51 −22.69 0.77 
HS-E 0.95 8.84 −4.26 0.96 
Lazi Shigatse 4,000 HS 0.72 17.47 −12.28 0.88 
HS-E 0.88 11.54 0.83 1.01 
Jiangzi Shigatse 4,040 HS 0.82 11.88 −7.45 0.93 
HS-E 0.86 11.63 11.05 1.11 
Damxung Lasa 4,200 HS 0.85 11.71 −15.06 0.85 
HS-E 0.89 9.94 13.88 1.14 
Shiquanhe Ali 4,278.6 HS  0.70 22.67 −28.95 0.77 
HS-E 0.96 9.42 2.78 1.03 
Tingri Shigatse 4,300 HS 0.80 14.86 −15.18 0.85 
HS-E 0.91 9.78 11.31 1.11 
Gerze Ali 4,414.9 HS 0.67 22.60 −30.09 0.70 
HS-E 0.98 5.60 2.00 1.02 
Nagqu Nagqu 4,507 HS 0.77 13.89 −18.58 0.81 
HS-E 0.77 14.00 18.27 1.18 
Shenzha Nagqu 4,672 HS 0.36 23.48 −30.55 0.69 
HS-E 0.93 7.86 2.22 1.02 
Ando Nagqu 4,800 HS 0.47 20.14 −32.33 0.68 
HS-E 0.84 10.90 12.10 1.12 
Average HS 0.69 15.51 15.82(ABS) 0.92 
   HS-E 0.89 9.52 7.21 1.07 
Station District Elevation (m) Method NSE RMSE (mm month−1MRE (%) ETHS/ET0 PM; ETHS-E/ET0 PM 
Chayu Linzhi 2,327.6 HS 0.65 15.66 16.85 1.17 
HS-E 0.90 8.21 −3.21 0.97 
Bomi Linzhi 2,736 HS 0.06 22.38 25.01 1.25 
HS-E 0.77 11.19 1.14 1.16 
Milin Linzhi 2,950 HS 0.33 16.45 13.50 1.14 
HS-E 0.83 8.72 9.06 1.09 
Linzhi Linzhi 2,991.8 HS 0.64 13.54 9.13 1.09 
HS-E 0.94 5.32 4.74 1.05 
Changdu Changdu 3,315 HS 0.71 17.11 16.86 1.17 
HS-E 0.78 14.84 22.10 1.22 
Gongga Shannan 3,555.3 HS 0.92 8.76 −0.23 1.00 
HS-E 0.94 8.66 5.98 1.09 
Zedang Shannan 3,560 HS 0.91 10.02 −4.41 0.96 
HS-E 0.94 8.03 0.02 1.00 
Luolong Changdu 3,640 HS 0.88 9.98 −5.92 0.94 
HS-E 0.95 7.24 6.18 1.06 
Lasa Lasa 3,648.9 HS 0.92 9.94 −8.29 0.92 
HS-E 0.94 8.23 −1.86 0.98 
shigatse Shigatse 3,836 HS 0.93 9.06 −2.95 0.97 
HS-E 0.91 10.43 11.25 1.11 
Pulan Ali 3,900 HS 0.76 18.51 −22.69 0.77 
HS-E 0.95 8.84 −4.26 0.96 
Lazi Shigatse 4,000 HS 0.72 17.47 −12.28 0.88 
HS-E 0.88 11.54 0.83 1.01 
Jiangzi Shigatse 4,040 HS 0.82 11.88 −7.45 0.93 
HS-E 0.86 11.63 11.05 1.11 
Damxung Lasa 4,200 HS 0.85 11.71 −15.06 0.85 
HS-E 0.89 9.94 13.88 1.14 
Shiquanhe Ali 4,278.6 HS  0.70 22.67 −28.95 0.77 
HS-E 0.96 9.42 2.78 1.03 
Tingri Shigatse 4,300 HS 0.80 14.86 −15.18 0.85 
HS-E 0.91 9.78 11.31 1.11 
Gerze Ali 4,414.9 HS 0.67 22.60 −30.09 0.70 
HS-E 0.98 5.60 2.00 1.02 
Nagqu Nagqu 4,507 HS 0.77 13.89 −18.58 0.81 
HS-E 0.77 14.00 18.27 1.18 
Shenzha Nagqu 4,672 HS 0.36 23.48 −30.55 0.69 
HS-E 0.93 7.86 2.22 1.02 
Ando Nagqu 4,800 HS 0.47 20.14 −32.33 0.68 
HS-E 0.84 10.90 12.10 1.12 
Average HS 0.69 15.51 15.82(ABS) 0.92 
   HS-E 0.89 9.52 7.21 1.07 

On the other hand, first, the NSE of 0.89 for the monthly ET0 HSE was greater than 0.79 at the daily scale. Second, the RMSE of 9.52 mm month−1 for the monthly ET0 HSE was lower than the 0.53 mm day−1 at the daily scale. Third, the MRE of 7.21% for the monthly ET0 HSE was less than the 7.44% at the daily scale in 20 different stations, meaning the monthly ET0 calculated by the HS-E model was more accurate than that at daily scale. Therefore, the HS-E model has smaller errors with the increasing calculation period and more obvious advantages over the HS model.

So, the simulation quality of the HS-E model was much better than that of the HS model at monthly time scale.

The performance on estimating seasonal (three-month) ET

Seasonal ET0 (n = 40 × 20) was estimated by the HS model and the HS-E model using the testing dataset. Then RMSE, MRE and NSE were calculated, using the results from PM equation as the true values, to evaluate the performance of the HS model and the HS-E model (Figure 5). First, The NSE of HS-E model for the seasonal ET0 was between 0.74 and 0.99 with an average of 0.90; while the NSE of HS model for seasonal ET0 ranged from −0.09 to 0.94, with an average of 0.69 (Figure 5). Second, the average RMSE of HS-E model with 25.47 mm season−1 was less than the average RMSE of HS model with 42.39 mm season−1. Third, the MRE of HS-E model varied from −4.74% (Lasa) to 20.44% (Changdu), with an average (the average of MRE absolute values) of 7.0%; while the HS model had an average MRE of 15.43%. Therefore, the HS-E model was superior to the HS model.

Figure 5

Performance on estimating seasonal ET0 by HS and HS-E models using data from 2000 to 2010.

Figure 5

Performance on estimating seasonal ET0 by HS and HS-E models using data from 2000 to 2010.

The ET0 values calculated by the HS-E model for seasonal data produced an average NSE of 0.90 at the seasonal scale, slightly greater than the 0.89 at the monthly scale and obviously greater than the 0.79 at the daily scale. An average RMSE of 25.47 mm season−1 was found for the seasonal scale, slightly lower than the 9.52 mm season−1 for the monthly scale and significantly lower than the 0.53 mm day−1 for the daily scale. The seasonal ET0 calculated by the HS-E model was slightly more reliable than calculated at the monthly scale, and both were much better than that at the daily scale. Therefore, the HS-E model has smaller errors with the increasing calculation period and with more obvious advantages over the HS model.

So, the simulation quality of the HS-E model was much better than that of the HS model at larger time scales.

DISCUSSION

Due to the difficulty in obtaining continuous climate data for ET0 calculation by PM method, many studies were conducted to simplify the calculation method of ET0 (Hargreaves & Allen 2003; Trajkovic 2005; Traore et al. 2010; Marti et al. 2011). It should be noted that the HS model used in this study is simpler compared with other existing studies. Many previous studies added continuous meteorological parameters, such as the relative humidity and sun hours, into the HS model to enhance its computational accuracy (Valiantzas 2015). However, increasing the complexity of the HS model also leads to higher demands for resources, including obtaining basic meteorological data and running a more complicated simulation (Ishak et al. 2010), which undermined the essence of HS model: using a simple temperature algorithm that only requires a few parameters.

The atmospheric temperature difference (Tmax–Tmin) was utilized as an indicator to express the amount of zenith radiation reaching the surface of the earth in the HS model. Ren et al. (2016a, 2016b) estimated the daily ET0 for arid to sub-humid areas in Inner Mongolia, China base on HS model. However, due to the unique nature of Tibet, the temperature at high elevation was not as high as expected under strong radiation (China Tibet Information Center 2005). The lower than expected atmospheric temperature resulted from the thinner atmosphere and lower absorption of long-wave radiation from the ground at higher elevation. Meanwhile, there is less cloud cover at higher elevation areas. And the less cloud cover means less long-wave radiation absorption from the ground in the day time and weaker counter radiation of the cloud against the ground at night, both leading to poor insulation (Allen et al. 1998). All of the above reasons weakened the sensitivity of temperature differences (Tmax–Tmin) to the radiation. Therefore, the elevation factor was included to enhance the ability of the model to estimate radiation and improve the accuracy of the HS model.

This paper mainly focused on the high-elevation areas of Tibet, so validation of the HS-E equation was also limited to that region. However, adjusting the temperature constant and using the elevation factor, the proposed HS-E model has the potential to easily estimate ET0 in a large variety of areas with easy to obtain parameters.

CONCLUSIONS

The standard PM method (adopted by FAO) requires a lot of continuous meteorological data that is difficult to obtain in areas of complex terrain and high elevation. And the popular HS model could not avoid negative ET0 values due to the inappropriate temperature constant. We adjusted the temperature constant in HS model by taking the minimum of the annual average temperature in 50 years; and we included an elevation factor in the HS model to better relate temperature difference to radiation. The elevation factor was derived by principle component analysis with geospatial data in Tibet. The improved HS model, HS-E model, was calibrated using meteorological data from 1981 to 2000 at 20 representative sites in Tibet, and validated using meteorological data from 2001 to 2010 at the same sites (using the results from PM equation as the true values).

The HS-E model not only avoided the occurrence of negative values in the ET0 estimation, but also enhanced the accuracy of the ET0 estimation. The accuracy of the HS-E model was shown to be superior to that of the HS model under different time scales (daily , monthly and seasonal) using our testing dataset. The NSE of the HS-E model reached 0.78 (daily), 0.89 (monthly) and 0.90 (seasonal), showing the high reliability of the model. Meanwhile, by analyzing the results of the NSE, RMSE and MRE at different time scales, the monthly ET0 was found to be closer to the PM standard ET0 compared to the daily ET0. The larger the time scale, the better the HS-E model performed. Overall, the HS-E model is feasible and convenient for use at high elevation areas. However, this study only used data from sites with elevations ranging from 2,000 to 5,000 m. In the future, it is recommended to consider the effect of elevation using more meteorological data at sites covering wider ranges of elevations.

ACKNOWLEDGEMENTS

This study was supported by the National Science Foundation of China (51579158, 51609154), State Key Laboratory of Simulation and Regulation of Water Cycle in River Basin (SKL2018ZY01) and scientific researching fund projects of China Institute of Water Resources and Hydropower Research (MK2017J02). The contributions of the editor and anonymous reviewers whose comments and suggestions significantly improved this article are also appreciated.

REFERENCES

REFERENCES
Allen
R. G.
,
Pereira
L. S.
,
Raes
D.
&
Smith
M.
1998
Crop Evapotranspiration Guidelines for Computing Crop Water Requirements
.
FAO Irrigation and Drainage. Paper 56
.
FAO
,
Rome, Italy
.
Annandale
J. G.
,
Jovanic
N. Z.
,
Benade
N.
&
Allen
R. G.
2002
Software for missing data error analysis of Penman-Monteith reference evaportransportation
.
Irrigation Science
21
,
57
67
.
Cai
J. B.
,
Liu
Y.
,
Lei
T. W.
&
Pereira
L. S.
2007
Estimating reference evapotranspiration with the FAO Penman-Monteith eqution using daily weather forecast messages
.
Agricultural Forest Meteorology
145
,
22
35
.
China Tibet Information Center
2005
Geography of Tibet, http://info.tibet.cn/newzt/rsxzzt/tyy/t20050310_15780.htm [accessed online 15 June 2010]. (In Chinese)
Er-Raki
S.
,
Chehbouni
A.
,
Khabba
S.
,
Simonneaux
V.
,
Jarlan
L.
,
Ouldbba
A.
,
Rodriguez
J. C.
&
Allen
R.
2010
Assessment of reference evapotranspiration methods in semi-arid regions: can weather forecast data be used as alternate of ground meteorological parameters
.
Journal of Arid Environment
74
,
1587
1596
.
Hargreaves
G. H.
&
Allen
R. G.
2003
History and evaluation of Hargreaves evapotranspiration equation
.
Journal of Irrigation and Drainage Engineering
129
,
53
63
.
Hargreaves
G. H.
&
Samani
Z. A.
1985
Reference Crop Evapotranspiration From Ambient air Temperature
.
Paper No. 85-2517
,
ASAE regional meeting
,
Grand Junction, Colorado
.
Hu
Q. F.
,
Yang
D. W.
,
Wang
Y. T.
&
Yang
H. B.
2011
Global calibration of Hargreaves equation and its applicability in China
.
Advance Water Science
22
,
160
167
.
Irmak
S.
,
Allen
R. G.
&
Whitty
E. B.
2003
Daily grass and alfalfa-reference evapotranspiration ratios in Florida
.
Journal of Irrigation and Drainage Engineering
129
,
360
370
.
Ishak
A. M.
,
Bray
M.
,
Remesan
R.
&
Han
D. W.
2010
Estimating reference evapotranspiration using numerical weather modelling
.
Hydrological Processes
24
,
3490
3509
.
Liu
Y.
&
Pereira
L. S.
2001
Calculation methods for reference evapotranspiration with limited weather data
.
Journal of Hydraulic Enginnering
3
,
11
17
.
Marti
P.
,
Gonzalez-Altozano
P.
&
Gasque
M.
2011
Reference evapotranspiration estimation without local data
.
Irrigation Science
29
,
479
495
.
Martinez
C. J.
&
Thepadia
M.
2010
Estimating reference evapotranspiration with minimum data in Florida
.
Journal of Irrigation and Drainage Engineering
136
(
7
)
484
501
.
Tabari
H.
,
Grismer
M. E.
&
Trajkovic
S.
2013
Comparative analysis of 31 reference evapotranspiration methods under humid conditions
.
Irrigation Science
31
,
107
117
.
Trajkovic
S.
2005
Temperature-based approaches for estimating reference evapotranspiration
.
Journal of Irrigation and Drainage Engineering
131
,
316
323
.
Wang
C.
2015
Escherichia Coli Transport Modeling at Soil Column Scale and Watershed Scale
.
Cornell University
,
Ithaca, NY
.
Wang
S. F.
,
Duan
A. W.
&
Zhang
Z. Y.
2008
Comparison and analysis of Hargreaves equation and Penman-Monteith equation during the different hydrological years in the semi-arid region
.
Transactions of the CSAE
24
,
29
33
.
Wang
C.
,
Schneider
R. L.
,
Parlange
J. Y.
,
Dahlke
H. E.
&
Walter
M. T.
2018
Explaining and modeling the concentration and loading of Escherichia coli in a stream – A case study
.
Science of the Total Environment
635
,
1426
1435
.
Yang
Y. H.
&
Zhang
Z. Y.
2009
Method for calculating Lhasa reference crop evapotranspiration by modifying Hargreaves
.
Advance Water Science
20
,
614
618
.