Evaluation of drought events in various climatic conditions using data-driven models and a reliability-based probabilistic model

Due to a wide range of socio-economic losses caused by drought over the past decades, having a reliable insight of drought properties plays a key role in monitoring and forecasting the drought situations, and finally generating robust methodologies for adapting to the various vulnerability of drought situations. The most important factor in causing drought is rainfall, but increasing or decreasing the temperature and consequently, evapotranspiration can intensify or moderate the severity of drought events. Standardized Precipitation Evaporation Index (SPEI), as one of the most well-known indices in the definition of the drought situation, is applied based on potential precipitation, evapotranspiration, and the water balance. In this study, values of SPEI are formulated for various climates by three robust Artificial Intelligence (AI) models: Gene Expression Programming (GEP), Model Tree (MT), and Multivariate Adaptive Regression Spline (MARS). Meteorological variables including maximum temperature (Tmax), minimum temperature (Tmin), average temperature (Tmean), relative humidity (RH), 24-h rainfall (P24) and wind speed (U2) were used to perform the AI models. Dataset reported from four synoptic stations through Iran, dating back to a 58-year period beginning in 1957. Each AI technique was run for all the climatic situations: Temperate-Warm (T-W), Wet-Warm (W-W), Arid-Cold (A-C), and Arid-Warm (A-W). Results of AI models development indicated that M5 version of MT provided the most accurate SPEI prediction for all the climatic situations in comparison with GEP and MARS techniques. SPEI values for four climatic conditions were evaluated in the reliability-based probabilistic framework to take into account the influence of any uncertainty and randomness associated with meteorological variables. In this way, the Monte-Carlo scenario sampling approach has been used to assess the limit state function from the AI models-based-SPEI. Based on the reliability analysis for all the synoptic stations, as the probability of exceedance values declined to below 75%, drought situations varied from “Normal” to “Very Extreme Humidity”.


