Abstract

Confronting drought and reducing its impacts requires modeling and forecasting of this phenomenon. In this research, the ability of different time series models (the ARIMA models with different structures) were evaluated to model and predict seasonal drought based on the RDI drought index in the south of Iran. For this purpose, the climatic data of 16 synoptic stations from 1980 to 2010 were used. Evaluation of time series models was based on trial and error. Results showed drought classes varied between ‘very wet’ to ‘severely dry’. The more occurrence frequency of ‘severely dry’ class compared to other drought classes represent the necessity of drought assessment and the importance of managing the effects of this phenomenon in the study area. Results showed that the highest severity of drought occurred at Abadeh, Shiraz, Fasa, Sirjan, Kerman, Shahre Babak and Saravan stations. According to selecting the best model fitted to the computed three-month RDI time series, results indicated that the MA model based on the Innovations method resulted in maximum cases with the best performance (37.5% of cases). The AR model based on the Yule–Walker method resulted in minimum cases with the best performance (6.3% of cases) in seasonal drought forecasting.

INTRODUCTION

Early indications of possible drought can help to plan drought mitigation strategies and measures. Modeling and predicting the behavior of the phenomenon in a standard time period is the best way to address damaging this natural phenomenon. Drought damages environmental factors and agriculture, vegetation, humans and wildlife, as well as local economies (Azarakhshi et al. 2011; Dastorani & Afkhami 2011; Nohegar et al. 2013; Zarei et al. 2013). This phenomenon begins naturally and ends slowly and its damage will be very serious and long-lasting, especially in arid and semi-arid regions. It can result in long-term harmful effects on water resources and the environment (Zarei et al. 2016; Zarei 2018).

In recent decades, many techniques have been used as suitable tools for modeling and forecasting meteorological information such as drought (Soltani et al. 2007; Shamshirband et al. 2015). Chun et al. (2013) evaluated and predicted the impact of climate change on drought in the UK using ARIMA models and the generalized linear model (GLM) approach. Results indicated that the drought pattern in the 2080s is less certain than for the 1961–1990 period, based on the Shannon entropy, but droughts are expected to be more clustered and intermittent. Jahandideh & Shirvani (2011) used time series models in Fars Province to forecast drought based on the standardized precipitation index. According to the results of this research, the SARIMA model with the minimum of Akaike's information criterion bias-corrected (AICC) was selected as the best model.

Mossad & Ali Alazba (2015) developed several ARIMA models for drought forecasting in a hyper-arid climate using the SPEI index. The results reveal that all developed ARIMA models demonstrate the potential ability to forecast drought over different time scales. This study recommends that ARIMA models can be very useful tools for drought forecasting. Djerbouai & Souag-Gamane (2016) investigated the ANN models for drought forecasting in the Algerois basin in Algeria in comparison with traditional stochastic models (ARIMA and SARIMA models). Results of this study showed that despite the linear property of SARIMA models, they yielded satisfactory performances, with respect to various model efficiencies, for SPI-12 forecasts at one-month lead time. Paul et al. (2017) used time series models (robust statistics) to assess the trend of rainfall data at Rajahmundry city located in lower Godavari basin, India. Phuong et al. (2018) evaluated the spatiotemporal variability of annual and seasonal rainfall time series in Ho Chi Minh City for the period 1980–2016. The results of trend estimation also indicated higher increasing rates of rainfall in the dry season compared to the rainy season at most stations. Other researchers have also used time series models to assess drought characteristics and conditions of climate parameters (Barua et al. 2010; Hong-yan et al. 2015; Wang et al. 2016; Zarei & Moghimi 2017; Bahrami et al. 2018; Zarei 2018).

Proper management of drought requires more precise monitoring, modeling and forecasting of this phenomenon using strong drought indices. In this research, the reconnaissance drought index (RDI) was used due to a number of advantages over other drought indices and because it can be used to quantify most types of drought events. This index was proposed by Tsakiris & Vangelis (2005) and for calculation of this index the ratio of cumulative values of precipitation to potential evapotranspiration was used. In addition, this index has certain advantages compared to indices based on precipitation, because it is more representative of the deficient water balance conditions. Since RDI resolved more climatic parameters, such as evapotranspiration, which had an important role in water resource losses in the Iranian basins, it was worthwhile to consider RDI in drought monitoring in Iran (Asadi Zarch et al. 2011). Therefore, for the regions with no sufficient data to calculate other indicators, the RDI index is the best and most comprehensive index.

The aims of this research are modeling and predicting of seasonal drought (calculated based on the RDI index, as a strong index for monitoring drought) in southern Iran using different time series models (ARIMA models). There is an integrated system from the model establishment, verification, and forecasting.

MATERIAL AND METHODS

Study area

Figure 1 shows the study area boundaries enclosing a 594,996.54 km2 area approximately between 25°17′N and 31°11′N latitudes and between 50°49′E and 62°20′E longitudes, in addition to the distribution of 16 synoptic stations used in this study area. These stations have a good distribution with sufficient length of meteorological data (31 years), which can spatially and temporally support drought modeling and forecasting studies. The central and northern areas are highlands and mountains, while the southern and western areas are mainly flat. The elevation of selected stations varied from 5 (at Jask) to 2,030 m (at Abadeh) from free sea surface level. The average annual precipitation varied from 53.1 (at Zabol) to 330.6 mm (at Shiraz) and the average annual potential evapotranspiration varied from 1,416.2 (at Jask) to 2,839.7 mm (at Bushehr). More detailed characteristics of the 16 surveyed stations are presented in Table 1.

Table 1

General characteristic of surveyed 16 stations

