Evaluation and Forecasting Meteorological Drought, Case Study: Kohgilooyeh and Boyer Ahmad

The SPI is the most widely used drought index to provide an acceptable estimation of drought characteristics. The objective of this study was to compare different threshold levels effect 27 on derived drought characteristics, assessment of the spatial variation of meteorological drought 28 properties as well as drought frequency, duration, and value in Kohgilooyeh and Boyer Ahmad 29 Province, Iran, using SPI for 1, 3, 6, 12, 24 and 48 months lead-times, and finally SPI forecasting 30 using Artificial Neural Networks (ANNs). For the first threshold level (scenario), drought 31 properties are extracted based on the standard level of zero, and for the second one, -1 is 32 considered. Results showed that the frequency of drought and wet periods decreased from SPI-1 33 to 48 for both scenarios in all stations. Max drought duration of stations had an increasing trend 34 from SPI-1 to 48. The average duration of dry periods changed as a function of the time scales; it 35 increased from SPI-1 to 48. Spatial variation of the drought average duration was considerable for 36 long-term drought. Max SPI value did not follow any spatial variation, as it was constant for all lead times in all stations. Average SPI values had a decreasing trend from SPI-1 to 9 but increased 38 from SPI-9 to 48 in all stations. Max average of SPI value observed in short-term drought and min 39 value in medium-term. SPI value general trend was similar in both scenarios, therefore drought 40 threshold level did not affect the results. The third objective was to develop neural network models 41 for drought forecasting. Different architectures are applied to find the best models to forecast SPI 42 over various lead times. The best forecasting results for SPI-3 and 6, obtained from the Quasi- 43 Newton training algorithm, when for SPI-1, 9, 12, 24, and 48, Levenberg-Marquardt was the best. 44 There was an increasing trend in performance measure R2 from SPI-1 to 48 and a decreasing trend 45 in Root Mean Square Error (RMSE). The best input lead-time for SPI-1 to 48 decreased from 11 46 to 1, the number of hidden layers decreased, but there was no significant trend in hidden neurons. 47 Drought properties could be considered in water resources management to supply water for various 48 demands. 49


