Introduction

The potential evapotranspiration is considered as an important element of the hydrological cycle, which plays an important role in agricultural studies, management plans of irrigation and drainage networks, and water structures. (Gundekar et al. 2008; López-Urrea et al. 2006; Shayannejad et al. 2022; Ostad-Ali-Askari 2022; Vanani et al. 2017; Snyder et al. 2005). One of the methods of reducing farm water losses is to estimate accurately the water requirement of plants. According to FAO standard, reference evapotranspiration is the amount of water consumed by a farm with reference crop cover (such as grass) during a determined period in order to prevent plants to suffer from water shortages during the growth period (Allen et al. 1998). Therefore, it is important to plan for the optimal use of water resources in order to achieve the sustainable development. Estimating potential evapotranspiration by using different time scales in particular climatic regions is considered as one of the most important atmospheric parameters. Therefore, it is highly important in the process of optimal usage of resources. The time series analysis method, GARCH model, is applied in order to investigate changes and estimate the potential evapotranspiration.

Constructing statistical models and extending them in order to be utilized for various processes has always been considered by experts and engineers due to its complexity and inexistence of sufficient knowledge about physical processes of hydrological cycle. However, the basis of making many decisions related to hydrological processes and water resources’ exploitation is time series’ estimation and analysis. Nowadays, applying time series is considered as a suitable tool to carry out various estimations. Time series is considered as a set of observations arranged based on time. If these observations are recorded or measured at regular and equal intervals, then a discrete time series will be obtained (Weigend, 2018). It is possible to carry out synthetic hydrological data generation, hydrological events forecasting, trend recognition and jump in data, and completion of statistical gaps by using time series models. Time series models are consisted of two main components including random and model components, which can be obtained through observational data and random variables by using various stochastic methods. Therefore, the structure of time models will be in accordance with structure of the hydrological series if there are correct selection and calculations (Salas, 1993). Different statistical models including Auto Regressive (AR) models, Moving Average (MA), Auto Regressive Moving Average (ARMA), and Auto Regressive Integrated Moving Average (ARIMA) contain a set of models with various parameters; however, they can be applied as possible choices of modeling. If a single random process has the mean of zero, a specified variance, and a multiple regression model, then it will be considered as an autoregressive process. Autoregressive processes of the moving average are time series with a stationary process, which is represented by ARMA (p, q), in which p shows the autoregressive order and q is the moving average. Estimating the potential evapotranspiration is considered as a tool which is utilized for the process of synthetic data generation. Many researches are carried out in the field of applying time series models in order to estimate evapotranspiration of reference crops, and other meteorological and hydrological phenomena. (Landeras et al. 2009) estimated the weekly evapotranspiration with Hargreaves-Samani equation by using ARIMA model in northern Spain. Then, they compared these models with estimations resulted from Artificial Neural Network (ANN). It was found that ARIMA models’ performance during September to November is better than ANN. Silva et al. (2010) and Tian and Martinez (2012) showed that weather forecasting by using numerical methods can be applied as an appropriate approach in order to estimate future evapotranspiration. Psilovikos and Elhag (2013) estimated daily evapotranspiration in Nile River Delta by utilizing seasonal ARIMA models, and presented an appropriate model in order to study the mentioned region. Luo et al. (2014) presented a model which was calibrated by Hargreaves model, and provided air temperature forecasting data for the future seven days. It was found that results of proposed model were highly in accordance with results of Penman–Monteith model. Also, the above model could estimate reference crop’s evapotranspiration with an absolute error of 1 mm per a day. Many researches in the field of estimating reference crop’s evapotranspiration (ET0) are also carried out in other areas of the world including researches conducted by Cohen et al. (2002); Hulme et al. (1994); Szilagyi (2001) and Thomas (2000).

Performing the process of estimating evapotranspiration during various periods is highly required due to the lack of future meteorological information in order to provide plans for water resources and manage farm’s irrigation. In the present study, the efficiency of GARCH series model related to processes of modeling and estimating potential evapotranspiration, which is estimated by FAO Penman–Monteith and Hargreaves methods, is investigated. Also, future values of potential evapotranspiration are modeled and estimated at the synoptic station of Tabriz. Also evaluating GARCH model’s efficiency in processes of modeling and potential evapotranspiration estimation is also considered in this study.

Materials and methods

Considered reference evapotranspiration equations are utilized in order to assess GARCH model’s efficiency related to two equations of FAO Penman–Monteith and Hargreaves. These two equations are defined as the following.

FAO Penman–Monteith equation

