Spatiotemporal change and attribution of potential evapotranspiration over China from 1901 to 2100

Global warming has accelerated surface water loss around the world. This study investigates in detail the change and attribution of potential evapotranspiration (PET) across China from 1901 to 2100 by the Hargreaves model, based on a 1-km temperature dataset downscaled from the low-spatial-resolution datasets using a Delta downscaling framework. Results showed that (1) relative to 1961–1990, PET increased by 0.62% in the historic period (1901–2017) and 6.43–12.89% for the future period (2018–2100), suggesting considerable future drying for China. Moreover, these increments had strong spatial variations and the largest increases were detected in high-elevation regions; (2) PET over entire China demonstrated a nonsignificant upward trend during the historic period and significant upward trends for the future period. For each period and GCM, significant upward PET trends occupied a much larger percent area than significant downward trends; (3) PET variations during the historic period were most sensitive to mean temperature (TMP), while in the future period it was more sensitive to maximum temperature (TMX), suggesting a change in the primary sensitivity factor due to global warming; (4) minimum temperature (TMN) made the largest contribution (45%) to PET variations during the historic period, while TMX had the largest contribution (36–40%) in the future period. Therefore, the primary contributing factor might transform from TMN to TMX under climate change; and (5) PET variations exhibited strong spatial heterogeneity, detected on fine geographic scales, due to the use of downscaled dataset. Overall, the results present a deep insight for planning coping strategies of global warming in China.


