Evaluation of temporal variability and stationarity of potential evapotranspiration in Peninsular Malaysia Open Access

Potential Evapotranspiration (PET) can be expressed as the maximum amount of water that can evaporate and transpire from a vegetated landscape without restrictions other than atmospheric demand (Thornthwaite 1948). Evapotranspiration occurs concurrently with the process of evaporation and transpiration. The evapotranspiration process is influenced by several climatic factors such as air temperature, wind speed, relative humidity and solar radiation (Lang et al. 2017). Over recent years, the capacity of PET has varied significantly due to human activities that are disastrous for both the agricultural sector as well as water resources management (Um et al. 2020). PET is one of the key components in the hydrological cycle as it is essential for proper water management practices, agricultural irrigation and scheduling.

There are six widely used PET estimation methods, which can be further classified into radiation-based methods and temperature-based methods. The Priestley–Taylor (PT) (Priestley & Taylor 1972), Makkink (MA) (Makkink 1957), and Turc (TU) (Turc 1961) methods are among the radiation-based methods whereas the Hargreaves–Samani (HS) (Hargreaves & Samani 1985), Thornthwaite (TH) (Thornthwaite 1948) and Blaney–Criddle (BC) (Blaney & Criddle 1950) methods are the temperature-based methods. Lu et al. (2005) compared the six estimation methods for long-term annual PET in the south-eastern United States. The results showed that the estimated PET values were highly correlated. However, the temperature-based methods showed bigger differences in the multivariate statistical tests compared with the radiation-based methods. In addition, Rao et al. (2011) compared the FAO-56 (Allen et al. 1998), Hamon (Hamon 1963) and PT methods at monthly and annual temporal time-scales in the Southern Appalachians. Based on the results, the PT method provided the best estimates of forest PET at two forested watersheds. It was also found that large underestimations of the forest PET were provided by the uncorrected FAO-56 and Hamon methods. Additionally, Tabari et al. (2013) adopted eight pan-evaporation-based methods, seven temperature-based methods, four radiation-based methods, and ten mass-transfer-based methods and compared them against the Penman–Monteith FAO-56 (PMF-56) model for PET estimation in Iran. The best estimation methods obtained from the study were the Snyder, BC and Romanenko for the pan-evaporation-based method, temperature-based method and mass-transfer-based method, respectively. For the radiation-based method, Tabari et al. (2013) developed a new equation using air temperature and solar radiation showed the best performance as compared with the existing available method in the estimation of PET. In Malaysia, Tukimat et al. (2012) adopted six methods, namely PT, MA, TU, HS, TH and BC, to estimate the historical PET and predict the future PET. The results showed that the radiation-based methods were capable of predicting future PET values more accurately than the temperature-based methods. Additionally, the radiation-based methods predicted closer future PET values as compared with those computed using the Penman–Monteith (PM) method (Penman 1948) with a statistical downscaling technique.

The trend analysis of PET is vital as it is one of the main variables contributing to climate change and temporal–spatial vector trends. According to Mansour et al. (2017), both parametric and non-parametric tests can be conducted to detect the patterns or trends in hydro-climatological time series. Shadmani et al. (2012) investigated the temporal trend of reference evapotranspiration (ET_o) using the Mann–Kendall (MK) and Spearman's Rho (SR) tests for various time-scales in Iran. Based on the findings, both increasing and decreasing trends were identified in the study region and the performances of the nonparametric trend tests (MK and SR) applied were consistent at the 5% significance level. Meanwhile, the slope severity was classified using the Theil–Sen test (Sen 1968; Theil 1992). Although the MK and SR tests were unique in their application, the values of significant trend demonstrated a minimal difference. On the other hand, Hu et al. (2019) adopted the MK test, modified MK test and Sen's slope estimator to analyse the annual and seasonal trends for various hydro-climatic variables in Kamo River Basin, Japan. Significant trends were detected for all hydro-climatic variables except for annual river discharge at the seasonal scale. Among the three methods, the modified MK test was capable of improving the trend analysis using a structure with full autocorrelation. Similarly, Ng et al. (2020) investigated the trends for monthly, seasonal and yearly rainfall series in Kelantan River Basin, Malaysia, using the MK test and Sen's slope estimator. From the results, it was evident that the northern part of the Kelantan River Basin exhibited significant increasing trends while the southwest region demonstrated decreasing trends.