In this method, the reference crop is a hypothetical grass covering with the height and reflection coefficient of 12 cm and 23%, respectively. It should be noted that plant resistance is the constant value of 70 s per meter. The following formulation is used in order to estimate the evapotranspiration of reference crop by using FAO Penman–Monteith method (Eq. 1) (Allen et al. 1998).

$${\text{ET}}_{{\text{o}}} = \frac{{0.408\Delta \left( {Rn - G} \right) + \gamma \left[ {900/\left( {T + 273} \right)} \right]U2\left( {ea - ed} \right)}}{{\Delta + \gamma \left( {1 + 0.34U2} \right)}}$$
(1)

where ETo: reference evapotranspiration [mm day−1], Rn: net radiation at the crop surface [MJ m−2 day−1], G: soil heat flux density [MJ m−2 day−1], T: mean daily air temperature at 2 m height [°C], u2: wind speed at 2 m height [m s−1], es: saturation vapor pressure [kPa], ea: actual vapor pressure [kPa], es—ea: saturation vapor pressure deficit [kPa], ∆: slope vapor pressure curve [kPa °C−1], γ: psychrometric constant [kPa °C−1].

Hargreaves equation

Hargreaves equation, which is with maximum and minimum temperatures, is utilized in order to calculate evapotranspiration during 24-h, weekly, 10-day, and monthly periods. The amount of reference crop’s evapotranspiration is also obtained by using Hargreaves method based on the following formulation as Eq. 2 (Hargreaves and Samani 1985):

$$ET_{0} = 0.0023 \times \left( {T_{mean} + 17.8} \right) \times \left( {T_{\max } - T_{\min } } \right)^{0/5} \times R_{a}$$
(2)

where ET0: reference evapotranspiration [mm day−1], Tmean: average daily air temperature [°C], Tmax: maximum daily air temperature [°C], Tmin: minimum daily air temperature [°C], Ra: radiation input at the top of the atmosphere [mm day−1].

Applied data

The amount of reference crop’s evapotranspiration is computed through FAO Penman–Monteith and Hargreaves methods by using daily data related to the airport synoptic station and Tabriz’s atmospheric conditions during 1993–2012. The maximum value per each month of the year was then determined as the reference crop’s potential value related to that month.

Introducing GARCH model

ARCH Model is presented by Engle (1982). According to this model, the random sentence has the average of zero and is serially non-correlated, but its variance will be variable if it has its past information. In this case, it is expected that the variance is not constant during the random sequence and follows the behavior of error sentences. In fact, ARCH model can explain the process of conditional variance according to its past information. Bollerslev (1986) suggested GARCH model, which is generalized to ARCH, for the first time as a solution for ARCH method’s problem. Also, GARCH restricts ARCH’s effect by using geometrical reduce of parameters’ number. This means that the effect of a shock on current fluctuations decreases over time. GARCH model’s (p, q) structure can be written as Eqs. 3 and 4:

$$\begin{aligned} (y_{t} |\varepsilon_{t} x_{t} + \varepsilon_{t - 1} ) = a^{ \circ } \sum\limits_{i = 1}^{s} {a_{i} y_{t - i} + \gamma } \hfill \\ (\varepsilon_{t} + \varepsilon_{t - 1} ) = N(0,\sigma_{t}^{2} ) \hfill \\ \end{aligned}$$
(3)
$$\sigma^{2} = E(\varepsilon_{t}^{2} |\varepsilon_{t - 1} ) = \beta^{ \circ } + \sum\limits_{i = 1}^{q} {\beta_{i} + \varepsilon_{t - i}^{2} } + \sum\limits_{j = 1}^{q} {\theta_{j} \sigma_{t - j}^{2} } + v_{t} \beta^{ \circ } \ge 0,\beta_{i} \ge 0,\theta_{j} \ge 0,\sum\limits_{i = 1}^{q} {\beta_{i} } + \sum\limits_{j = 1}^{q} {\theta_{j} < 1} ,\tilde{v}_{t} = IIN(0,\sigma_{v}^{2} )$$
(4)