Introduction
Global warming has become an indisputable fact and has large long-term effects on water resources, plant growth, crop yields, and economic behavior (Naumann et al. 2018). Its impacts on water resources are the most obvious and extensive, because it can directly and quickly increase surface evaporation (Cook et al. 2014). Furthermore, global temperatures are projected to continuously increase (IPCC 2013), which will accelerate water loss from Earth's surface (Aouissi et al. 2016). Potential evapotranspiration (PET), one of crucial hydrological variables that can quantify the water loss of a region, could be used to calculate actual evapotranspiration (Aouissi et al. 2016), plan crop irrigation (Fan et al. 2018), project vegetation dynamics (Peng and Li 2018), and could act as input for ecosystem or hydrological models (Odusanya et al. 2019). Therefore, detailed spatiotemporal PET changes and attribution should be investigated to help develop sustainable strategies that cope with negative effect of warming in the globe.
Numerous approaches have been employed to estimate PET, including Budyko (Budyko 1974), Priestley-Taylor (Priestley and Taylor 1972), Hargreaves (Hargreaves and Samani 1985), and FAO Penman-Monteith (Allen et al. 1998) models. These models have been widely evaluated and employed to study PET variations around the world (Hargreaves and Allen 2003;Zhao et al. 2004;Yang et al. 2017;Odusanya et al. 2019). Furthermore, many studies have examined PET attribution using the FAO Penman-Monteith model and discussed the relative contributions of wind speed, solar radiation, mean temperature, and vapor pressure (Zheng et al. 2009;Han et al. 2018). However, few studies have investigated the specific influences of temperatures on PET changes under climate warming. The Hargreaves model adopts monthly minimum, maximum, and mean temperatures (TMN, TMX, and TMP, respectively) to calculate PET and thus can better attribute long-term PET variations to changes of these temperature variables, which in turn would enhance our understanding of how global warming induces surface water loss. Therefore, the Hargreaves model could be employed to investigate long-term PET changes and the attribution of such changes to a warmer climate.
Currently, climate data used to calculate PET typically comes from either station observations (Zheng et al. 2009; or proxy climate data (Cook et al. 2014;Yang et al. 2017). Although these data can produce reliable and long-term PET results, respectively, they cannot elucidate the detailed spatial patterns of PET because of poor weather station density and climatic proxy data spatial resolution (e.g., >50 km) (Harris et al. 2014;Gao et al. 2018). In contrast, high-spatial-resolution PET results would help develop suitable adaptation and mitigation strategies at subregional and even local scales (Giorgi et al. 2009;Peng et al. 2017;Peng et al. 2019). Therefore, to generate highspatial-resolution climate datasets containing a long PET time series is very urgent.
Several studies reported that the Delta downscaling method can be used to generate the high-spatial-resolution climate datasets (Brekke et al. 2013;Wang and Chen 2014;Peng et al. 2017). This method employs a low-spatial-resolution monthly time series and high-spatial-resolution reference climatology as inputs; the latter helped draw a fine-scale pattern of the climate that considered topographical effects (Mosier et al. 2014;Peng et al. 2017). The resulting downscaled dataset exhibited higher accuracy than raw datasets because of its ability to resolve topographical effects .
To date, many investigations regarding PET have been completed for China. However, these studies were conducted using datasets based on low-density weather stations with short time series (Chang et al. 2019;Jiang et al. 2019;Li et al. 2019;Liu et al. 2019;Xu et al. 2019) or climatic proxy data with low spatial resolution (Chai et al. 2018;Tian et al. 2018;Wang et al. 2019). Limited knowledge on the detailed spatial patterns of PET over China impedes the development of sustainable strategies to adapt and mitigate the climatic warming of China. Furthermore, little work concern to the detailed contributions of temperature changes to long-term PET variations, leaving a gap how climatic warming influences potential water loss across China.
Current work studied the changes and attribution of longterm PET in China. Specifically, the Delta downscaling framework and Hargreaves model were employed to generate high-spatial-resolution temperature time series and calculate PET for China, respectively. First, monthly TMN, TMX, and TMP with 0.5′ from January 1901 to December 2100 were generated based on climatic proxy data with 30′. Next, the downscaled temperatures and calculated PET were compared to observations. Spatiotemporal changes and trends in annual PET during 1901-2100 were then analyzed over China. Finally, the sensitivities and contributions of TMN, TMX, and TMP changes on annual PET variations over China were assessed.

Data acquisition and process
This study used the Climatic Research Unit (CRU) TS v. 4.02 datasets (Harris et al. 2014) and CMIP5 datasets simulated by 27 general circulation models (GCMs) (Reclamation 2013) with 30′ (approximately 55 km) as the climatic proxy data, including monthly TMN, TMX, and TMP. These two datasets are global and have been widely employed to study climate change. Currently, the CRU covers a period of 1901-2017, while CMIP5 has data for 1950-2100 under several representative concentration pathway (RCP) scenarios. For generating high-spatial-resolution climate dataset based on these two datasets, the dataset with 0.5′ (approximately 1 km) obtained from the WorldClim v. 2.0 was seen as the reference data. The reference dataset has mean monthly data for 1970-2000 and includes orographic effects and observed monthly climate information (Fick and Hijmans 2017).
The Delta method was employed to downscale the CRU and CMIP5 temperature data to 1-km temperature dataset over China. The produces of Delta spatial downscaling have been described in Peng et al. (2019). The new temperature dataset covers two periods, i.e., the historical period  and future period (2018-2100) under the RCP2.6, RCP4.5, and RCP8.5 scenarios. The temperature from 27 GCMs during 1950-2005 periods was also downscaled, which was used to select the best GCMs for temperature projections through comparing with the observations. This study employed the three best GCMs to downscale future temperatures and calculated PET, which may reduce the uncertainty of future PET changes in GCMs.
Monthly observed TMN, TMX, and TMP data during 1951-2016 and from 745 national weather stations ( Fig. 1) were downloaded from the National Meteorological Information Center of China (http://data.cma.cn/en), for evaluating the accuracy of downscaled data. In addition, the monthly mean observed pan evaporation from 1981 to 2010 was obtained from 1795 national weather stations across China (Fig. 1); this dataset was used to validate the calculated PET. Specifically, the evaporation pan in this study has a caliber of 20 cm and precision of 0.1 mm.

Calculating PET
This study employed the Hargreaves model (Hargreaves and Samani 1985;Hargreaves and Allen 2003) to calculate monthly PET based on downscaled 0.5′ monthly temperature data. Monthly PET (mm) is calculated as: where R a is monthly theoretical insolation (mm) (Allen et al. 1998), which is calculated using the solar constant, latitude, solar decimation, and day number. The unit of temperature is degrees Celsius (°C).

Evaluating downscaled temperatures and PET
To evaluate the 0.5′ downscaled temperatures, this study employed mean absolute error (MAE), the Nash-Sutcliffe efficiency coefficient (NSE), and root mean square error (RMSE) as followings: where P i is the raw/downscaled value, O i is the observed value, and n is the number of observations. After evaluating the downscaled values for the 27 GCMs during 1950-2005, the three best GCMs were chosen to project downscaled temperatures over China. In addition, this study used mean monthly observations of pan evaporation during the 1981-2010 period to compare to the calculated PET values. Linear statistic relationships were established for each season to reduce the influence of seasonality on PET.

Analyzing PET trends
To detect PET trends, the linear regression relationship between annual PET and time was established at each grid point over China, using: where b represents the annual trend of PET during a period. We used 10b for long-term PET trends (mm/decade). Trends were deemed significant at the 95% confidence level.

Attributing PET changes
Following the Hargreaves model (Eq. 1), PET is a function of R a , TMP, TMX, and TMN, and can be expressed as: The sensitivity of PET to an independent variable x was calculated from the derivative of PET with respect to x (∂PET/ ∂x). To compare the sensitivity of factors in Eq. 6, the partial derivative was transformed into a non-dimensional form (McCuen 1974): where S(x i ) is the sensitivity coefficient of climatic factor x i . According to Eq. (1), the sensitivity coefficient for R a is equal to 1. A positive sensitivity coefficient indicated that PET increases with the magnitude of the associated climatic factor, and vice versa. PET variations stimulated by a climatic factor were calculated through multiplying the long-term trend of the factor by its partial derivative . Accordingly, the contribution of climatic factors on PET changes was calculated using the following differential equation (Zheng et al. 2009): More simply: where TR PET is the long-term PET trend; C(R a ), C(TMP), C(TMX), and C(TMN) are the contributions of R a , TMP, TMX, and TMN changes on the long-term PET trend, respectively. In addition, the long-term trend of each climatic factor was calculated similarly to Eq. (5). Because R a in Eq.
(1) is a theoretical solar radiation, it is a constant on annual timescales for each grid point over China. Thus, the contribution of solar radiation on annual PET variations = 0, according to the Hargreaves model. To evaluate the robustness of the calculation, we compared the total contributions in Eq. (9) with the PET trend calculated in Eq. (5). Furthermore, the proportional absolute contribution of each climate factor was estimated as: where p(x) is the relative percent contribution of climatic factor x (i.e., TMP, TMX, and TMN), and ranges from 0 to 1.  Figure 3 shows annual PET change from 1901 to 2100 overall of China. Compared to the base period , PET exhibited a slight (0.62%) increase for the historic period . For the future period (2018-2100), PET increased by 6.43, 8.47, and 12.89% for RCP2.6, 4.5, and 8.5, respectively. The increase in the historic period was non-significant, while significant upward trends were detected for the future period with 2.45, 8.64, and 22.0 mm/decade under RCP2.6, 4.5, and 8.5, respectively. Figures 4 and 5 show spatial patterns of PET changes during the historic and future periods, and standard errors introduced by the three best GCMs are presented in Figs. S1 and S2. Compared to the historic period, PET exhibited large projected changes in the future period and its sub-  Table 2 presents the statistical information for each spatial pattern in Figs. 4 and 5. The results indicated that (1) PET changes in the historic and future periods increased by 0.66% and 7.1-14.24%, respectively, relative to the base period. This is similar to the PET changes over all of China (Fig. 3); (2) under each RCP, PET increased with time; and (3) PET increased more with each higher emission scenario in the future periods except 2018-2040. The coefficient of variation indicated that although the historic period had lower standard deviation (0.70%) than the future periods (2.81-8.78%), the historic period had larger spatial variation (106.06%) than future periods (37.58-50.87%). However, the difference between maximum and minimum showed that the future period had more extreme spatial variation (380.55-803.66%) than the historic period (38.38%), especially in 2071-2100 under RCP8.5 (1216.47%).   Table 3 lists the statistical information of the zones with significant trends across China. Because future PET trends in GCMs had large spatial differences, results are independently presented for BNU-ESM, CESM1-CAM5, and GFDL-ESM2M. Overall, the annual PET trends varied with the future RCP scenarios, GCMs, and study periods. Moreover, for each period and GCM, the percent area of significant increase trend was greater than that of decrease trend. Specifically, during historical period, significantly upward and downward trends were detected, with mean values of 3.48 and 2.33 mm/decade and percent areas of 15.14 and 40.67%, respectively. Based on average of three best GCMs, future period had both significantly upward and downward trends, with mean values of 3.89 ± 0.79 mm/decade and 4.57 ± 0.23 mm/decade under RCP2.6, 6.04 mm/decade and 8.91 ± 1.38 mm/decade under RCP4.5, and 7.56 ± 2.13 mm/decade and 21.99 ± 1.95 mm/decade under RCP8.5. Besides, they had percent areas of 1.55 ± 0.25% and 45.62 ± 13.90% under RCP2.6, 0.94% and 94.83 ± 2.32% under RCP4.5, and 0.20% and 99.85 ± 0.09% under RCP8.5. Moreover, the significant positive trend in the future period was larger for each higher emission scenario in each GCM.

Spatiotemporal changes of annual PET
The coefficient of variation showed that under RCP2.6, the BNU-ESM trend had the largest spatial variation (58.72%), while the smallest was detected from GFDL-ESM2M under RCP8.5 (12.52%). The BNU-ESM and GFDL-ESM2M trends under RCP8.5 showed the most extreme (41.66 mm/decade) and moderate (1.98 mm/decade) spatial variations, respectively. Moreover, spatial distributions of annual PET trends and their significance for future sub-periods were showed in Figs. S3-S8 and Tables S1-S3.    Fig. S9. Results indicated that (1) for each temperature, sensitivity coefficients in historic and future periods were quite spatially similar; and (2) positive and negative TMN sensitivity coefficients were mainly distributed over high-elevation regions (e.g., Tibetan Plateau, Tianshan Mountains, Aertai Mountains, Qilian Mountains, and northeast China) and low-elevation regions (e.g., Tarim Basin, center and south of China), respectively. The opposite distribution was observed for TMX and TMP. Table 4 lists corresponding geographic statistical information for the sensitivity coefficients. Results indicated that (1) the Mean of sensitivity coefficients for TMN and TMP were negative, ranging from −0.11 to −0.24 and − 0.26 to −0.82, respectively. This suggested decreasing PET with increasing TMN and TMP. However, TMX exhibited the opposite results, because its mean ranged from 0.61 to 0.74; (2) the absolute Mean for TMN and TMX increased from historic to future periods, while it decreased for TMP. Meanwhile, the absolute Means for TMP (− 0.82) and TMX (0.67 under RCP2.6,0.69 under RCP4.5,and 0.74 under RCP8.5) were the largest in the historic and future periods, respectively, This suggests a transformation in the primary sensitivity factor on PET variations under future warming; and (3) in the future period, the absolute Mean for TMN and TMX increased with higher emission scenarios, but decreased for TMP. According to Eqs. (1) and (7), the sensitivity coefficient of R a = 1 at each grid point over China, suggesting that a 10% increase in R a would induce a 10% PET increase. Compared to temperature, PET was overall more sensitive to solar radiation. However, the Min, Max, Std, and CV results in Table 4 suggested that PET was more sensitive to temperatures than solar radiation in some regions, such as southern China, for TMN and TMX, and the Tibetan Plateau, Tianshan Mountains, Aertai Mountains, Qilian Mountains, and northeast of China for TMP (Fig. 8).

Attribution of annual PET change
To further investigate the influences of temperatures on annual PET variations, the contribution of each temperature was analyzed. To validate the rationality of this method, we compared the total contributions of temperatures and the actual PET trend (Fig. 9); results showed good statistical consistency between them during the historic and future periods. Figure 10 maps spatial patterns of the contributions of TMN, TMX, and TMP on annual PET variations and their standard errors are showed in Fig. S10, while Table 5 lists the corresponding statistical information. Overall, the contribution of each temperature on annual PET variations was different from its sensitivity coefficient due to the incorporation of its trend. Statistics indicated that (1) a positive contribution was detected in the Mean for TMX and TMP for historic and future periods, ranging from 2.57 to 27.42 mm/decade and from 1.60 to 20.00 mm/decade, respectively. In contrast, a negative contribution was observed for TMN, ranging from −1.91 to −26.00 mm/decade; (2) TMN had the largest contribution (− 4.83 mm/ decade, 45%) in the historic period, followed by TMP (3.33 mm/decade, 32%) and TMX (2.57 mm/decade, 23%). However, TMX had the largest contribution (ranging from 2.70 to 27.42 mm/decade, from 40 to 36%) in the  Note: Min, Max, Mean, Std, and CV represent the minimum, maximum, mean, standard deviation, and coefficient of variation, respectively future under each RCP, followed by TMN (ranging from −1.91 to −26.00 mm/decade, from 31 to 35%) and TMP (ranging from 1.60 to 20.00 mm/decade, from 29 to 30%); and (3) with higher emission scenarios, the negative contribution for TMN and positive contribution for TMX and TMP both increased. The Min, Max, Std, and CV in Table 5 suggested that there were large spatial variations for the contributions of these temperatures. However, the CV in Tables 4 and 5 revealed that the sensitivity of each temperature was greater than its contribution. In addition, because R a in Eq. 1 is constant on annual timescales. Therefore, the actual annual PET variation was completely controlled by the aforementioned temperature variations. Figure 11 maps spatial patterns of the dominant temperature contribution on annual PET variations during historic and future periods, based on Eq. (10). Results indicated that (1) for both historic and future periods, high-elevation PET variations were dominated by TMP; and (2) low-elevation PET changes were primarily associated with TMX and TMN for the historic period, but TMX in the future period.

Conclusion and discussion
Although numerous works have studied PET changes, trends, and attribution in China, they typically used lowdensity weather station or low-spatial-resolution climatic proxy data (Chai et al. 2018;Tian et al. 2018;Chang et al. 2019;Jiang et al. 2019;Li et al. 2019;Liu et al. 2019;Wang et al. 2019;Xu et al. 2019). Moreover, they often employed the FAO Penman-Monteith model to attribute PET variations to wind speed, solar radiation, TMP, and vapor pressure (Zheng et al. 2009;Han et al. 2018). These approaches may limit understanding of PET variations at fine geographic scales and the detailed attribution of temperatures on PET changes as a result of climate change. This study investigated the detailed changes, trends, and attribution (including TMN, TMX, and TMP) of PET across China from 1901 to 2100, based on a 1-km temperature dataset downscaled from CRU and CMIP5 GCMs using the Delta downscaling framework and Hargreaves model. The results can help plan coping strategies of the global warming in China. The PET change over China in each period relative to the base period exhibits strong spatiotemporal variability (Figs. 3,4,and 5 and Table 2). Although PET increased only slightly (0.62%) for the 1901-2017 historical period, there were much larger increases (6.43-12.89%) for the 2018-2100 future period. These PET increases over all of China were induced by ongoing global warming (Cook et al. 2014). In general, PET increased more in high-elevation regions than low-elevation regions. This is because temperature increases in highelevation regions were greater than those in low-elevation regions (Peng et al. 2017).
PET trends and associated significance exhibited strong spatiotemporal variability over China during each time period (Figs. 6 and 7 and S3-S8, Tables 3 and S1-S3). Overall, significant increasing and decreasing trends were observed across different regions of China. However, the PA of the former was much greater than the latter, for each period and GCM. Chen and Frauenfeld (2014) reported that precipitation by the end of the century (2071-2100) would increase by 6% (5 mm/decade) under RCP2.6, 12% (11 mm/decade) under RCP4.5, and 16% (15 mm/decade) under RCP8.5. For the same period, we projected PET increases of 7.69% (12.83 mm/decade) under RCP2.6, 11.69% (17.13 mm/decade) under RCP4.5, and 21.51% (26.65 mm/decade) under RCP8.5 (Tables 2 and S3). This suggested that water loss over China may be accelerated in the future, which would threaten food and ecological security of this region. Moreover, the PET changes and trends showed deep spatial heterogeneity, suggesting that coping strategies should be planned at small geographic scale.
Current work attributed PET variations to temperatures (Figs. 8,9,10,and 11 and Tables 4 and 5), based on the Hargreaves model. The sensitivity analysis indicated that PET variations were more sensitive to TMP (− 0.82), followed by TMX (0.61) and TMN (− 0.11) in the historic period. In the future period, PET was more sensitive to TMX (from 0.67 to 0.74), followed by TMP (from −0.52 to −0.26) and TMN (from −0.24 to −0.17). Therefore, the main sensitivity factor on PET variations over China would transform from TMP to TMX and from negative to positive under global warming. However, the spatial pattern of each temperature's sensitivity would remain stable between the historic and future periods. The contribution analysis indicated that TMN had the largest contribution (45%) on PET variations during the historic period, followed by TMP (32%) and TMX (23%). In the future period, TMX had the largest contribution (36-40%), followed by TMN (31-35%) and TMP (29-30%). Therefore, the main contributing factor on PET variations over China will transform from TMN to TMX under climate change.
The above results allow us to infer that (1) in the historic period, although the PET was most sensitive to TMP over China, TMN had the largest contribution to PET variations; Note: Min, Max, Mean, Std, and CV represent the minimum, maximum, mean, standard deviation, and coefficient of variation, respectively and (2) future PET was most sensitive to TMX, which also had the largest contribution to future PET variations. However, the primary temperature contributor varied widely spatially, as TMP and TMX were dominant in high-and lowelevation regions, respectively (Fig. 11). In addition, although the total contributions of temperatures were close to the actual PET trend for each period, they were not an exact match (Fig.  9). Such differences were also detected by other studies (Zheng et al. 2009;, and may be due to the fact that only a first-order differential equation was used to calculate contributions. In other words, non-linear interactions among climatic factors may also impact PET variations ). Our future work will involve reducing uncertainties and improving the accuracy of resolved temperature contributions.
To cut down uncertainty of future temperature data (Hawkins and Sutton 2011;Ning and Bradley 2015;Ning and Bradley 2016), this study evaluated 27 GCMs by comparing historic downscaled temperatures from GCM and observations. The three best GCMs (i.e., GFDL-ESM2M, CESM1-CAM5, and BNU-ESM) were selected to downscale future temperatures. However, Wang and Chen (2014) downscaled raw GCM dataset to a higher spatial resolution over China and reported that EC-EARTH, with the smallest MAE (> 2.0°C), performed best to generate future monthly TMP across China. Moreover, the aforementioned four GCMs were evaluated in this and their studies; our study represented lower MAE (1.66°C) than Wang and Chen (2014). These differences are because our study adopted high-spatial-resolution reference climatology to generate the temperature and surface Fig. 10 Contributions of minimum, maximum, and mean temperatures on annual PET variations observations to evaluate the downscaled results. Overall, spatial uncertainty for climate change investigations can reduce result accuracy. Accordingly, this study highlighted the spatial heterogeneity of climate change and the associated development of sustainable strategies on fine spatial scales.

Supplementary Information
The online version contains supplementary material available at https://doi.org/10.1007/s00704-021-03625-w. Fig. 11 Dominant contribution factors on annual PET variations during historic and future periods. The red, green, and blue synthesis was performed on the stretched relative percent contributions of maximum, mean, and minimum temperatures (TMX, TMP, and TMN)