Several studies on the stationarity of hydro-climatic variables have been carried out by researchers across the world using different methods. There are five commonly used methods: Augmented Dickey–Fuller (ADF) (Dickey & Fuller 1979), Dickey–Fuller test (Fuller 1976) statistic using generalised least squares (DF-GLS), Elliot–Rothenberg–Stock (ERS) (Elliott et al. 1996), Kwiatkowski–Phillips–Schmidt–Shin (KPSS) (Kwiatkowski et al. 1992) and Phillips–Perron (PP) (Phillips & Perron 1988) tests. Khalili et al. (2013) adopted these five methods to evaluate the stationarity of streamflow for the daily, monthly and yearly time series in West Azarbaijan of Iran. The results revealed that the daily, monthly and yearly time series were stationary upon the removal of streamflow with trend components from the time series studied. Moreover, Adewole & Serifat (2015) examined the stationarity of the rainfall series from the year 1971 to 2010 in Nigeria using the ADF and PP tests. The findings indicated that all rainfall series analysed were stationary. Additionally, Rutkowska & Ptak (2017) applied the ADF and KPSS tests to analyse the stationarity for average monthly and annual flow of three rivers (two Polish and one American river). It was found that the ADF and KPSS tests provided similar results without the presence of seasonality. However, with the presence of seasonality, the results of the KPSS test became distorted. Besides that, Ukkola et al. (2019) investigated the statistical significance of trends for precipitation, temperature and pan evaporation using the covariance stationarity method at annual and seasonal time-scales in Australia. It was found that stationarity was not detected for temperature at both annual and seasonal time-scales. Ng et al. (2020) adopted the KPSS test to evaluate the stationarity of several rainfall series in Kelantan River Basin, Malaysia. The results indicated that the rainfall series were stationary except at two rainfall stations during Inter-monsoon 1, Inter-monsoon 2 and yearly series.

The implications due to the changes in the PET of a tropical region with high humidity are much higher for the regional and global climate. Besides that, the Southeast Asia region which constitutes the maritime continent plays a major role in convective activities around the globe. Thus, the analysis of the trend of PET in a tropical region like Malaysia is highly important for understanding both local and global climate systems (Pour et al. 2020). In addition to that, the most diverse ecosystem is commonly found in the tropical region and such an ecosystem is very sensitive to climate variability. Since PET is an important element that controls the climatic water balance, a huge impact may fall on the biodiversity of the ecosystem even if a small change occurs in the PET (Hui-Mean et al. 2018). On the other hand, Malaysia has encountered significant climate change in recent years with an increase in daily temperature and climatic models have projected that this situation will persist for the coming years. As such, the impacts of global warming in Malaysia have resulted mainly due to the rise in temperature and other hydrological variables such as the PET. Although various studies have been carried out to comprehend the driving forces behind the changes in rainfall and temperature, studies related to the trend analysis of the PET are still very limited in tropical regions like Malaysia. Therefore, this study aims to analyse the trend and conduct a stationarity test on the PET within Peninsular Malaysia. Historical data of seven stations located in Peninsular Malaysia from year 2009 to year 2018 were used for the estimation of PET. The results of the study can be beneficial in providing an appropriate scientific framework for water resource planning and management strategies.

Malaysia is a country located in Southeast Asia that comprises two distinct parts which are peninsular and East Malaysia. Peninsular Malaysia was formed by 11 states with a total area of 131,794 km². Its weather is hot and humid throughout the year. The main focus was on seven stations in Peninsular Malaysia, namely Bayan Lepas station, Ipoh station, Subang station, KLIA Sepang station, Kuantan station, Muadzam Shah station and Pulau Langkawi station.

The historical meteorological data used were obtained from the Malaysia Meteorological Department (MMD) for the period of year 2009 to year 2018. The daily climate data obtained were the 24 hours mean, maximum and minimum temperature, mean surface wind speed, relative humidity, solar radiation and rainfall. For the monthly climate data, the mean, maximum and minimum temperature were obtained. The following information for the studied stations such as station code, longitude, latitude and record period are shown in Table 1. The geographical locations of the seven stations are shown in Figure 1.