Introduction
Drought events, known as one of the most severe climatic phenomena, can significantly affect agricultural production, groundwater quality resources, and water resources management (Saadat et al. 2011;Maghrebi et al. 2020Maghrebi et al. , 2021Noori et al. 2021). Due to the difficulty of determining the drought duration, there is no reliable relationship for drought events which is globally accepted (Bartlett and Singh 2018). Different definitions of drought situations have been proposed by meteorologists with various perspectives based on different variables. One of the most comprehensive definitions provided by Palmer (1965). In accordance with this, drought is a period in which the amount of moisture or any other indicator has a negative anomaly compared to the average conditions of the region (Palmer 1965). Due to various influential factors of drought occurrence, drought types are generally divided into various groups: meteorological drought, agricultural drought, hydrological drought, and socio-economic drought (Heim Jr 2002). Numerous indicators have been presented to assess the types of drought, the most common of which can be Deciles Index (DI) (Gibbs and Maher 1967), Percent of Normal Precipitation Index (PNPI) (Willeke et al. 1994), Palmer Drought Severity Index (PDSI) (Palmer 1965), Surface Water Storage Index (SWSI) (Shafer and Dezman 1982), Standardized Precipitation Index (SPI) (McKee et al. 1993) and Standardized Precipitation Evaporation Index (SPEI) (Vicente-Serrano et al. 2010;Bayissa et al. 2018). Among these indices, SPI and SPEI are known as the most widely used indices. The SPI has been widely applied in previous investigations due to its straightforward computation and its robust capability to find drought situations at various time scales (Li et al. 2020). One of the major restrictions of the SPI is the lack of water balance based on the rate of evapotranspiration. SPEI index is a type of SPI which considers the water balance for calculating the drought index based on rainfall and potential evapotranspiration (Mouatadid et al. 2018). SPEI can integrate multi-timescales aspects of the SPI with climatic data about evapotranspiration process which makes it high applicable for climatic investigations. Statistically, SPEI needs only climatological data without presumptions related to the properties of the underlying system.
Estimation of drought indicators has several applications in agricultural research, irrigation, hydrological engineering, and water resources management. Over the last decade, quite a few researches have been conducted to detect drought situations through the world due to the importance of the SPEI prediction.
In this regard, Potop and Mozny (2011) evaluated the new SPEI in the Czech Republic and proved the high potential of the SPEI to find various situations of the drought event. Tornros and Manzel (2014) analyzed drought conditions in Jordan and found that the sixmonth SPEI indicated the highest dependency on estimated soil moisture. In addition, Normalized Difference Vegetation Index (NDVI) was selected as the most efficient index for determining annual variations of drought. Chunping et al. (2015) examined the temporal and spatial changes of drought using the two indicators SPI and SPEI in the situation of climate changes in China. From their study, it was found that SPEI has a higher potential of visualizing drought in climate change conditions than SPI (Tan et al. 2015). Le et al. (2016) did research on the climatic signals to predict SPEI in Vietnam's Hangu Province using an Artificial Neural Network (ANN). They found that the developed model could benefit from the development of long-term policies for reservoirs and irrigation operations and alternative plant plans in the field of drought risk. Mumtaz et al. (2018) conducted a study with the aim of developing a robust drought modeling strategy to predict SPI at different time scales in Pakistan. In their research, Adaptive Neuro Fuzzy Inference System (ANFIS) model was used to estimate the SPI. Ultimately, the ANFIS model developed for 6 and 12-month time scales had more accuracy than three-month SPI. Given the superiority of the ANFIS model and its potential of providing reliable information about uncertainty, the proposed model can be used to predict drought risk for strategic water and agricultural resource management decision-making (Ali et al. 2018). Additionally, Khan et al. (2020) explored for the first time the potential of developing drought forecasting models in Pakistan using three advanced Machine Learning (ML) methods: Random Forest (RF) and Support Vector Machine (SVM). From their investigations, it was found that SPEI is directly related to the relative humidity of the Mediterranean Sea and the northern district of the Caspian Sea. While in the temperate season, SPEI has a direct relationship with the humid district southeast of the Bay of Bengal and the northern regions of the Mediterranean Sea and the Caspian Sea. Liu et al. (2020) introduced a new drought index, known as the Integrated Agricultural Drought Index (IDI), based on ANN and Remote Sensing (RS) data in northern China. They found that IDI based on an ML method could reduce the hypothesis used in many existing indicators that input and output data (or meteorological information) are linearly correlated. Additionally, Dikshit et al. (2021) applied Deep Learning (DL) methodology to predict SPEI for New South Wales and Australia using data from the period of 1901-2018. From their study, they found that the DL technique is capable of forecasting accurately at short-term and long-term lead times.
One of the methods used to predict the drought event situations is the Artificial Intelligence (AI) techniques which could detect efficiently general patterns/variations between drought indices and meteorological data (e.g., Ozger et al. 2011;Belayneh et al. 2014;Abbot and Marohasy 2014;Salcedo-Sanz et al. 2016). For instance, ANN (Djerbouai and Souag-Gamane 2016;Khan et al. 2018;Mulualem and Liou 2020), ANFIS (Kisi et al. 2019;Rahmati et al. 2020), Group Method of Data Handling (GMDH) (Aghelpour et al. 2021), and Support Vector Regression (SVR) Dikshit et al. 2020) are the cases of Artificial Intelligence (AI) techniques which were employed to estimate time series modeling.
From the literature review, it was found that some robust AI techniques such as Gene Expression Programming (GEP), Multivariate Adaptive Regression Spline (MARS), and Model Tree (MT) have not been yet employed to forecast drought index. These AI techniques have two major advantages. Firstly, these AI models are capable of presenting reliable mathematical expressions for recognizing behavioral patterns of complex systems. Secondly, these formulation-based AI approaches can automatically reduce the number of input variables that have a low level of contribution to the performance of AI models. In this way, these types of AI techniques can be practically used for forecasting drought index because there is a wide range of meteorological variables (e.g., temperature, dew point, evapotranspiration, precipitation) affecting this drought index.
On the contrary, the most frequently used AI techniques (e.g., SVM, ANN, and ANFIS) to predict drought index which acted as a black box, were not capable of obtaining mathematical expressions based on meteorological data. Additionally, ignorance of this issue might be likely to reduce the applicability of black box-based AI models in comparison with formulation-based-AI models.
The current research work includes two major contributions. First, in order to create a mathematical model that incorporates the meteorological parameters affecting the SPEI, it is necessary to have empirical relationships based on AI models which can provide an accurate estimation of SPEI. In the case of the second contribution, according to the random nature of the climatic data, SPEI is unavoidably sensitive to any variation in all meteorological data; then, the evaluation of reliability analysis for the SPI categories is efficiently applied by experts. More importantly, it is clear that experts in different fields of water science have not paid yet attention to the reliability analysis of SPEI during the last decade. In this way, it can be said that there is a ferocious demand for obtaining AI modelsbased equations to accurately estimate SPEI and further reliability analysis.
In the current research, the performance of artificial intelligence models of Multivariate Adaptive Regression Splines (MARS), Gene Expression Programming (GEP), and Model Tree (MT) is examined for four various climatic conditions of Iran. After that, AI models are developed and additionally, the efficiency levels in the prediction of SPEI for all the climatic situations is evaluated by various statistical benchmarks. Ultimately, the reliabilitybased probabilistic framework is introduced for the proposed AI technique with efficient performance when the randomness and uncertainty of each value of SPEI is considered.