Station Elevation from free sea surface level (m) Average of annual precipitation (mm) Average of annual PET (mm) Aridity index Climate 
Abadeh 2,030 133.44 1,752 0.076 Arid 
Shiraz 1,484 330.57 1,985.6 0.169 Arid 
Fasa 1,288 284.99 1,978.3 0.156 Arid 
Sirjan 1,739 134.77 1,602.35 0.092 Arid 
Kerman 1,754 130.13 1,306.7 0.074 Arid 
Shahre Babak 1,834 145.69 1,346.85 0.103 Arid 
Bam 1,067 54.32 1,832.3 0.027 Hyper-arid 
Bushehr 262.88 2,839.7 0.201 Semi-arid 
Bandar Abbas 10 160.18 1,474.6 0.081 Arid 
Bandar Lengeh 23 122.72 1,762.95 0.077 Arid 
Iranshahr 591 105.15 2,606.1 0.037 Hyper-arid 
Jask 118.64 1,416.2 0.08 Arid 
Chabahar 121.2 1,960.05 0.09 Arid 
Saravan 1,195 109.74 1,460 0.042 Hyper-arid 
Zabol 489 53.12 2,584.2 0.021 Hyper-arid 
Zahedan 1,370 75.7 1,788.5 0.042 Hyper-arid 
Station Elevation from free sea surface level (m) Average of annual precipitation (mm) Average of annual PET (mm) Aridity index Climate 
Abadeh 2,030 133.44 1,752 0.076 Arid 
Shiraz 1,484 330.57 1,985.6 0.169 Arid 
Fasa 1,288 284.99 1,978.3 0.156 Arid 
Sirjan 1,739 134.77 1,602.35 0.092 Arid 
Kerman 1,754 130.13 1,306.7 0.074 Arid 
Shahre Babak 1,834 145.69 1,346.85 0.103 Arid 
Bam 1,067 54.32 1,832.3 0.027 Hyper-arid 
Bushehr 262.88 2,839.7 0.201 Semi-arid 
Bandar Abbas 10 160.18 1,474.6 0.081 Arid 
Bandar Lengeh 23 122.72 1,762.95 0.077 Arid 
Iranshahr 591 105.15 2,606.1 0.037 Hyper-arid 
Jask 118.64 1,416.2 0.08 Arid 
Chabahar 121.2 1,960.05 0.09 Arid 
Saravan 1,195 109.74 1,460 0.042 Hyper-arid 
Zabol 489 53.12 2,584.2 0.021 Hyper-arid 
Zahedan 1,370 75.7 1,788.5 0.042 Hyper-arid 
Figure 1

Geographic position and digital elevation model of the study area and spatial distribution of selected synoptic meteorological stations.

Figure 1

Geographic position and digital elevation model of the study area and spatial distribution of selected synoptic meteorological stations.

Method

Data collection

To assess the drought characteristics based on the RDI index, meteorological data of 16 synoptic stations (www.irrimo.ir) were used. Variability in climate conditions of synoptic stations was one of the main criteria in selecting the synoptic stations. The climatic conditions of selected stations was calculated using the UNEP aridity index (UNEP 1992). 
formula
(1)
where P is the average of annual precipitation (mm) and PET is the average of annual potential evapotranspiration. Table 2 was used for the classification of the degree of climate dryness according to the UNEP index. For PET calculation, the FAO-56 Penman–Monteith (FAO-56 PM) equation (Allen et al. 1998) was used as follows: 
formula
(2)
where PET is potential evapotranspiration, Δ is the slope of the saturation vapor pressure function, Rn is net radiation, G is soil heat flux density, is psychometric constant, T is mean air temperature, U2 is average 24-h wind speed at 2 m height and VPD is vapor pressure deficit.
Table 2

Classification of climate according to UNEP aridity index (UNEP 1992)

P/PET1 Climate 
Lower than 0.05 Hyper-arid 
0.05–0.20 Arid 
0.2–0.5 Semi-arid 
0.5–0.65 Sub-humid 
Greater than 0.65 Humid 
P/PET1 Climate 
Lower than 0.05 Hyper-arid 
0.05–0.20 Arid 
0.2–0.5 Semi-arid 
0.5–0.65 Sub-humid 
Greater than 0.65 Humid 

P: Average of annual precipitation (mm) and PET: Average of annual potential evapotranspiration.

Since the RDI considers the proportion of precipitation to PET, at first, precipitation data and the data that was required to compute PET (maximum and minimum temperature, maximum and minimum relative humidity, wind speed, and sunshine) were gathered, and the missing data was substituted with the corresponding long-term mean (Dinpashoh et al. 2011). ET0 was estimated using the PM-56 equation (Allen et al. 1998) at the monthly time scale for each synoptic station.

The RDI drought index

The reconnaissance drought identification and assessment index was used in this research (RDI) because of the possible role of ET0 in the detection of drought events. This index was proposed by Tsakiris & Vangelis (2005), utilizing the ratios of precipitation over Reference Crop Evapotranspiration (ET0) for different time scales to be representative of the desired region. Calculating RDI was carried out using equations as follows:

At the first step, αk is calculated as the coefficient of the ith year in aggregated form as follows, using a monthly time step: 
formula
(3)
where Pij and PETij are precipitation and potential evapotranspiration of the jth month of the ith year that usually starts from October (K = 1) in the study region. N is the total number of years of the available data. Equation (3) can be calculated for any period of the year. It could also be recorded starting from any month of the year except October. As previously mentioned, ET0 was used to represent PETij, which was estimated using the Penman–Monteith equation (Allen et al. 1998). The stations selected in this study are the official meteorological synoptic stations in each province, with complete data required for estimating ET0 using the Penman–Monteith method. Computation of normalized RDI (RDIn) was the next step, using as the arithmetic mean of αk values: 
formula
(4)
Finally, the standardized RDI (RDIst) is computed as follows: 
formula
(5)
where yk is the , is the arithmetic mean of yk and is the standard deviation. According to several studies on various data from several locations and different time scales, it was concluded that αk values follow the gamma distributions, which has been found to be more successful and it has been proven that the calculation of RDIst could be performed better by fitting the gamma probability density function (pdf) at the given frequency distribution of αk, following the method described below (Tsakiris et al. 2008; Asadi Zarch et al. 2011; Kousari et al. 2014). This method tends to solve the problem of calculating RDIst for the small time scales, such as monthly, which may include zero precipitation values (αk = 0), for which Equation (4) could not be applied (Tsakiris et al. 2008). The gamma distribution is defined by its frequency or probability density function: 
formula
(6)
where α and β are shape and scale parameters, x is the precipitation amount and is the gamma function. Maximum likelihood solutions were used to estimate α and β. 
formula
(7)
 
formula
(8)
 
formula
(9)
In some cases, precipitation distribution contains zeros and the gamma function is undefined for x = 0, thus, the cumulative probability for x = 0 becomes: 
formula
(10)
where q is the probability of zero precipitation and G(x) is the cumulative probability of the incomplete gamma function. If m is the number of zeros in an αk time scale/s, then q could be estimated by m/n. The cumulative probability H(x) is then transformed to the standard normal random variable z with mean zero and the variance of one (Abramowitz & Stegun 1965), which is the value of RDIst (Tsakiris et al. 2008). In the present study, three-month RDIst was used to represent three-month RDI. Drought category classification suggested for the RDI is illustrated in Table 3 (Tsakiris 2004; Tsakiris et al. 2007; Banimahd & Khalili 2013; Zarei et al. 2016).
Table 3