The temporal patterns of trends in the meteorological and time series data of the study area were analysed. The Mann–Kendall (MK) test was used to identify the significance of the trend for the time series while the Sen's slope estimator test was used to determine the trend magnitude.

The Z value is used to determine the significance of a trend at the chosen α significance level. In a two-sided test, the null hypothesis is accepted if |Z| is larger than and the rejection of the null hypothesis defines the trend as significant.

For the stationarity test of PET, the Augmented Dickey–Fuller (ADF) test and Kwiatkowski–Phillips–Schmidt–Shin (KPSS) test were used.

According to Dickey & Fuller (1979), if the ρ value is less than ⁠, there is a presence of unit root in that time series.

The PET estimated using the Turc method and the Makkink method were compared with the Penman method for the daily, monthly and annual PET time series. The Penman method has been selected as the reference model for comparison since it is considered as the most physical and reliable method in the verification of other PET methods (Lang et al. 2017; Lakatos et al. 2020; McColl 2020). As shown in Tables 2–4, the Turc method achieved better statistical performances as compared with the Makkink method for all the stations and the time series studied. This is obvious and clear where the Turc method obtained a lower percentage of bias (PBIAS), normalised root mean square (NRMSE) and mean bias error (MBE) values. Similarly, the Turc method also acquired higher coefficient of determination (R²) values as compared with the Makkink method, which shows that the Turc method is able to produce a closer estimation of PET than the Penman method. A huge difference in the PET values obtained by the Makkink method could be due to increasing evaporative conditions that lead to underestimation of PET values as shown in Figure 2. The Turc method was also recommended to be adequate as a radiation-based method compared with the Makkink method by Liu et al. (2017), who also obtained similar findings. Overall, the Ipoh station acquired the lowest PBIAS, NRMSE and MBE values among the other six stations studied for the daily, monthly and annual PET time series. For the daily PET time series, the Ipoh station acquired PBIAS, NRMSE and MBE values of −0.050, 0.056 and −0.002, respectively. As for the monthly PET time series, the Ipoh station acquired PBIAS, NRMSE and MBE values of −0.050, 0.035 and −0.062, respectively. For the annual PET time series, the Ipoh station acquired PBIAS, NRMSE and MBE values of −0.050, 0.020 and −0.743, respectively. As such, it can be concluded that the Turc method is able to estimate PET values with higher accuracy. On the other hand, it is to be noted that the Muadzam Shah station obtained the highest R² values using the Turc method in the estimation of PET for the daily, monthly and annual PET time series, which were 0.985, 0.976 and 0.971, respectively. The reason that the Turc method is able to achieve such performance is due to its simplicity. The Turc method only requires three parameters for PET estimation, thus making it easy to be applied. Hence, the Turc method was selected as the best PET estimation method since it acquired PET values that were closer to the PET values of the Penman method. Similar findings were also obtained by Tukimat et al. (2012) where the Turc method was found to be more sensitive towards changes in temperature due to climate change among the radiation-based methods applied.

The Mann–Kendall (MK) test was used to determine the existence of an increasing or decreasing monotonic trend whereas Sen's slope estimator was used to estimate the magnitude of the trend. From the MK test, the trend can be identified as increasing if the Z value obtained is positive and vice versa. Meanwhile, the trend is significant at the 5% significance level if the value obtained is less than 0.05. Based on Table 5, Pulau Langkawi, KLIA Sepang and Kuantan stations exhibited increasing trends since the Z values are positive while Bayan Lepas, Ipoh, Subang and Muadzam Shah stations exhibited decreasing trends since the Z values are negative for the daily PET time series. Among these stations, only Pulau Langkawi, Subang and Kuantan stations exhibited significant trends with values of 0.001, 0.036 and <0.0001, respectively. These stations also acquired greater Sen's slope values of 0.00005, −0.00001 and 0.00008, respectively, which indicates a higher rate of change in the daily PET trend.