Overview of case study
Iran is a high plateau located in the northern hemisphere in the warm region and is one of the arid and semi-arid regions of the world. As in many parts of the world in Iran, especially in very arid regions with uneven spatial and temporal distribution of rainfall throughout the year is limited. The total annual rainfall in Iran is about 413 billion cubic meters, of which 280 billion cubic meters are returned to the atmosphere through evapotranspiration. As a result, the total amount of renewable water per year is less than 2000 cubic meters for arid and very dry areas. However, sometimes the balance of water resources in these areas is often negative (Kiafar et al. 2017). According to the various climatic conditions of Iran in this study, we divided Iran into four temperate climates: Temperate-Warm (Ramsar Station), Wet-Warm (Ahvaz Station), Arid-Cold (Tabriz Station), and Arid-Warm (Kerman Station). As seen in the Supplementary Materials file, Fig. S1 shows the geographical location of the areas of the case study. Additionally, Table 1 indicates the characteristics of the synoptic stations in each case study district.
The present study was conducted in four meteorological stations in the different climates of Iran. In this regard, maximum temperature (T max ), minimum temperature (T min ), average temperature (T mean ), relative humidity (R H ), 24-h rainfall (P 24 ), and wind speed (U 2 ) have been obtained from the Meteorological Organization of Iran. The climatic information was prepared daily with a statistical period of 58 years (1/1/1957 to 12/31/2014). In order to calculate the SPEI drought index on a 12-month time scale, Excel2013 software was used to convert daily data into monthly and then annual data. Some statistical characteristics of the meteorological variables for each climate are shown in Table 2. In each synoptic station, 685 series of climatic variables were accumulated. In general, 75% of the data (514 observations) were allocated for the implementation of the training phase and 25% of the data were allocated for the implementation of the testing phase.

Standardized precipitation evaporation index
At the end of the previous decade, Vicente-Serrano et al. (2010) used the SPEI to calculate the annual drought. Generally, SPEI is a multi-quantitative index in which precipitation and temperature data are used. Although the calculation methodology of the SPEI is the same as the SPI, there is a subtle difference between two methods. The difference is that SPEI methodology uses the values of the difference between precipitation and evapotranspiration potential (PET). The results of Mavromatis (2007) investigations showed that simple and complex methods for calculating potential evapotranspiration are similar in terms of accuracy level. Therefore, Vicente-Serrano et al. (2010) used the model of Thornthwaite (1948) to estimate potential evapotranspiration. In the present study, a programming code was provided in MATLAB R2013a software to calculate the SPEI. Rainfall and temperature data were defined as input as monthly, as well as latitude of synoptic stations and year of onset and time scale for all the climates. For this reason, PET is firstly calculated as follows: in which, T is the average monthly temperature (°C), m is the coefficient of dependence on I, I is the heat index or the sum of 12 months i, K is the correction factor in terms of month and latitude, NDM is the number of days in a month and N denotes the maximum number of hours of radiation. After calculating the evapotranspiration potential (PET), the difference between precipitation (P) and PET for i month was obtained.
To calculate the drought index, three-parameter distribution is needed to cover the negative values in the D data. The results of selecting the most appropriate distribution function show that the logarithmic logistic function fits well with the time series of data at various time scales. Thus, the cumulative probability function of the D data series based on the logarithmic logistic function is calculated as follows: where is the scale parameter, is the shape parameter and is the main parameter for D values in the range γ > D < ∞.
Finally, the SPEI is calculated as standardized values of F(x): in which, W = √ −2 ln (P) which is valid for P ≤ 0.5 and P variable denotes the probability of increasing the set values d. The values of L 0 , L 1 and L 2 as well as d 1 , d 2 and d 3 are constant. The classification of this index indicates that the SPEI situation changes from very extreme humidity (SPEI < −2) to very severe drought (SPEI > 2), as seen in Table S1 (See Supplementary Materials).