Drought classification of RDI (Asadi Zarch et al. 2011; Zarei et al. 2016)

Drought class RDI value 
Extremely wet RDI ≥ 2 
Very wet 1.5 < RDI < 1.99 
Moderately wet 1 < RDI < 1.49 
Normal 0 < RDI < 0.99 
Near normal –0.99 < RDI < 0 
Moderately dry –1.49 < RDI < –1 
Severely dry –1.99 < SPI < –1.5 
Extremely dry RDI ≤ –2 
Drought class RDI value 
Extremely wet RDI ≥ 2 
Very wet 1.5 < RDI < 1.99 
Moderately wet 1 < RDI < 1.49 
Normal 0 < RDI < 0.99 
Near normal –0.99 < RDI < 0 
Moderately dry –1.49 < RDI < –1 
Severely dry –1.99 < SPI < –1.5 
Extremely dry RDI ≤ –2 

Modeling and forecasting

Stochastic models

Usually, the stochastic models known as time series models (ARIMA) have been used in scientific applications for the analysis of time series. Autoregressive integrated moving average (ARIMA) models are mathematical models of persistence, or autocorrelation, in a time series introduced by Box & Jenkins (1976). These models are expressed by a series of equations. An AR (autoregressive) model is a subset of ARIMA models and describes a time series as a linear function of its past values plus a noise term . The order of the AR model shows the number of past values involved. The equation of the simplest AR model (the first-order autoregressive (AR (1))) is given by: 
formula
(11)
where is a stationary mean zero time series and is the first-order autoregressive coefficient.
The moving average (MA) model is another form of ARIMA model in which the time series is described as a linear function of its prior errors plus a current error . The equation of first-order moving average (MA (1)) is given by: 
formula
(12)
where is a stationary mean zero time series, and are the error terms at time t and t − 1, and is the coefficient of the first-order moving average.
A general autoregressive moving average (ARMA) model, ARMA (p, q), is given by: 
formula
(13)

The ARIMA model (the integrated ARMA) is a broadening of the class of ARMA that includes differencing (an important technique in data transformation; it attempts to de-trend to control autocorrelation and achieve stationary time series).

The first order differencing is denoted by the following equation: 
formula
(14)
where B is the backshift operator. The differencing of order d is defined by: 
formula
(15)
Usually, the linear trend is removed by single differencing and the quadratic trend is removed by double differencing. Seasonality and trend of period d can be removed by introducing the lag-d differencing operator : 
formula
(16)

ARIMA modeling generally involves three stages as follows. First stage: model identification by specifying the type of the model (AR, MA, ARMA, or ARIMA) and its order. This identification is sometimes undertaken by looking at plots of the sample autocorrelation function (ACF) and sample partial autocorrelation function (PACF) and sometimes it is employed by an autofit procedure fitting many different possible model structures and orders and using a goodness-of-fit statistic to select the best model. Second stage: estimate the coefficients of the model by minimizing the sum of squared residuals. Third stage: model diagnostics. In this stage, it is very important to check that the residuals of the candidate model are random and normally distributed and the estimated parameters are statistically significant. It should be noted that the best model of all which fit the data is the one which has the fewest parameters.

Cross-validation and transformation

In the cross-validation procedure, the data set split into two sets: training sample and prediction. The training sample is used to develop a model for prediction and the prediction set is used to evaluate the rationality and predictive ability of the selected model. This validation procedure is the statistical practice of splitting a sample of data into two subsets so that the analysis is initially performed on one subset and the other subset is retained for subsequent use in confirming and validating the initial analysis. For fitting the ARMA model, it must be at least likely that the data are in fact a realization of an ARMA process and, in particular, a realization of a stationary process. In the stationary time series, the statistical properties such as mean, variance, autocorrelation and so on are all constant over time. In order to obtain a stationary time series, a sequence of mathematical transformations (Box–Cox transformation, mean subtraction, and the differencing) was used.

The Box–Cox transformation , for a sequence of observations (), is given by the following equation: 
formula
(17)

When the variability of the data increases or decreases with the level, this transformation is useful. The variability can be made nearly constant by a suitable choice of . For instance, the variability of a set of positive data whose standard deviation increases linearly can be stabilized by choosing (Brockwell & Davis 2002).

The lowest order of differencing is the correct amount of differencing that resulted in time series which fluctuate around a well-defined mean value and whose ACF plot decays rapidly to zero, either from above or below. Thus, at every stage of differencing, the plots of the sample ACF and the sample PACF were checked to see where the ACF/PACF were out of the bounds ±1.96/. The stationary series is the series with a sample ACF that decays fairly rapidly. If the ACF plot has a polynomial trend, it shows that the series still has some trends. The periodicity of ACF shows that the series has seasonality and some more differencing for the data should be applied.

Model selection

The residual ACF/PACF of the models and the randomness of the residuals should be checked. For the sample with large n, the sample autocorrelations of an independent and identically distributed (iid) sequence () with zero mean and finite variance are approximately iid with normal distribution N (0,1/n) (Amei et al. 2012). The consistency of the observed residuals with iid noise was considered by examining the sample correlations of the residuals and rejecting the iid noise hypothesis if more than two or three out of 40 falls outside the bounds ±1.96/ or if one falls far outside the bounds (Brockwell & Davis 2002).

The Ljung–Box test (Ljung & Box 1978) was used to check whether the residuals of a fitted model are iid in ARIMA modeling or not. This test was based on the autocorrelation plot and tests the overall independence based on a few of the time lags. The Ljung–Box test is defined as follows:

  • H0: The sequence data are iid

  • Ha: The sequence data are not iid

The test statistic is , where is the estimated autocorrelation at lag-k and is equal to , n is sample size, m is the number of lags being tested and are the residuals after a model has been fitted to a series Z1, …, Zn. For large sample sizes (n), the distribution of is approximately under the null hypothesis, where p + q is the number of parameters of the fitted model. The hypothesis of iid is rejected if at level and the sequence data do have autocorrelations significantly different from zero and a new search for a fitted ARMA model for a mean-corrected data set should be followed.

Model comparison

In this research, the AICC statistic (Akaike 1974), the bias-corrected version of AIC statistic, as an information criterion to select candidate models using the ITSM2000 package (Brockwell & Davis 2002), was used. As a rough guide, the small value of AICC is an indication of a good model. Maximum likelihood estimation is the base of final decisions between models. ITSM2000 contains some other model-selection statistics, such as a BIC statistic. Bayesian modification of the AIC statistic is the BIC statistic (Schwarz 1978) which is evaluated at the same time and used in the same way as the AICC. The definition of each information statistic is as follows: 
formula
(18)
 