For the monthly time series, Pulau Langkawi, Ipoh, KLIA Sepang and Kuantan stations exhibited increasing trends since the obtained Z values were positive. The Bayan Lepas, Subang and Muadzam Shah stations exhibited decreasing trends since the obtained Z values were negative, as shown in Table 6. None of the stations showed significant increasing trend at the significance level of 5% since their values were more than 0.05. However, Pulau Langkawi, Bayan Lepas and Kuantan stations acquired greater Sen's slope values of 0.03, −0.049 and 0.058, respectively, which indicates a higher rate of change in the monthly PET trend. A previous study by Hui-Mean et al. (2018) has also revealed similar results in which most of the stations in Peninsular Malaysia exhibited an increasing trend for the monthly PET time series especially at the Ipoh station. The impact of the increasing trend for the monthly PET time series also indicates that the increasing temperature may lead to more severe drought especially in Ipoh.

As for the annual time series, Pulau Langkawi, Ipoh, Subang, Muadzam Shah and Kuantan stations exhibited increasing trends since the obtained Z values were positive while Bayan Lepas and KLIA Sepang stations exhibited decreasing trends since the obtained Z values were negative, as shown in Table 7. Among these stations, only Kuantan station showed a significant increasing trend at the significance level of 5% with a value of 0.047. In contrast, Pulau Langkawi and Kuantan stations showed greater Sen's slope values of 11.254 and 10.978, respectively, which showed a higher rate of change in the annual PET trend. A possible reason that the PET values obtained were much lower in most of the central mountainous stations such as Bayan Lepas, Ipoh, Subang and Muadzam Shah was the higher relative humidity. In comparison, most of the coastal region stations such as Pulau Langkawi, KLIA Sepang and Kuantan achieved higher PET values due to the lower relative humidity. These results were consistent with the findings of Pour et al. (2020) where there is an increase of the annual evapotranspiration in Peninsular Malaysia. Additionally, the increase in the trend of evapotranspiration in Peninsular Malaysia after the 1990s was also contributed by the large increase in global temperature as revealed by Di Liberto et al. (2018). As such, this has also proved the effect of global warming in the rise in the PET trend in Peninsular Malaysia.

Overall, a combination of increasing and decreasing trends was detected across daily, monthly, and annual PET time series. It is worth mentioning that the Pulau Langkawi and Kuantan stations exhibited increasing trends throughout the PET time series analysed. According to Suhaila & Yusop (2018), the large increase in mean temperature may contribute to the rapid change in the PET trend. Apart from that, a larger loss of water through evapotranspiration may be recorded due to the gradual rise in the overall temperature and the effects of global warming (Roshan et al. 2013). Thus, the increasing trend of PET is likely to reduce water catchment yield, water supply and water availability assuming that there is no change in the rainfall patterns for the particular region.

Trend analysis may be able to detect the monotonic trend in a PET time series, however, it cannot determine the stationarity of the daily, monthly and annual PET time series. Therefore, two tests were conducted to check for the stationarity of the PET time series. Firstly, the Augmented Dickey–Fuller (ADF) test was applied to check the existence of a unit root for the time series at the significance level of 5%. The unit root is said to be present when the value is less than = 0.05. Secondly, the Kwiatkowski–Phillips–Schmidt–Shin (KPSS) test was used to check the stationarity of the time series at the significance level of 5%. The PET time series is said to be stationary when the null hypothesis is accepted with a value greater than 0.05. In contrast, the PET time series is said to be non-stationary when the null hypothesis is rejected with a value less than 0.05.

For the daily PET time series, there was the presence of a unit root for all stations. These stations acquired values that were less than in the ADF tests. Meanwhile, the values of the KPSS test at Pulau Langkawi, KLIA Sepang and Kuantan stations were less than indicating that the daily PET time series were non-stationary, as shown in Table 8. According to Tramblay et al. (2013), global warming may influence the change in the behaviours of a hydro-climatic time series. The changes in the behaviours of hydro-climatic parameters can induce non-stationarity in the time series.

For the monthly PET time series, the presence of a unit root was identified at all stations except for the KLIA Sepang station using the ADF test since the values obtained for KLIA Sepang station were greater than ⁠, as shown in Table 9. Besides that, the KPSS test applied on the monthly PET time series also found stationarity at all stations except for the KLIA Sepang station with a value less than 0.05. Therefore, it was observed that the means and variances of the PET time series for the KLIA Sepang station varied significantly over time. This can be explained by the occurrence of an uneven cycle or trend that has altered the means and variances gradually over time, as explained by Nashwan et al. (2019).