where yt: Dependent variable in period t, xt: Explanatory variable in period t, εt Remaining in the period t, σ2: A conditional variance that is interpreted to predict the time series fluctuations in t period, εt-1: Includes a set of information up to time (t–1), in addition to εt. Equation (3), which is a criterion to determine the conditional mean, is considered as a function of exogenous variables with disruption component of (εt). If (εt) noise follows a normal distribution with the average of zero and conditional variance of (σ2t) in this equation, then it will be possible to apply Eq. (4). Equation (4) specifies the conditional variance. The equation of conditional variance consists of three parts of mean oscillation (β), (εt-i), and GARCH component (σ2t-j); where ARCH is the index of the previous period’s news and appears as the second power of conditional equation’s residue. Also, GARCH component represents the estimation of previous periods’ fluctuations. It should be noted that conditional variance’s residue existed in Eq. (4) has a normal distribution with the variance of (σ2). In other words, residue will be the white noise here. This condition is in accordance with all variance models of conditional heterogeneity.

Modeling processes using GARCH method

There were four modeling steps with the purpose of estimating evapotranspiration:

According to basic assumptions of time series modeling (static), the static nature of data series has to be ensured during the first step. Augmented Dickey–Fuller (ADF) method was applied in this research by using “Eviews” software in order to carry out the static test. “Eviews” is known as the application software utilized in the field of economics, which is designed and presented by the World Bank. This software is used in hydrologic time series in order to investigate the static; however, there are two common methods of Augmented Dickey Fuller (ADF) and Kwiatkowski Phillips Schmidt Shin (KPSS) for stationary hydrologic time series test. ADF method is a root test, which is developed by Dickey and Fuller (1979) and corrected by Said and Dickey (1984).

AR (1) model, which is a correlated model with a delay of 1, is as the following Eq. 5:

$$Z_{t} = \varphi_{1} Z_{t - 1} + \varepsilon_{t}$$
(5)

where “φ1” is the model’s coefficient, while “εt” is the standard random series. It can be said according to ADF method that if the absolute value of “φ1” is less than one, then time series will be static; however, time series will be astatic if the mentioned value is equal to one. The second step is determining the type or model, and selecting the best model. It should be noted that Autocorrelation Function (ACF), Partial Autocorrelation Function (PACF), or Akaike Information Criterion (AIC test) methods are applied for this purpose. The curve of partial autocorrelation function of different delays will be plotted after calculating partial autocorrelation function; however, the last consecutive point, which is located out of the probability range, will be determined by plotting the confidence interval. A line vertically lies on X-axis from the mentioned point, and related delay will be read. This delay determines the model’s order.

Statistics are produced during the third stage after selecting the best model. Model’s random parameters were calculated by “Excel” software in this research. Therefore, artificial data can be generated according to the shape by having model’s coefficients and normal random data. Determining the model error value and Goodness of Fit Test are carried out in the fourth step. In this study, the mean error of root squares was used in order to examine model’s error as well as the amount of explanatory factor in order to express model’s accuracy.

Evaluation indicators

The statistical mean of root mean square error (RMSE), absolute error percentage (AE), and coefficient of determination (R2) were utilized in this study (Eqs. 6–8):

$$RMSE = \sqrt {\left( {\frac{1}{n}\sum\limits_{i = 1}^{n} {\left( {O_{i} - P_{i} } \right)^{2} } } \right)}$$
(6)
$$AE = \frac{100}{n}\sum\limits_{i = 1}^{n} {\frac{{\left| {O_{i} - P_{i} } \right|}}{{O_{i} }}}$$
(7)
$$R^{2} = \frac{{\left[ {\sum\nolimits_{i = 1}^{n} {(P_{i} - P_{ave} )(O_{i} - O_{ave} )} } \right]^{2} }}{{\sum\nolimits_{i = 1}^{n} {(P_{i} - P_{ave} )^{2} } \cdot \sum\nolimits_{i = 1}^{n} {(O_{i} - O_{ave} )^{2} } }}$$
(8)

where in “Pi” is the evapotranspiration value, “Oi” is the observational evapotranspiration value, “Oava” is the average of modeled evapotranspiration values, “Pava” is the average of observational evapotranspiration values and finally, “n” is the number of investigated months.

The root means square of error is always positive; also, the best mode of operation is when this value is close to zero. The coefficient of determination is an index that shows the linearity of relationship between measured values; however, it should be noted that the linear relationship will be more obvious if mentioned value is closer to one. To have a reliable estimation of evapotranspiration, the absolute value of this statistic (for instance, the absolute error statistics (AE)), has to be applied. Also, model’s accuracy will increase by decrease of this statistic.

Results and Discussions

Monthly values of potential evapotranspiration were calculated by using two methods of FAO Penman- Monteith and Hargreaves. Results showed that FAO Penman–Monteith method estimated higher values of evapotranspiration compared to Hargreaves method based on the root mean square error and absolute error (respectively, 1.07 mm/day and 20%).