formula
(19)
 
formula
(20)
where is the maximum likelihood estimator of , and r is the number of parameters estimated in the model, including a constant term, that is equal to (p+q+ 1). In all of the latter equations, the second term is a penalty for increasing r. Therefore, the best model is the model that adequately describes data and has the fewest parameters.

The best ARIMA models that fitted to the data were used to forecast future values of the time series from the observed values.

Application

After computation of three-month RDI for 16 synoptic stations of the study area (southern Iran), ARIMA models were fitted to these data and then optimized by eliminating non-significant coefficients. Then, using various indicators, the goodness of fitted models was evaluated. For this purpose, heterogeneity and randomness of the residuals, model validity in the forecast, and comparison between ACF/PACF of data and the fitted model were considered.

Model validation

Evaluation of the forecast performance of all developed models was carried out using different measures of goodness of fit such as the Nash–Sutcliffe efficiency coefficient (NSE), root mean squared error (RMSE), mean absolute error (MAE) and correlation coefficient (R-squared value), calculated as follows: 
formula
(21)
 
formula
(22)
 
formula
(23)
 
formula
(24)

where RDIi is computed RDI for the subset of observed data (2008–2010), is the arithmetic mean of RDIi, is the forecasted RDI for the subset of observed data (2008–2010), and is the arithmetic mean of .

In general, high values of NSE (up to 100%) and correlation coefficient (R-squared value) and small values for RMSE and MAE indicate a good model.

RESULTS AND DISCUSSION

The results of drought class probabilities relative to the three-month RDI for 16 synoptic stations of this study region for observed data (time period of 1980–2010) are shown in Table 4. Generally, the drought classes varied between ‘very wet’ to ‘severely dry’, and the number of non-zero probabilities related to the ‘severely dry’ class is considerably higher than other drought classes that indicated the drought occurrence in this time period at this study area. In the Bam, Bandar Abbas, Bandar Lengeh, Jask, Chabahar and Zabol stations, drought class probabilities varied between ‘very wet’ to ‘moderately dry’, which indicated the lowest severity of drought in these stations. The Abadeh, Shiraz, Fasa, Sirjan, Kerman, Shahre Babak and Saravan stations showed the drought class probabilities between ‘moderately wet’ to ‘severely dry’. Therefore, in the latter stations, the highest severity of drought occurred. In Bushehr and Zahedan stations, drought class probabilities varied between ‘moderately wet’ and ‘moderately dry’ which indicated the normal conditions related to drought.

Table 4

Drought class probabilities relative to the three-monthly RDI

  Analytical steady class probabilities
 
Station Extremely wet Very wet Moderately wet Near normal Moderately dry Severely dry Extremely dry 
Abadeh 0.00 0.00 0.19 0.57 0.09 0.15 0.00 
Shiraz 0.00 0.00 0.23 0.53 0.08 0.15 0.00 
Fasa 0.00 0.00 0.23 0.57 0.07 0.13 0.00 
Sirjan 0.00 0.00 0.19 0.63 0.06 0.13 0.00 
Kerman 0.00 0.00 0.22 0.53 0.10 0.15 0.00 
Shahre Babak 0.00 0.00 0.19 0.59 0.11 0.11 0.00 
Bam 0.00 0.03 0.17 0.51 0.29 0.00 0.00 
Bushehr 0.00 0.00 0.27 0.45 0.27 0.00 0.00 
Bandar Abbas 0.00 0.02 0.23 0.52 0.24 0.00 0.00 
Bandar Lengeh 0.00 0.08 0.15 0.43 0.35 0.00 0.00 
Iranshahr 0.00 0.00 0.22 0.55 0.06 0.17 0.00 
Jask 0.00 0.12 0.13 0.39 0.36 0.00 0.00 
Chabahar 0.00 0.09 0.15 0.41 0.35 0.00 0.00 
Saravan 0.00 0.00 0.21 0.58 0.09 0.12 0.00 
Zabol 0.00 0.12 0.13 0.40 0.35 0.00 0.00 
Zahedan 0.00 0.00 0.23 0.47 0.30 0.00 0.00 
  Analytical steady class probabilities
 
Station Extremely wet Very wet Moderately wet Near normal Moderately dry Severely dry Extremely dry 
Abadeh 0.00 0.00 0.19 0.57 0.09 0.15 0.00 
Shiraz 0.00 0.00 0.23 0.53 0.08 0.15 0.00 
Fasa 0.00 0.00 0.23 0.57 0.07 0.13 0.00 
Sirjan 0.00 0.00 0.19 0.63 0.06 0.13 0.00 
Kerman 0.00 0.00 0.22 0.53 0.10 0.15 0.00 
Shahre Babak 0.00 0.00 0.19 0.59 0.11 0.11 0.00 
Bam 0.00 0.03 0.17 0.51 0.29 0.00 0.00 
Bushehr 0.00 0.00 0.27 0.45 0.27 0.00 0.00 
Bandar Abbas 0.00 0.02 0.23 0.52 0.24 0.00 0.00 
Bandar Lengeh 0.00 0.08 0.15 0.43 0.35 0.00 0.00 
Iranshahr 0.00 0.00 0.22 0.55 0.06 0.17 0.00 
Jask 0.00 0.12 0.13 0.39 0.36 0.00 0.00 
Chabahar 0.00 0.09 0.15 0.41 0.35 0.00 0.00 
Saravan 0.00 0.00 0.21 0.58 0.09 0.12 0.00 
Zabol 0.00 0.12 0.13 0.40 0.35 0.00 0.00 
Zahedan 0.00 0.00 0.23 0.47 0.30 0.00 0.00 

Modeling and forecasting

Initially, in order to obtain stationary time series, Box–Cox transformation (to achieve constant variance), mean subtraction, and the differencing (to eliminate trend (k = 1) and periodicity (d = 4)) were used. Figure 2 shows the computed three-month RDI time series for observed data (1980–2010) for the synoptic meteorological stations of Abadeh, Bushehr and Zabol (results of other stations are not shown due to space). Figure 3 shows the stationary time series for the mentioned stations.

Figure 2

Computed three-monthly RDI time series for observed data (1980–2010) for some of the synoptic meteorological stations (for example).

Figure 2

Computed three-monthly RDI time series for observed data (1980–2010) for some of the synoptic meteorological stations (for example).

Figure 3

Stationary three-monthly RDI time series for observed data (1980–2010) for some of the synoptic meteorological stations (for example).

Figure 3