For the annual PET time series, there was no presence of unit root identified using the ADF test at all stations, as shown in Table 10. As for the KPSS test, non-stationarity was detected at the Kuantan station with a value of 0.019 for the annual PET time series. To the best knowledge of the authors, there is very limited literature on the stationarity test of PET especially in Peninsular Malaysia. As such, the results of this study were not able to be compared with the findings from other studies. However, several studies on the stationarity of reference evapotranspiration have been carried out in other countries. For example, Kazemi et al. (2020) applied both the ADF test and KPSS test in Iran. The results revealed that the reference evapotranspiration time series in Iran contained a unit root and were not stationary. Similarly, Mirdashtvan et al. (2020) applied the KPSS test on various hydro-climatic variables such as precipitation, mean air temperature and pan evaporation in Iran.

The focus in this study was to carry out the trend analysis and stationarity test for daily, monthly and annual PET times series at seven stations in Peninsular Malaysia. A comparison was made for the different estimation methods of the PET time series. From the results, the Turc method was found to show better accuracy in the estimation of PET values in a warm and humid country like Malaysia. After the PET estimation, the Mann–Kendall (MK) test and Sen's slope estimator were applied for the trend analysis of the PET time series obtained from the Turc method. As for the stationarity test of the PET time series, the Augmented Dickey–Fuller (ADF) test and Kwiatkowski–Phillips–Schmidt–Shin (KPSS) test were applied.

The trend analysis using the MK test and Sen's slope estimator were applied on the PET trends for the daily, monthly and annual time series. Both increasing (positive) and decreasing (negative) trends were identified across all the time series analysed. In short, Pulau Langkawi and Kuantan stations were found to exhibit increasing trends while Bayan Lepas station exhibited decreasing trends for the all the three daily, monthly and annual PET time series. Additionally, significant trends were identified at Pulau Langkawi, Subang and Kuantan stations for the daily PET as well as Kuantan station for the annual PET time series. According to Hui-Mean et al. (2018), the increasing trends of the PET could be influenced by rising temperature that was attributed to the effect of global warming. The increasing trends of PET may cause water deficiency and subsequently lead to drought. In contrast to that, the decreasing trends of PET at Bayan Lepas station could be due to its geographical location since it is situated close to the coastal region. Based on the statistical characteristics of meteorological data by Binti Ahmad et al. (2019), the assumption is that high relative humidity may cause additional soil moisture to be present and result in the increase of water availability. However, further investigation on how well relative humidity is correlated to the trends of PET is still required. A comprehensive understanding of the trends of PET is essential for better planning and management of water resources especially for a tropical country with huge biodiversity.

The stationarity of PET for the time series was also checked using the ADF and KPSS tests. The stationarity tests were adopted to identify the statistical characteristics of the PET time series over a varying period. For the ADF test, a unit root was present at all stations for all PET time series analysed except at KLIA Sepang station for the monthly PET time series and Kuantan Station for the annual PET time series. Among the seven stations studied using the KPSS test, the results showed that only Bayan Lepas, Ipoh, Subang and Muadzam Shah stations were found to be stationary for all the three daily, monthly and annual PET time series. It can be concluded that the non-stationarity in PET time series may be influenced by extreme rainfall events or even other meteorological variables such as temperature and relative humidity.

One of the suggested recommendations is to carry out further studies with a longer length of input data (more than ten years) and check if this could improve the accuracy of results. Besides that, exploration of other PET estimation methods such as the Priestley–Taylor and Hargreaves–Samani methods are suggested for comparison as well as the use of the Innovative Trend Analysis (ITA) method for trend analysis and Phillips–Perron for the stationarity test. The adoption of different methods with different characteristics may be able to improve the accuracy of the results. Moreover, the empirical methods that were applied such as the Mann–Kendall (MK) test, Sen's slope estimator, Augmented Dickey–Fuller (ADF) test and Kwiatkowski–Phillips–Schmidt–Shin (KPSS) test could have certain limitations and their results are not always consistent.

All relevant data are included in the paper or its Supplementary Information.