Gene-Expression Programming (GEP)
Gene expression programming (GEP) is a generalization of the two algorithms: Genetic Algorithm (GA) and Genetic Programming (GP) (Koza and Koza 1992). GEP is a search model which develops computer programs with various configurations such as mathematical expressions, polynomial expressions, and logical expressions. Different shapes, similar to decomposition trees, are combined in GP. The performance of the GA algorithm is based on nature and based on the survival of the fittest. In this way, it first forms a population of organisms and by applying a series of operators on this population, creates an optimal set that has specific characteristics of the population (Ferreira 2006;Gholampour et al. 2017;Iqbal et al. 2020). GEP was firstly used by Ferreira (2001). In GEP, a relationship is established by chromosomes which have uncomplicated, certain lengths, and linear configurations. The GEP approach has five steps of development. The first step is to fix the fitness function of an individual program. Mean Squared Error (MSE) is applied for the evaluation of fitness function. The second stage is determining the terminals and the mathematical function to generate the chromosomes. The third stage of development is identifying the structure of chromosomes. And more, the well-matching algebraic operators are essentially fixed to establish a general equation among genes. Finally, a broad range of genetic operators are fixed (Ebtehaj et al. 2015).
In this study, in order to estimate the drought index using the GEP model, GenXpro-Tools5 software was used. In GEP modeling, maximum and minimum temperature, average temperature, 24-h relative humidity, and wind speed were used as input parameters of the SPEI as an output parameter of the model. In the present research, GEP model uses a wide range of fitness functions (such as Relative Hits with Root Mean Square Error, Mean Square Error, Relative Absolute Error, Mean Absolute Error, positive correlation). Through the previous investigations, MSE criterion was successfully considered for evaluation of chromosome in each generation (Maroufpoor et al. 2019;Shiri et al. 2020;Kazemi et al. 2021). In this way, MSE criterion was used to assess the goodness of training and testing performance. The settings used in GEP modeling to estimate the SPEI are given in Table S2.
The best relationships extracted from 10,000 generations and 3 genes for estimating SPEI in the Ramsar, Ahvaz, Tabriz, and Kerman stations are, respectively, as follows: To obtain robust GEP models, there are three setting parameters that require to be regulated. Firstly, stop conditions of GEP models are fundamentally on some criteria: number of generation, random R-square threshold, and random generation number. However, we applied an empirically applicable alternative (R-square). In this investigation, this alternative is a good indication of GEP performance since GEP performance will continue till the goodness of performance corresponded to the training and testing phases is met. Another possibility is to select an addition linking function ( +) for the combination of subtrees in the GEP which was found to be a more powerful alternative. The final possibility of GEP settings is the selection of arc tangential, exponential, and tangent hyperbolic functions (introduced as inner functions) for the SPEI simulation in all synoptic stations.