Stationary three-monthly RDI time series for observed data (1980–2010) for some of the synoptic meteorological stations (for example).

Model diagnostics and fitting

Determination of the order of p and q in the AR, MA, and ARMA models was examined by using the ACF/PACF plots (Figure 4). According to Figure 4, up to 5% of residual ACF/PACF/s are out of zero range (the dotted lines) and the sample and model ACF/PACF are close together indicating the appropriateness of the model. After selecting the best values of p and q, different models were fitted to the computed three-month RDI time series for observed data (1980–2010) in all of the stations and, according to the AICC statistic, the model with the minimum value of AICC was selected as the best model. The results of model diagnostics and fitting shown in Table 5 indicate the best models with the minimum value of AICC at each station. Results indicated that the best performances in seasonal drought forecasting were as follows: the AR model based on the Yule–Walker method in 6.3% of cases; the AR model based on the Burg method in 12.5% of cases; the MA model based on the Hannan–Rissanen method in 25% of cases; the MA model based on the Innovations method in 37.5% of cases; the ARMA model based on the Hannan–Rissanen method in 25% of cases; and the ARMA model based on the Innovations method in 12.5% of cases.

Table 5

The best models after difference 4 and mean-correction in different stations in south of Iran

Station Best model Method AICC statistic BIC statistic Ljung–Box statistic P-value 
Abadeh AR(16) Yule-Walker 194.098 197.542 0.396 
Shiraz MA(4) Hannan-Rissanen (or Innovations) 170.349 135.424 0.134 
Fasa ARMA(8,12) Hannan-Rissanen 227.068 203.502 0.259 
Sirjan MA(18) Innovations 238.921 215.767 0.065 
Kerman MA(18) Hannan-Rissanen 191.927 181.489 0.826 
Shahre Babak MA(6) Hannan-Rissanen (or Innovations) 207.631 191.071 0.521 
Bam MA(8) Innovations 238.960 218.163 0.260 
Bushehr ARMA(1,4) Innovations 136.307 135.985 0.615 
Bandar Abbas ARMA(16,4) Hannan-Rissanen 230.502 219.198 0.342 
Bandar Lengeh MA(8) Hannan-Rissanen 227.956 202.856 0.402 
Iranshahr ARMA(5,4) Hannan-Rissanen 303.331 296.684 0.576 
Jask MA(8) Innovations 247.128 232.586 0.465 
Chabahar ARMA(1,4) Innovations (or Hannan-Rissanen) 289.153 285.894 0.570 
Saravan MA(25) Innovations 301.312 269.031 0.315 
Zabol AR(16) Burg 194.282 197.680 0.571 
Zahedan AR(24) Burg 252.766 254.663 0.446 
Station Best model Method AICC statistic BIC statistic Ljung–Box statistic P-value 
Abadeh AR(16) Yule-Walker 194.098 197.542 0.396 
Shiraz MA(4) Hannan-Rissanen (or Innovations) 170.349 135.424 0.134 
Fasa ARMA(8,12) Hannan-Rissanen 227.068 203.502 0.259 
Sirjan MA(18) Innovations 238.921 215.767 0.065 
Kerman MA(18) Hannan-Rissanen 191.927 181.489 0.826 
Shahre Babak MA(6) Hannan-Rissanen (or Innovations) 207.631 191.071 0.521 
Bam MA(8) Innovations 238.960 218.163 0.260 
Bushehr ARMA(1,4) Innovations 136.307 135.985 0.615 
Bandar Abbas ARMA(16,4) Hannan-Rissanen 230.502 219.198 0.342 
Bandar Lengeh MA(8) Hannan-Rissanen 227.956 202.856 0.402 
Iranshahr ARMA(5,4) Hannan-Rissanen 303.331 296.684 0.576 
Jask MA(8) Innovations 247.128 232.586 0.465 
Chabahar ARMA(1,4) Innovations (or Hannan-Rissanen) 289.153 285.894 0.570 
Saravan MA(25) Innovations 301.312 269.031 0.315 
Zabol AR(16) Burg 194.282 197.680 0.571 
Zahedan AR(24) Burg 252.766 254.663 0.446 
Figure 4

Plots of sample and model autocorrelation function (ACF) and the sample and model partial autocorrelation function (PACF) for some of the synoptic meteorological stations (for example).

Figure 4

Plots of sample and model autocorrelation function (ACF) and the sample and model partial autocorrelation function (PACF) for some of the synoptic meteorological stations (for example).

According to the results (Table 5), at the Abadeh, Zabol, and Zahedan stations, the AR model was diagnosed as the best model and the general form of this model is shown in the following equation: 
formula
(25)

The AR coefficients (Øi) are presented in Table 6.

Table 6

The AR coefficients (Øi) obtained for autoregressive model

Station AR coefficients
 
Ø1 Ø2 Ø3 Ø4 Ø5 Ø6 Ø7 Ø8 
Abadeh 0.0000 0.1151 0.0000 −0.7702 0.0000 0.0000 0.0000 −0.3291 
Zabol 0.0000 0.1560 0.0000 −0.8214 0.0000 0.0000 0.0000 −0.3164 
Zahedan 0.0000 0.0000 0.0000 −0.8305 0.0000 0.0000 0.0000 −0.5666 
 Ø9 Ø10 Ø11 Ø12 Ø13 Ø14 Ø15 Ø16 
Abadeh 0.0000 0.0000 0.0000 −0.2209 0.0000 0.0000 0.0000 −0.3342 
Zabol 0.0000 0.0000 0.0000 −0.4456 0.0000 0.0000 0.0000 −0.4289 
Zahedan 0.0000 0.0000 0.0000 −0.4663 0.0000 0.0000 0.0000 −0.3032 
 Ø17 Ø18 Ø19 Ø20 Ø21 Ø22 Ø23 Ø24 
Abadeh – – – – – – – – 
Zabol – – – – – – – – 
Zahedan 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 −0.2431 
Station AR coefficients
 
Ø1 Ø2 Ø3 Ø4 Ø5 Ø6 Ø7 Ø8 
Abadeh 0.0000 0.1151 0.0000 −0.7702 0.0000 0.0000 0.0000 −0.3291 
Zabol 0.0000 0.1560 0.0000 −0.8214 0.0000 0.0000 0.0000 −0.3164 
Zahedan 0.0000 0.0000 0.0000 −0.8305 0.0000 0.0000 0.0000 −0.5666 
 Ø9 Ø10 Ø11 Ø12 Ø13 Ø14 Ø15 Ø16 