Model estimation and hypothesis testing

Data graphical observation is considered as the first step of time series analysis. The process of evapotranspiration changes of FAO Penman–Monteith and Hargreaves methods are shown in Figs. 1 and 2. Uncertainty has to be quantitated by using the conditional variance heterogeneity model (GARCH) in order to estimate the uncertainty effect resulted from volatilities. First, the reliability of evapotranspiration’s time series has to be investigated. Results of the reliability test, which is generalized by Generalized Dickey-Fuller test (ADF), indicate that the level of evapotranspiration stability obtained by FAO Penman-Montith method is determined to be 90%, while the level of stability of evaporation’s time series obtained by Hargreaves method is determined to be 95%. It should be noted that evaporation’s time series is investigated at first.

Fig. 1
figure 1

Evolution-Evapotranspiration Potential from the first month of the year 1993 to the last month of 2012 by the FAO Penman–Monteith method

Full size image
Fig. 2
figure 2

Evolution-Evapotranspiration Potential from the first month of the year 1993 to the last month of 2012 by the Hargreaves method

Full size image

Selecting an appropriate model to form the mean equation of evaporation according to two following methods is the main purpose of the second step. ARMA (p, q) model is selected for the mean equation of series by using the correlation graph, “R2”, Akaike-Schwartz criterion, Akaike information criterion (AIC), and Bayesian information criterion (BIC). According to obtained information, ARMA (2, 3) with BIC = 2.81, AIC = 2.73, and R2 = 80% was selected for FAO Penman–Monteith method. Also, ARMA with BIC = 1.99, AIC = 1.92, and R2 = 92% was selected based on acquired information. After selecting an appropriate pattern for the conditional mean equation and considering minimum values of AIC and BIC criteria, GRACH (1.1) model was considered as the optimal model of FAO Penman–Monteith AND Hargreaves methods. Results of the fitness of GARCH (1.1) model related to two series are summarized in the following Tables 1, 2, 3 and 4.

Table 1 Conditional mean equation evapotranspiration by the FAO Penman–Monteith method
Full size table
Table 2 Equation of conditional variance of evapotranspiration by the FAO Penman–Monteith method
Full size table
Table 3 Conditional mean equation evapotranspiration by the Hargreaves method
Full size table
Table 4 Equation of conditional variance of evapotranspiration by the Hargreaves method
Full size table

The square of standardized residues indicates that there is no sequential autocorrelation in standardized residues, which means the absence of ARCH effects on residues. However, it indicates the proper fitting of the GARCH model for both FAO-Penman-Monteith and Hargreaves methods (Table 5 and 6).

Table 5 Correlation charts (AC) and partial autocorrelation function (PAC) of the FAO Penman–Monteith residual series with 20 interruptions
Full size table
Table 6 Correlation charts (AC) and partial autocorrelation function (PAC) of the Hargreaves residual series with 20 interruptions
Full size table

In the third stage, estimating evapotranspiration was studied by using selected models. However, estimating evapotranspiration of both methods is investigated by applying GARCH (1.1) model (Figs. 3 and 4).

Fig. 3
figure 3

Observed and predicted evapotranspiration charts (FAO Penman Monteith method)

Full size image
Fig. 4
figure 4

Observed and predicted evapotranspiration charts (Hargreaves method)

Full size image

As it can be observed by both methods’ diagrams, observed evapotranspiration are highly close to the lower limit of predicted evapotranspiration. After performing the process of calibrating results obtained from the observatory evapotranspiration evaluation that was carried out by FAO Penman–Monteith, Hogwarts, and simulated methods, it was found that coefficients of determination (R2) are determined to be 0.84 and 0.96, while the root mean square errors (RMSE) are 0.76743 and 5097/0 mm per a day, respectively.

Conclusions

Results showed that Time Series is considered as a precise tool to estimate evapotranspiration values. It was found that GARCH (1.1) time series has better results for FAO Penman-Monteith and Hargreaves methods compared to other models. Also, it simulates the process of time series changes with less error. After performing the process of calibrating results of observatory evapotranspiration evaluation that was carried out by FAO Penman-Monteith, Hargreaves, and simulated methods, it was found that estimated GARCH model had better results than FAO Penman-Monteith method and finally, simulates the process of time series’ changes with less error. Observed and predicted evapotranspiration charts of both methods indicated that observational evapotranspiration was highly close to the lower limit of estimated evapotranspiration. Therefore, applying lower limit estimation as a prediction value was suggested.