Multivariate adaptive regression spline (MARS)
The MARS model, as a non-parametric regression method, was firstly proposed by Friedman (1991). The MARS model is capable of providing flexible, precise, and fast regression models for predicting continuous and numerical output variables (Yilmaz et al. 2018). Additionally, MARS has the capability to process a large volume of high-dimensional data. The main merit of MARS is the capability to express the complex and nonlinear (8) relationship between the predictor and response variables. The MARS model can model the equations between the predicted and observed independent variables. In the modeling process, the input variables (or vectors) are grouped into intervals (subsets) and the base functions which are proportional to each interval. The base function denotes information about input variables. A minus function is determined in a limited range, and its start and endpoints are named nodes. The node represents the point at which performance behavior changes (Mehdizadeh et al. 2017;Kisi and Parmar 2016;Fernández et al. 2013;Balshi et al. 2009). The basic functions of the relationship among input variables are shown as follows: where K is a node constant called (threshold limit), x is the input variable and y is the target variable. Also, the general relationship is described by the following formula: where i , C i , C o ′ and e are the base functions, fixed functions, initial constants, and the number of functions, respectively. The MARS model consists of two development stages. In the first stage (step forward), MARS predicts the dependent variable (or target values) with a constant value. This is the average value of the data of the target variable. This step will continue until all possible basic functions are added to the model. In step one, a complex and overly complex model are generated with a fair number of nodes. In the next step (step back), less significant basic functions are detected by using MARS model; therefore, these basic functions are eliminated by the pruning process. This process continues to check all the basic functions. Ultimately, an optimal model is obtained (Mehdizadeh et al. 2017).
In this research, the MARS model has been implemented using the programmed codes of MATLABR2013a software, so that the number of basic functions (NBF) for the Ramsar, Ahvaz, Tabriz, and Kerman stations are equal to 19, 37, 32 and 33, respectively; additionally, the initial constant values in the models are fixed −1.387, 0.16632, 0.81793 and −1.4038, respectively. Base functions and corresponding constant coefficients used in modeling drought estimates for the different climates are presented in Tables S3-S6. All the climatic variables had contributions to the creation of basis functions in order to predict SPEI values in all the climatic conditions (i.e., Temperate-Warm, Arid-Cold, Wet-Warm, and Arid-Warm).

Model Tree (MT)
Machine Learning (ML) techniques, especially decision tree modeling, have shown successful applications in water resources and flood management (Bhattacharya and Solomatine 2003;Quinlan 1992;Granata et al. 2018;Granata 2019). Although M5MT-based modeling has no popularity as well as GA, ANN, and GP, it was found to be significantly robust in many scientific and engineering applications (Singh 2007). The decision tree is a powerful and common tool for categorizing and representing a set of rules that lead to a category or value and has a (12) 1 3 tree structure similar to a flowchart (Berson et al. 1999). The M5 model tree, as the most commonly used classification in the family of decision tree, was proposed by Quinlan (1992) for predicting continuous data.
M5MT is one of the data mining methods that is inspired by the nature of trees, i.e., it has a structure similar to a tree (leaf branch root) and similar processes (growth and pruning) (Quinlan 1992). M5MT unlike other decision tree models that represent discrete classes or classes as output, it creates a multivariate linear model for the data in each node of the tree model. In order to generate M5MT, a division criterion is used. The standard deviation is the class values that reach a node as a quantity of error and calculate the expected reduction in this error as the test result of each attribute in that node. The standard deviation (SDR) reduction is obtained from the following equation (Witten et al. 2005;Quinlan 1992).
where SDR is the standard deviation reduction, SD is the standard deviation, S is the input of each node and M i is the input of the i th node.
Weka3.9 software was used to implement the M5MT model. The tree structure of the M5 algorithm for estimating drought is shown in Table S7. As seen in Table S7, SPEI values were estimated in the Ramsar station with Temperate-Wet climate using three model tree rules and additionally P 24 variable is only splitting variable for the creation of tree structure. In the Wet-Warm climate, MT created 13 rules with two splitting variables: P 24 and T mean . In Tabriz station with Wet-Cold condition, linear equations given by MT have six rules with two splitting variables P 24 and T mean . In the Arid-Warm climate, SPEI values were predicted by 15 rules and four splitting variables (i.e., P 24 and T mean , T max , and R H ). Also, the relationships extracted from the M5 models are shown for all the synoptic stations in Table S8.

Results and discussions
Results of this study are presented in two sections: (i) statistical measures of AI models to predict SPEI values for various climatic situations and (ii) reliability analysis of AI results.