Abadeh 0.0000 0.0000 0.0000 −0.2209 0.0000 0.0000 0.0000 −0.3342 
Zabol 0.0000 0.0000 0.0000 −0.4456 0.0000 0.0000 0.0000 −0.4289 
Zahedan 0.0000 0.0000 0.0000 −0.4663 0.0000 0.0000 0.0000 −0.3032 
 Ø17 Ø18 Ø19 Ø20 Ø21 Ø22 Ø23 Ø24 
Abadeh – – – – – – – – 
Zabol – – – – – – – – 
Zahedan 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 −0.2431 

The zero values indicate the non-significant coefficients.

The dashed lines indicate the absence of coefficient in model at the related station.

At Shiraz, Sirjan, Kerman, Shahre Babak, Bam, Bandar Lengeh, Jask and Saravan stations, the MA model was diagnosed as the best model and the general form of this model is shown in the following equation: 
formula
(26)
The MA coefficients (θi) are presented in Table 7.
Table 7

The MA coefficients (θi) obtained for moving average model

Station MA coefficients
 
θ1 θ2 θ3 θ4 θ5 θ6 θ7 θ8 θ9 
Shiraz 0.0000 0.0000 0.0000 –1.0000 – – – – – 
Sirjan 0.0000 0.0000 0.0000 –1.0740 0.0000 0.0000 0.0000 0.0000 0.0000 
Kerman 0.0000 0.0000 0.0000 –0.9647 0.0000 0.0000 0.0000 0.0000 0.0000 
Shahre Babak 0.0000 0.4378 0.0000 –0.9397 0.0000 –0.3769 – – – 
Bam 0.0000 0.2328 0.0000 –1.2899 0.0000 –0.2070 0.0000 0.3155 – 
Bandar Lengeh 0.0000 0.0000 0.0000 –1.2538 0.0000 0.0000 0.0000 0.2537 – 
Jask 0.0000 0.0000 0.0000 –1.2335 0.0000 0.0000 0.0000 0.3187 – 
Saravan 0.0000 0.0000 –0.1767 –1.0604 0.0000 0.0000 0.3871 0.0000 0.0000 
 θ10 θ11 θ12 θ13 θ14 θ15 θ16 θ17 θ18 
Shiraz – – – – – – – – – 
Sirjan 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 –0.1365 
Kerman 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 –0.0306 0.0000 –0.0310 
Shahre Babak – – – – – – – – – 
Bam – – – – – – – – – 
Bandar Lengeh – – – – – – – – – 
Jask – – – – – – – – – 
Saravan 0.0000 0.0000 0.0000 0.0000 –0.2883 0.0000 0.0000 0.0000 0.0000 
 θ19 θ20 θ21 θ22 θ23 θ24 θ25   
Shiraz – – – – – – –   
Sirjan – – – – – – –   
Kerman – – – – – – –   
Shahre Babak – – – – – – –   
Bam – – – – – – –   
Bandar Lengeh – – – – – – –   
Jask – – – – – – –   
Saravan 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.3351   
Station MA coefficients
 
θ1 θ2 θ3 θ4 θ5 θ6 θ7 θ8 θ9 
Shiraz 0.0000 0.0000 0.0000 –1.0000 – – – – – 
Sirjan 0.0000 0.0000 0.0000 –1.0740 0.0000 0.0000 0.0000 0.0000 0.0000 
Kerman 0.0000 0.0000 0.0000 –0.9647 0.0000 0.0000 0.0000 0.0000 0.0000 
Shahre Babak 0.0000 0.4378 0.0000 –0.9397 0.0000 –0.3769 – – – 
Bam 0.0000 0.2328 0.0000 –1.2899 0.0000 –0.2070 0.0000 0.3155 – 
Bandar Lengeh 0.0000 0.0000 0.0000 –1.2538 0.0000 0.0000 0.0000 0.2537 – 
Jask 0.0000 0.0000 0.0000 –1.2335 0.0000 0.0000 0.0000 0.3187 – 
Saravan 0.0000 0.0000 –0.1767 –1.0604 0.0000 0.0000 0.3871 0.0000 0.0000 
 θ10 θ11 θ12 θ13 θ14 θ15 θ16 θ17 θ18 
Shiraz – – – – – – – – – 
Sirjan 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 –0.1365 
Kerman 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 –0.0306 0.0000 –0.0310 
Shahre Babak – – – – – – – – – 
Bam – – – – – – – – – 
Bandar Lengeh – – – – – – – – – 
Jask – – – – – – – – – 
Saravan 0.0000 0.0000 0.0000 0.0000 –0.2883 0.0000 0.0000 0.0000 0.0000 
 θ19 θ20 θ21 θ22 θ23 θ24 θ25   
Shiraz – – – – – – –   
Sirjan – – – – – – –   
Kerman – – – – – – –   
Shahre Babak – – – – – – –   
Bam – – – – – – –   
Bandar Lengeh – – – – – – –   
Jask – – – – – – –   
Saravan 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.3351   

The zero values indicate the non-significant coefficients.

The dash lines indicate the absence of coefficient in model at the related station.

At Fasa, Bushehr, Bandar Abbas, Iranshahr and Chabahar stations, the ARMA model was diagnosed as the best model and the general form of this model is illustrated in the following equation: 
formula
(27)
The ARMA coefficients (Øi and θi) are presented in Table 7.

We checked the residual ACF/PACF of the models and the randomness of the residuals. According to the P-values obtained for the Ljung–Box statistic in various lags (Table 5), homogeneity and randomness of residuals (p > 0.05) and the reliability of forecasts at the level of 95% can be deduced.

Model validation

To validate the models with minimum AICC statistic that were obtained for different stations (Table 5), using the subset of observed data (1980–2007), the values for the 2008–2010 time period (approximately 10% of the time period used for modeling) were forecast and compared with values observed in this time period. The results of this comparison are presented in Table 8. According to the Pearson test, a significant correlation between observed and forecasted values is noticed at the significance level of 0.01. This comparison is also shown in Figure 5, including all of the synoptic meteorological stations of this study area. Also, Table 9 shows the results of NSE, RMSE and MAE measures to evaluate the goodness of fitted models for forecasting. High values of NSE and correlation coefficient (R-squared value) and small values for RMSE and MAE indicate a good model.

Table 8

The ARMA coefficients (θi) obtained for autoregressive moving average model

Stations ARMA coefficients
 