Introduction 53
Drought has become more important in water resources management. It is recognized as an 54 environmental disaster, occurs in all climate zones. Drought is related to the precipitation, 55 temperature, relative humidity, rainfall intensity, duration, and temporal distribution (Mishra and 56 Sing, 2010). Droughts impact surface and ground waters and can lead to reduced water supply, 57 power generation and recreation activities, low water quality, crop failure, and economic and social 58 activities (Riebsame et al., 1990). Due to population growth, agricultural, energy, and industrial 59 5 2000), but fail to give helpful insights into the system, like white-box models (Vos and Rientjes, 129 2005). For example, ERs and Auto-Regressive Moving Average (ARMA) time series have a 130 limited ability to capture non-stationarities and non-linearities in hydrologic data, so alternative 131 models developed. Among various DDMs, ANNs are the most popular (Solomatine et al., 2008). 132 In recent decades ANNs have shown considerable ability in hydrologic and water resources. An 133 application of ANNs to solve civil engineering problems began in the 1980 decade (Flood and 134 Kartam, 1994 a,b). Primary concepts and capabilities of ANNs in hydrologic modeling are 135 described in ASCE (2000 a, b) and Govindaraju and Rao (2000). Because of the rehabilitation of 136 complex non-linear dependencies, ANNs are more accurate than models like ARMA and linear 137 regression. 138 This paper focuses on meteorological drought using SPI for 1, 3, 6, 12, 24, and 48 months 139 lead-times (moving average) to evaluate the effect of different threshold levels on drought 140 characteristics and assessment of the spatial variation of derived drought properties such as 141 frequency, duration, and value. Different return periods are considered because of two reasons. 142 First: SPI-1 and 3 represent short-term, SPI-6 and 9 represent medium-term, and SPI-12, 24, and 143 48 represent long-term drought. Second: SPI series are strongly dependent on adjacent SPI series. 144 The ANN models were used to forecast SPI in different lead times. With our knowledge, 145 this kind of research has not been performed in the Kohgilooye and Boyer Ahmad province yet. 146 The novelty of the present study is to calculate the time series of SPI in multiple time scales, 147 comparison of different threshold levels effect on drought characteristics, and assessment of ANN 148 models to forecast meteorological drought over different lead-times in this region. 149 150 2. Methodology 151  2.1. Standardized Precipitation Index (SPI) 152 SPI used to study different aspects of droughts such as forecasting (Mishra and Desai, 2006; 153 The SPI is computed based on the long-term precipitation record. Records fitted to a 160 probability distribution function (PDF), then transformed to a normal PDF (Mishra and Desai,161 2006), so the mean SPI for the desired period is zero (Mckee et al., 1993). The length of record 162 and PDF nature could be the SPI limitations (Mishra and Sing, 2010). It is performed separately 163 for each temporal basis. On the SPI standard classification, a drought event occurs when the 164 volume of SPI is continuously negative and ends when it becomes positive (Table 1) The gamma distribution PDF is defined as: 169 Where > 0 is a shape factor, > 0 is a scale factor, and > 0 is the amount of precipitation. 171 ( ) is the gamma function which is defined as: 172 Edwards and Mckee (1997) suggested Thom (1958) equations for maximum likelihood α 174 and β parameters approximation as follows: 175 The resulting parameters used to find the cumulative probability of observed precipitation 180 for the given month and time scale: 181 Substituting t for x β reduces the equation to the incomplete gamma function. Since the gamma 183 function is undefined for x = 0 and a precipitation distribution may contain zeros, the cumulative 184 probability becomes: 185 Where q is the probability of zero precipitation. 187 The cumulative probability H(x) is transformed to the standard normal random variable Z 188 with mean zero and variance one (value of SPI

203
Through calibration (training), the internal pattern of neuron connectivity (weight) is 204 determined, based on the data given to the network (Vos and Rientjes, 2005). The weights (Wi) 205 reflect the importance of the correlation between inputs and outputs. The initial selection of the 206 weights is performed randomly when weight adjustment is completed during training (learning). 207 The bias is considered as a weight of an input variable, increasing the network flexibility (Farmaki 208 et al., 2010). ANNs often trained using algorithms that minimize a performance measure such as 209 the RMSE. The relations between input and output are determined by training. There are two types of 225 training ANNs: supervised and unsupervised. In supervised, a set of training data, containing 226 inputs and correlated outputs, used. During the iteration, the network adjusts the weights, in a way 227 that calculated and desired outputs are as close as possible. If the network is properly trained, it is 228 learned and can predict the outputs with an unknown function as a black-box model. A network 229 with more weights models a complex function but proved to be over-fitting, when a network with 230 few weights may not be powerful. The solution is to use another data set for model validation, to 231 check the trained model. Training and validation errors calculated. If the first decreased and the 232 second increased, the network starts to over-fit, but if both decreased, the network was trained 233 correctly Figure 4. In Unsupervised training, the network is provided with a data set, trying to learn 9 Gradient Descent such as BP, Newtonian optimization, and Levenberg-Marquardt (L-M) are the 236 most popular training algorithms, documented by Haykin (1999). Alternative algorithms are linear 237 least squares and simplex optimization (Hsu et al., 1995) and global optimization methods such as 238 Simulated Annealing and Genetic Algorithm (Goldberg, 2000). The performance of the predictions, resulting from many models such as ANNs, evaluated 258 by the following goodness of fit measures: 259 Where the subscript m and s represent the observed and simulated SPI values, respectively; 261 p = total number of events considered, and 262 Where ̅ is the mean value of the x. 264 265

Case Study 266
The study was carried out in the Kohgilooye and Boyer Ahmad province, Iran, including 267 major parts of the three important river basins of Karoon, Maroon-Jarahi, and Zohreh-Hendijan, 268 located in the Sought West of Iran, Figure 5. The watershed boundaries are between 30°, 9 ′ to 269 31°, 32 ′ N and 49°, 57 ′ to 50°, 42 ′ E with a 16249 km 2 area. It is a mountainous area, with a wide 270 range of altitudes, from 4409 m in Dena to less than 500 m in Lishtar, with a complex topography 271 and dominated by a steep slope. In the low elevated areas the mean annual precipitation is 350 mm 272 and in the elevated areas more than 800 mm. Temperature variation range in low lands is between 273 10 ℃ to 47 ℃ and in elevated areas between −10 ℃ to 37 ℃ . 274

Results 292
The precipitation series was examined (Peterson et al., 1998) using the Mann-Whitney test 293 to assess the homogeneity of the data. The temporal gaps (<10%) in the meteorological stations; 294 data completed using grid fit upon the reference series.

Drought characteristics in different threshold levels 312
For extracting SPI characteristics, two threshold levels (scenarios) were considered. In the 313 first scenario, drought properties were extracted based on the standard level of zero (Table 1). 314 Thus, when SPI values were less than zero (negative), near normal to extreme droughts were 315 determined. In the second scenario, the drought properties were extracted based on -1 as the level 316 and compared with the first level. 317 The derived SPI characteristics of all studied stations are presented in Figure 7(a-e); the 318 number of drought events in Figure (a); max drought duration in (b); average drought duration in 319 (c); max SPI value in (d) and average SPI value in (e). It could be seen that the frequency of 320 drought events decreased from SPI-1 to 48 for both scenarios in all of the stations. It was also 321 observed that the spatial variation of SPI-1 to 12 (short and medium) was more than the others, 322 Figure 7(a). 323 In the second threshold level, the frequency of events of all lead times increased, which was 324 due to the splitting of events that merged with a mild drought (between 0 and -1). It depends on 325 the severity of drought events, as in the stations with mild drought, the number of events decreased. 326 It was also dependent on the lead times because the number of droughts with a low lead-time 327 increased but number of high lead-time events decreased. SPI with low lead-time often have a 328 lower severity than high lead-time drought, so they are divided into some events. 329 Max drought duration had an increasing trend from SPI-1 to 48 in both scenarios, figure  330 7(b). The difference of max drought duration varied in different SPI lead-times, as for SPI-1, 3 331 (short), and 48 (long) there was not any significant difference among stations, but for SPI-6 and 9 332 (medium), the variation was not neglectable. In the second level, there was seen a decreasing trend 333 from SPI-1 to 6 and an increasing trend from SPI-6 to 48 (with a different rate from the first level). 334 The average SPI value has a decreasing trend from SPI-1 to 9 but increased from SPI-9 to 347 48 in both scenarios. This trend could be seen in all stations, Figure 7(e). The min average SPI 348 value was observed for SPI-1 and 3, and max in SPI-6 and 9. The difference of this drought 349 character was considerable for all SPI lead times, in all of the stations.

ANN 364
NNs are a class of models discover pattern from the data. Although many types of NNs 365 developed, MLFF is the widest architecture for the prediction and forecasting of hydrologic 366 variables (Gupta et al., 2000). In an FF-BP, the weights connections feed activities only in the 367 forward direction, input to the output layer. In this paper, an MLFF-BP for forecasting time steps 368 is discussed. Hidden nodes with appreciate non-linear transfer functions are used to process the 369 input information. The nodes of the hidden layer allow to capture of the data pattern and perform 370 the non-linear mapping between inputs and outputs. The hidden nodes also allow modeling the 371 trend and seasonal variations as non-stationarities (Maier and Dandy, 1996). Increasing the number 372 of parameters by adding hidden neurons and layers, complicates network training. To find the 373 optimal size of hidden layer neurons, the number of layers systematically increases, till the network 374 performance no longer improved (ASCE, 2000 a). In this study, the number of hidden nodes 375 increased from 2 to 10, and the performance of the model was tested. 376 The hidden layers play crucial roles in most successful training. To obtain the optimal 377 network architecture, the number of layers is determined by iterations. The layers increased from 378 training algorithm (derived from BP) tested here to find the best model. The activation function 381 determines the relationship between the input and outputs. The FF-NNs adopted in this study used 382 the tan-sigmoid activation function as the most popular choice (Elshorbaghy et al., 2009). 383 The ANN is used for time series forecasting when the input nodes are reconnected to some 384 past observed values to predict the variable at future time-steps, so the number of input nodes 385 corresponds to the number of lagged observations (Mishra and Desai, 2006). In this study, the 386 number of lagged observations increased from 1 to 14 to select the best model. 387 Data sets normalized before training, using: 388 Where and 0 represent the normalized and original data, and and the 390 minimum and maximum values of original data. The network trained for 1000 epochs, using the 391 BP algorithm with a learning rate of 0.01, and a momentum coefficient of 0.9. The available data 392 was split into three parts, 70% for calibration, 20% validation, and 10% for testing. 393 The optimal architecture for SPI-1 was the L-M training algorithm with 3 hidden layers and 394 10 hidden neurons. Comparison of model performance measures presented in Table 3, for 1 to 14 395 months lead-times. The model with the best performance for SPI-1 resulted from 11 months lead-396 time (input neurons), with total R2=0.7328 and RMSE=0.188. Figure 8 shows a scatter plot of 397 observed versus forecasted SPI-1 for training, validation, testing, and all periods. The scattering 398 around the 45-degree line supports the conclusion made earlier. Performance measures of 399 forecasted SPI-1, 3, 6, 9, 12, 24, and 48, presented in Table 4. The best results of SPI-3 and 6 were 400 obtained from the Quasi-Newton training algorithm when for SPI-1, 9, 12, 24, and 48, L-M was 401 the best one. 402 Comparing all models, there was an increasing trend in the measure R2 from SPI-1 to 48; 403 0.7328 to 0.9934 and decreasing trend in RMSE; 0.188 to 0.1194. It can be concluded that long-404 term SPIs are forecasted more precisely. The best input lead-time of SPI-1 to 48 decreased from 405 11 to 1, and the number of hidden layers also decreased. It could have resulted from the data 406 volume used for the calculation of SPIs, with different lead times, but no significant trend was 407 observed in the hidden neurons. 408 409

Conclusion 417
Since precipitation is the base of the water supply component, an analysis of precipitation 418 deficit characteristics is critical in drought risk assessment. SPI is based only on the precipitation 419 data. It is standardized and can be computed on different time scales, allowing to monitor the 420 various kinds of drought. In this study, the SPI drought indices for 1 to 48 months lead-times were 421 computed, and the characteristics were extracted in all meteorological stations, at two threshold 422 levels. In the first scenario, drought properties were extracted based on the standard level of zero, 423 and in the second one, derived with -1 as the threshold level. Then the properties were compared. 424 Evaluation of SPI characteristics highlighted the spatial variations of the meteorological 425 drought frequency, duration and, value in the province. The frequency of drought events decreased 426 from SPI-1 to 48 in all stations, at both thresholds. It could be from the moving-average effect in 427 different SPI lead times. Max drought duration of stations had an increasing trend from SPI-1 to 428 48. The average duration of dry periods changed noticeably as a function of the time scales, as it 429 increased from SPI-1 to 48. Spatial variation of the average drought duration was considerable in 430 a long-term drought. Max SPI value did not follow any spatial variation, as it was constant for all 431 lead times in all of the stations. Average SPI value had a decreasing trend from SPI-1 to 9 but 432 increased from SPI-9 to 48. A max average of SPI was observed in short drought and the minimum 433 value in medium SPI. The general trend was the same in both scenarios, so the threshold level did 434 not affect the results. The spatial variation of short, medium, and long-term drought could be used 435 in water resource management strategies to supply water for various demands. 436 One more objective of the study was to generate SPI time series, over multiple durations. It 437 was observed that at shorter time scales the dry and wet periods showed a high temporal frequency. 19 duration of dry periods changed noticeably as a function of time scale, as for SPI-3 the average 441 duration was 5.1 months and for SPI-9 was 7.2. Comparing all the time scales, there was a 442 decreasing trend of severity from SPI-1 to 48. 443 Another objective was to develop NN models for drought SPI forecasting. The application 444 of ANN in drought forecasting was successfully demonstrated. Different NN architectures were 445 applied to forecast SPI series and the best models selected, comparing observed and forecasted 446 SPIs over various lead times, using model performance measures. The hidden layers varied from 447 1 to 3, and corresponding hidden nodes increased from 2 to 10. The L-M, BP with Steepest 448 Descent, Conjugate Gradient Descent, and Quasi-Newton were tested as training algorithms. The 449 input lead times increased from 1 to 14 months, and the model with the best performance, selected 450 for SPI forecasting. 451 Evaluation of model performances showed that the Quasi-Newton training algorithm had the 452 best performance in forecasting SPI-3 and 6, when L-M was the best in forecasting SPI-1, 9, 12, 453 24, and 48. So the best training algorithm was a little different depending on the SPI time scale. 454 Comparing all models, there was an increasing trend in performance measure R2 for SPI-1 to SPI-455 48: 0.7328 to 0.9934 and decreasing one in RMSE from 0.188 to 0.1194. It could be concluded 456 that the ANN model was more accurate in forecasting long-term drought (correlated with 457 groundwater level and reservoir storage) than short-terms (related to the precipitation and soil 458 moisture). 459 The best input lead-time (neurons) of SPI-1 to 48 decreased from 11 to 1, and the number of 460 hidden layers also decreased, but there was no significant trend in hidden neurons. It could be from 461 the volume of information contained in long-term SPIs such as SPI-48 in comparison with the 462 medium and short-term SPIs. Thus the quality of input information (long time average) leads to 463 the simpler architecture of ANN. This paper highlights the importance of ANN models in 464 forecasting and comparison of drought properties over short, medium, and long lead times. The 465 results could be used in drought forecasting for water resources management plans.