Definition of statistical measures
To evaluate the performance of under proposed approaches and identify the best ones for drought index estimation in each climate, several statistical measures were employed. In this way, Index of Agreement (IA), Root Mean Square Error (RMSE), Mean Absolute Percentage of Error (MAPE), Bias Index (BIAS), and Scatter Index (SI) have been applied: where N is the number of samples, SPEI i Obs is the measured SPEI, SPEI i Pre is the predicted SPEI, SPEI Obs is the mean of measured SPEI, and SPEI Pre is the mean of estimated SPEI.
IA is a dimensionless parameter that is used as a standard criterion for evaluating the error of models and has a value in the range of [0, 1]. This index is indicative of the ratio of the mean square error and the potential error. Additionally, RMSE is used to measure the prediction accuracy of the numerical model, which, taking the square root of the errors, indicates a non-negative value and includes a value in the range of [0, '+ ∞) MAPE calculates the average absolute percentage of error in a set of predicted values. The values of this index are expressed as a percentage and, like RMSE, which low values of MAPE indicate good performance of the model. The BIAS represents the difference between the means of the observed and predicted data. The closer the values are to 0, the better the model performance. If the BIAS > 0 over prediction for numerical model (the mean of SPEI Pre values are greater than the mean of SPEI Obs values) is met, otherwise, under prediction takes place ( SPEI Obs > SPEI Pre ). Ultimately, SI value indicates the percentage of RMS difference with respect to the mean of SPEI Obs . Above-mentioned statistical measures were applied frequently to evaluate problems related to the various water resources management (Noori et al. 2010(Noori et al. , 2011Park et al. 2016;Deo et al. 2017;Rezaie-Balf et al. 2018;Kisi et al. 2019;Aghelpour et al. 2021).
This section evaluates the performances of the GEP, MARS and M5MT for the SPEI estimation. Additionally, strengths and weaknesses associated with models performance in the prediction of SPEI in each climate are discussed. Table 3 indicates the quantitative results of AI techniques in the SPEI prediction for the training stages. According to the statistical measures, MARS model has the best performance in all the climatic conditions, followed by M5MT and GEP models. In the Temperate-Warm conditions, MARS technique had more satisfying performance in terms of IA (0.  Similarly, for the Kerman station, the statistical measures associated with the MARS model indicated comparatively better performance than other AI models. Figure 1 illustrates the performance of AI models for various climates in the training stage. As seen in Fig. 1a, it was found that all the values of SPEI predicted by the MARS model had the best performance whereas, in Fig. 1b, the GEP model had the highest level of scattering data (SI = −20.569) in the Temperate-Wet climate. In the Wet-Warm climate, Fig. 1e illustrated the graphical performance of the GEP model had significant over/under prediction. In Fig. 1h, GEP model under predicted SPEI between −3 and −2 for the Tabriz station, whereas over prediction was seen for SPEI = 2-3. In the Kerman station with Arid-Warm climate, Fig. 1k indicated that GEP model under/over predicted partially for the SPEI = −3-0 and additionally, for SPEI = 2-3, under prediction was observed. Scatter plot observed in Fig. 1j illustrated the satisfying performance of MARS model and Fig. 1 indicated the relative good performance of M5MT.

Performance of data-driven models
As indicated in Table 4, the results of AI techniques in the SPEI estimation were provided in the testing stages. Statistical measures given in Table 4    For the Arid-Warm climate in the Kerman station, the statistical results related to the M5MT provided relatively better efficiency in the SPEI estimation than GEP and MARS models. Figure 2 presented the performance of AI models for different climatic situations in the testing stage. As illustrated in Fig. 2a, MARS technique had partial over prediction for SPEI between −3 and −2 in the Temperate-Wet climate whereas, in Fig. 3b, the GEP model [Eq. (9)] indicated both over/under predictions (SI = 4.851) for the SPEI values between 1 and 2. Additionally, Fig. 2c depicted insignificant over/under predictions for the performance of M5MT for SPEI between −3 and −1.5. In Fig. 2d, MARS model indicated significant overestimation for the SPEI values varying from −3 to −1 for the Ahvaz station, whereas a large amount of error prediction was illustrated in Fig. 2e for the SPEI values between −2 and −1. In addition, Fig. 2f illustrated good performance by M5MT. In the Tabriz station with Arid-Cold climate, Fig. 2g-i indicated relatively satisfying efficiency in the prediction of SPEI values. In the case of the Arid-Warm climate, as illustrated in Fig. 2l, the graphical performance of the GEP model indicated remarkable over prediction for the SPEI values between −3 and −1, whereas, for the SPEI = 1-3, relatively significant underestimation was seen. 1 3

Reliability-based probabilistic assessment of drought event
For a reliability-based probabilistic assessment of drought events, it is essential to address the problem with a limit state function (LSF). In this procedure, the LSF is used to describe the exceedance of Standardized Precipitation Evaporation Index (SPEI) from presumed values. Moreover, the uncertainty from the random nature of the effective parameters (or meteorological information) in the problem is addressed by introducing a series of random variables. In this regard, the probabilistic distribution for random variables has to be developed based on the database. These parameters are used to assess the SPEI by considering the randomness of each climatic variable. Hence, the next section deals with developing the probabilistic distribution of input variables with random nature.

Probabilistic distribution of random parameters
The measured data for climate parameters have been seen to change randomly. For such statistical data, results were reported in the form of the mean, maximum, and minimum values. However, to describe the yearly variation of these parameters, it would be much appropriate to develop a probabilistic distribution for these parameters. In this form, the effect of variation in each parameter can be assessed in probabilistic methods. To this end, the most suitable probabilistic distribution has to be selected for each random parameter. In this paper, the Kolmogorov-Smirnov test (also known as the K-S test) (Ang and Tang 2007) was used for this purpose. This method is a well-known approach for developing probabilistic models for random parameters. In this procedure, the cumulative distribution of the presumed distribution is compared to the numerical cumulative frequency of the random parameter. Then the maximum error is measured as: where G n (r) and F x (r) present the numerical and theoretical cumulative distribution of the parameter r, correspondingly. κ shows the significance level and is presented by typical engineering probability textbooks (Ang and Tang 2007). Generally, the K-S test is conducted for a group of presumed theoretical distributions and the one with the lowest D n that also satisfies the D n < D n κ condition is selected as the well-suited theoretical distribution. Fig. 3 Schema of the probabilistic model for assessing the confidence level for SPE Many probabilistic distributions can be used for capturing the statistical properties of a random parameter. A group of six well-known probabilistic distributions was considered in this study to see which of them is most suitable for any random parameter. These probabilistic distributions are Poisson, Exponential, Normal, Lognormal, Uniform, and Rayleigh distributions. Fig. S2 depicts results for the theoretical and the empirical cumulative distribution for two parameters in this study (i.e., T max and U 2 ) for the city of Kerman. As shown, the "Normal" and the "Lognormal" distributions are the most suitable distributions for these parameters. Such a procedure was also conducted for each parameter in every studied city. Results for D n were summarized in Fig. S3.

Probabilistic model for SPEI
As mentioned before, a limit state function (LSF) is employed in reliability methods for defining the boundary between the failure and safety of the problem. When it comes to drought assessment, such an LSF expresses the exceedance of Standardized Precipitation Evaporation Index (SPEI) from a specific value. The LSF can be expressed in the form of a mathematical relation which is a function of random variables: where SPEI A indicates the "Acceptable" threshold for the Standardized Precipitation Evaporation Index and SPEI M is the "Measured" value from the realization of random input parameters: In this equation, h is the function or algorithm and ε is the model error in the form of a random variable. ε is also known as the "uncertainty" parameter and shows the difference between the observed SPEI and the predicted ones from M5MT (Tables S7 and S8) for all the climatic situations. This parameter is usually considered as a random parameter with a normal distribution (Mahsuli and Haukaas 2013). As shown in Fig. S4, the residual quantile of the error against the normal quantiles was distributed along the 45° line which indicates the normality of errors.
Depending on the result of LSF, the failure or safety status is specified:

Reliability assessment of drought events through SPEI
Since the effective parameters in drought events have a random nature, results are sensitive to the variation of the input parameters (climatic variables). This is more outstanding when random parameters have a large dispersion. In this case, probabilistic methods are more suitable approaches for assessing the event since they consider both randomness and uncertainty of parameters, and hence, more realistic results are obtained. Such results can be employed for decision-making. By employing the reliability-based probabilistic analysis, the probability of exceedance (PE) (or probability of failure (p f )) for SPEI is measured. To this end, the violation of the LSF is assessed for different realizations of input parameters. As shown in Fig. 3, the realization of seven input variables (including the error) is generated through their probabilistic distribution. These values are used to assess SPEI M . Then, the violation of the LSF is evaluated for estimating the failure probability: where X 1 , X 2 ,…, X n are the limit state function input variables.
In this study, the Monte-Carlo scenario sampling was used for estimating p f for each SPEI A . To this end, realization for LSF is obtained from the scenario sampling of input parameters. Then, the probability of exceedance is assessed by computing the ratio of the number of samples with failure LSF to the total number of sampling: Figure S5 depicts the exceedance probabilities for different values of SPEI given by Eq. (23) where the randomness and uncertainty in each climatic variable were considered. The exceedance probability (PE) measures the chance of experiencing an event and how probable is it for an event to pass a specific limit (i.e., SPEI A ). Figure S5 illustrated a larger value of PE for the smaller values of SPEI in all the synoptic stations. In Fig. S5, all the curves follow a descending trend which is indicative of less PE for larger values of SPEI. To characterize the situation of drought for the case study areas, specific ranges of SPEI are introduced by quite a few descriptive terms, as seen in Table 3.
In the Ahvaz station with Arid-Warm climate, it was found that there are 85-99% chances for a climatic situation to exceed the "Very Severe Drought", indicating SPEI varies from −5 to −2. As chance declined to 80%, the situation of "Severe Drought" was met when SPEI values changed from −1.5 to −1.99. SPEI values in the "Moderate Drought" situation, the descending slope associated with the curve rises when the chance for that SPEI values varies from 75 to 80%. As seen in Fig. S5a, for "Normal" situation, the chance values declined to 30%. In the three situations of humidity, the chance values (PE) changed from 33% in the "Medium Humidity" to below 15% in the "Very Extreme Humidity". For other synoptic stations, variations of PE values versus SPEI ranges were similar to those reported by Ahvaz station.

Conclusion
The current investigation aimed to design a probabilistic framework for drought situations by using three AI models (i.e., MARS, GEP, and M5MT) by considering the SPEI evaluation in various climatic situations: Temperate-Warm, Wet-Warm, Arid-Cold, and Arid-Warm. To develop the AI models, the six meteorological variables were considered. (25) ... � f X 1 , X 2 , ..., X n dX 1 dX 2 , ..., dX n (26) I X 1 , X 2 , ..., X n , I X 1 , X 2 , ..., X n = 1 if LSF X 1 , X 2 , ..., X n ≤ 0 0 if LSF X 1 , X 2 , ..., X n > 0 MARS models presented basis functions-based-equations for all the climates when interaction among basis functions was the second-order polynomial. For each synoptic station, GEP model-based-formulation, in which all the climatic variables had contributions, was obtained by summation of three genes with satisfying accuracy level. In addition, M5MT provided a set of multivariable linear regression equations, rooting from if-then rules, to evaluate the SPEI for all the climatic situations. According to the performance of AI models, MARS models developed in the training stages of all the climates had the best efficiency and followed by the M5MT and GEP models. In the case of testing results, M5MT had satisfying performance in the prediction of SPEI. Therefore, the results of M5MT, as the most efficient equations in all the climates, were used to propose probabilistic models in order to evaluate the drought situations of the climates in the four case study areas.
Through results of reliability analysis in the Ahvaz station, in the case of humidity situations, probability of exceedance (PE) values obtained below 33% indicated humid climate with various severities from "Medium Humidity" for 25-33% to "Very Extreme Humidity" for PE < 15%. For other climatic situations, the analysis of PE values had the same trend as it was found for the Wet-Warm climate in the Ahvaz station.