Ø1 Ø2 Ø3 Ø4 Ø5 Ø6 Ø7 
Fasa 0.0000 0.0000 0.0000 –1.0425 0.0000 0.0000 0.0000 
Bushehr –0.2459 – – – – – – 
Bandar Abbas 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 
Iranshahr 0.0000 0.1505 –0.1454 0.0000 0.2894 – – 
Chabahar –0.0864 – – – – – – 
 Ø8 Ø9 Ø10 Ø11 Ø12 Ø13 Ø14 
Fasa –0.9844 – – – – – – 
Bushehr – – – – – – – 
Bandar Abbas 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.1563 
Iranshahr – – – – – – – 
Chabahar – – – – – – – 
 Ø15 Ø16 θ1 θ2 θ3 θ4 θ5 
Fasa – – 0.0000 0.0000 0.0000 0.0000 0.0000 
Bushehr – – 0.1146 0.0000 0.0000 –0.8314 – 
Bandar Abbas 0.0000 –0.0209 0.0000 0.0000 0.0000 –1.0000 – 
Iranshahr – – 0.0000 0.0000 0.0000 –1.0002 – 
Chabahar – – 0.2072 0.0000 0.0000 –0.8485 – 
 θ6 θ7 θ8 θ9 θ10 θ11 θ12 
Fasa 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 –1.0001 
Bushehr – – – – – – – 
Bandar Abbas – – – – – – – 
Iranshahr – – – – – – – 
Chabahar – – – – – – – 
Stations ARMA coefficients
 
Ø1 Ø2 Ø3 Ø4 Ø5 Ø6 Ø7 
Fasa 0.0000 0.0000 0.0000 –1.0425 0.0000 0.0000 0.0000 
Bushehr –0.2459 – – – – – – 
Bandar Abbas 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 
Iranshahr 0.0000 0.1505 –0.1454 0.0000 0.2894 – – 
Chabahar –0.0864 – – – – – – 
 Ø8 Ø9 Ø10 Ø11 Ø12 Ø13 Ø14 
Fasa –0.9844 – – – – – – 
Bushehr – – – – – – – 
Bandar Abbas 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.1563 
Iranshahr – – – – – – – 
Chabahar – – – – – – – 
 Ø15 Ø16 θ1 θ2 θ3 θ4 θ5 
Fasa – – 0.0000 0.0000 0.0000 0.0000 0.0000 
Bushehr – – 0.1146 0.0000 0.0000 –0.8314 – 
Bandar Abbas 0.0000 –0.0209 0.0000 0.0000 0.0000 –1.0000 – 
Iranshahr – – 0.0000 0.0000 0.0000 –1.0002 – 
Chabahar – – 0.2072 0.0000 0.0000 –0.8485 – 
 θ6 θ7 θ8 θ9 θ10 θ11 θ12 
Fasa 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 –1.0001 
Bushehr – – – – – – – 
Bandar Abbas – – – – – – – 
Iranshahr – – – – – – – 
Chabahar – – – – – – – 

The zero values indicate the non-significant coefficients.

The dash lines indicate the absence of coefficient in model at the related station.

Table 9

Forecast performance of the best models at different stations

Station NSE RMSE MAE R2 
Abadeh 62.0 0.590 0.443 0.770 
Shiraz 48.6 0.692 0.467 0.619 
Fasa 36.8 0.691 0.419 0.510 
Sirjan 41.6 0.750 0.526 0.492 
Kerman 78.3 0.437 0.324 0.785 
Shahre Babak 50.0 0.722 0.536 0.504 
Bam 49.4 0.676 0.493 0.500 
Bushehr 53.3 0.588 0.375 0.686 
Bandar Abbas 72.9 0.466 0.400 0.735 
Bandar Lengeh 57.6 0.604 0.487 0.593 
Iranshahr 61.5 0.636 0.584 0.738 
Jask 78.2 0.434 0.357 0.793 
Chabahar 33.7 0.802 0.655 0.410 
Saravan 31.1 1.044 0.874 0.405 
Zabol 31.0 0.598 0.468 0.520 
Zahedan 69.7 0.565 0.433 0.745 
Station NSE RMSE MAE R2 
Abadeh 62.0 0.590 0.443 0.770 
Shiraz 48.6 0.692 0.467 0.619 
Fasa 36.8 0.691 0.419 0.510 
Sirjan 41.6 0.750 0.526 0.492 
Kerman 78.3 0.437 0.324 0.785 
Shahre Babak 50.0 0.722 0.536 0.504 
Bam 49.4 0.676 0.493 0.500 
Bushehr 53.3 0.588 0.375 0.686 
Bandar Abbas 72.9 0.466 0.400 0.735 
Bandar Lengeh 57.6 0.604 0.487 0.593 
Iranshahr 61.5 0.636 0.584 0.738 
Jask 78.2 0.434 0.357 0.793 
Chabahar 33.7 0.802 0.655 0.410 
Saravan 31.1 1.044 0.874 0.405 
Zabol 31.0 0.598 0.468 0.520 
Zahedan 69.7 0.565 0.433 0.745 

CONCLUSIONS

Drought predicting plays an important role in the planning and management of water resource systems by significantly reducing drought impacts. In this study, different statistical models were evaluated to select the best time series models for forecasting droughts based on the RDI index in southern Iran. For comparison and evaluation of forecasting models, various performance measures (NSE, RMSE, MAE and R2) were used. Drought classes of observed data (1980–2010) varied between ‘very wet’ to ‘severely dry’. The occurrence frequency of the ‘severely dry’ class was considerably higher compared to other drought classes. This result indicates the drought occurrence in this time period in this study area. Drought class probabilities indicated the lowest severity of drought in Bam, Bandar Abbas, Bandar Lengeh, Jask, Chabahar and Zabol stations and the highest severity of drought in Abadeh, Shiraz, Fasa, Sirjan, Kerman, Shahre Babak and Saravan stations.

Between different models (AR, MA, ARMA and ARIMA) fitted to the computed three-month RDI time series for observed data (1980–2010) in all of the stations, the model with the minimum value of AICC was selected as the best model. Results indicated that in 18.8% of cases for the AR models, in 50.0% of cases for the MA models and in 31.3% of cases, the ARMA models had the best performance in seasonal drought forecasting. For evaluating the forecast performance of all fitted models, different measures of goodness of fit (NSE, RMSE and MAE) were computed. The results verified the goodness of fitted models.

The results of this study are applicable for water resources managers, decision makers and governments to face this phenomenon well prepared. They can achieve this by proper management of water resources consumption, especially in the agricultural sector, which is the main consumer of water resources in Iran, in order to control or prevent the effects of drought.

REFERENCES

REFERENCES
Abramowitz
M.
&
Stegun
I. A.
1965
Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables
.
Applied Mathematics Series – 55
,
Washington, D.C
.
Akaike
H.
1974
A new look at the statistical model identification
.
IEEE Trans. Autom. Contr.
19
,
203
217
.
Allen
R. G.
,
Pereira
L. S.
,
Raes
D.
&
Smith
M.
1998
Crop Evapotranspiration
.
FAO Irrigation and Drainage Paper 56
.
Food and Agriculture Organization
,
Rome
.
Amei
A.
,
Fu
W.
&
Ho
C. H.
2012
Time series analysis for predicting the occurrences of large scale earthquakes
.
Int. J. Appl. Sci. Technol.
7
(
2
),
64
75
.
Asadi Zarch
M.
,
Malekinezhad
H.
,
Mobin
M.
,
Dastorani
M.
&
Kousari
M.
2011
Drought monitoring by reconnaissance drought index (RDI) in Iran
.
Water Resour. Manage.
25
,
3485
3504
.
Azarakhshi
M.
,
Mahdavi
M.
,
Arzani
H.
&
Ahmadi
H.
2011
Assessment of the Palmer drought severity index in arid and semi-arid rangeland: case study: Qom province
.
Iran. J. Range Watershed
16
,
77
85
.
Bahrami
M.
,
Bazrkar
S.
&
Zarei
A. R.
2018
Modeling, prediction and trend assessment of drought in Iran using standardized precipitation index
.
J. Water Clim. Change
.
https://doi.org/10.2166/wcc.2018.174 (in press).
Barua
S.
,
Perera
B. J. C.
,
Ng
A. W. M.
&
Tran
D.
2010
Drought forecasting using an aggregated drought index and artificial neural network
.
J. Water Clim. Change
1
(
3
),
193
206
.
Box
G. E. P.
&
Jenkins
G. M.
1976
Time Series Analysis, Forecasting and Control
.
Holden-Day
,
San Francisco
.
Brockwell
P. J.
&
Davis
R. A.
2002
Introduction to Time Series and Forecasting
,
2nd edn
.
Springer-Verlag
,
New York
.
Dastorani
M. T.
&
Afkhami
H.
2011
Application of artificial neural networks on drought prediction in Yazd (Central Iran)
.
Desert
16
,
39
48
.
Dinpashoh
Y.
,
Jhajharia
D.
,
Fakheri-Fard
A.
,
Singh
V. P.
&
Kahya
E.
2011
Trends in reference crop evapotranspiration over Iran
.
J. Hydrol.
399
,
422
433
.
Hong-yan
L.
,
Li-jun
X.
&
Xiao-jun
W.
2015
Relationship between solar activity and flood/drought disasters of the Second Songhua river basin
.
J. Water Clim. Change
6
(
3
),
578
585
.
Jahandideh
M.
&
Shirvani
A.
2011
Drought forecasting based on the standardized precipitation index using time series models in Fars Province
.
Iran. Water Res. J.
5
(
9
),
19
28
.
Ljung
G. M.
&
Box
G. E. P.
1978
On a measure of lack of fit in time series models
.
Biometrika
65
,
297
303
.
Nohegar
N.
,
Heydarzadeh
M.
&
Malekian
A.
2013
Assessment of severity of droughts using geostatistics method (case study: southern Iran)
.
Desert
18
,
79
87
.
Paul
A.
,
Bhowmik
R.
,
Chowdary
V. M.
,
Dutta
D.
,
Sreedhar
U.
&
Sankar
H. R.
2017
Trend analysis of time series rainfall data using robust statistics
.
J. Water Clim. Change
8
(
4
),
691
700
.
Phuong
D. N. D.
,
Linh
V. T.
,
Nhat
T. T.
,
Dung
H. M.
&
Loi
N. K.
2018
Spatiotemporal variability of annual and seasonal rainfall time series in Ho Chi Minh city, Vietnam
.
J. Water Clim. Change
.
https://doi.org/10.2166/wcc.2018.115 (in press).
Schwarz
G.
1978
Estimating the dimensions of a model
.
Ann. Stat.
6
(
2
),
461
464
.
Shamshirband
S.
,
Gocic
M.
,
Petkovi
D.
,
Saboohi
H.
,
Herawan
T.
,
Mat Kiah
L.
&
Akib
S.
2015
Soft-computing methodologies for precipitation estimation: a case study
.
IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens.
8
(
3
),
1353
1358
.
Soltani
S.
,
Modarres
R.
&
Eslamian
S.
2007
The use of time series modeling for the determination of rainfall climates of Iran
.
Int. J. Climatol.
27
,
819
829
.
Tsakiris
G.
2004
Meteorological Drought Assessment
.
Paper prepared for the needs of the European Research Program MEDROPLAN (Mediterranean Drought Preparedness and Mitigation Planning)
,
Zaragoza
,
Spain
.
Tsakiris
G.
&
Vangelis
H.
2005
Establishing a drought index incorporating evapotranspiration
.
Eur. Water
9
(
10
),
3
11
.
Tsakiris
G.
,
Pangalou
D.
&
Vangelis
H.
2007
Regional drought assessment based on the reconnaissance drought index (RDI)
.
Water Resour. Manage.
21
(
5
),
821
833
.
Tsakiris
G.
,
Nalbantis
I.
,
Pangalou
D.
,
Tigkas
D.
&
Vangelis
H.
2008
Drought meteorological monitoring network design for the Reconnaissance Drought Index (RDI)
. In:
1st International Conference ‘Drought Management: Scientific and Technological Innovations’
,
Zaragoza, Spain
.
UNEP
.
1992
World Atlas of Desertification
.
Edward Arnold
,
London
.
Wang
Y.
,
Li
J.
,
Feng
P.
&
Hu
R.
2016
Analysis of drought characteristics over Luanhe River basin using the joint deficit index
.
J. Water Clim. Change
7
(
2
),
340
352
.
Zarei
A. R.
&
Moghimi
M. M.
2017
Environmental assessment of semi-humid and humid regions based on modelling and forecasting of changes in monthly temperature
.
Int. J. Environ. Sci. Technol.
DOI 10.1007/s13762-017-1600-z
.
(in press)
.
Zarei
R.
,
Sarajian
M.
&
Bazgeer
S.
2013
Monitoring meteorological drought in Iran using remote sensing and drought indices
.
Desert
18
,
89
97
.
Zarei
A. R.
,
Moghimi
M. M.
&
Mahmoudi
M. R.
2016
Analysis of changes in spatial pattern of drought using RDI index in south of Iran
.
Water Resour. Manage.
30
(
11
),
3723
3743
.