The parameter estimation of statistical distributions is important for regional frequency analysis (RFA). The accuracy of different parts of RFA such as estimating the regional quantiles of the selected statistical distribution, determining the heterogeneity measure, and choosing the best distribution based on the Monte Carlo simulation, may be influenced by using the different values of regional parameters. To fulfill this aim, in the present study, a new model is developed for regional drought frequency analysis. This model utilizes the L-moments approach and the adjusted charged system search as an advanced meta-heuristic algorithm, in which some modifications on the equations of the algorithm are performed to improve its standard variant. The verification of the regional parameters estimated by the new methodology yields accurate results compared to other models. Furthermore, this study illustrates the usefulness of the robust discordancy measure against the classic one. For this purpose, different values of the subset factors (α) are utilized in the robust discordancy measure, and finally, the best value of subset factor is found equal to 0.8, which can accurately recognize discordant sites within the region.
INTRODUCTION
Drought is a normal part of any climate (Keshavarz et al. 2013). Generally, droughts take place when a region receives consistently below average precipitation (Hossein et al. 2012). Compared with other natural hazards, the damaging effects of droughts are non-structural, and much larger in terms of the spatial extent (Reddy & Ganguli 2012), where 22% of economic damages and 33% of the number of affected persons are related to droughts (Belayneh et al. 2014). Therefore, utilizing regional frequency analysis (RFA) seems necessary to recognize the properties of this phenomenon. RFA can circumvent the limitation of statistical estimates, such as the absence of lengthy records or too short records (Zhang et al. 2012). RFA leads to more accurate quantile estimations than at-site analysis (Santillán et al. 2014).
To perform RFA, the considered region should be homogeneous, otherwise, it should be subdivided into some homogeneous regions (Seckin et al. 2013). Among many homogeneity tests in the hydrologic literature, the L-moments (LMOM)-based tests proposed by Hosking & Wallis (1993, 1997) are known as the most powerful ones (Viglione et al. 2007). Thus, numerous studies have utilized the LMOM approach to the RFA of hydrologic phenomena such as rainfall, stream flow, drought, and other fields of water engineering (Modarres 2010; Chérif & Bargaoui 2013; Dikbas et al. 2013; Rahman et al. 2013; Rajsekhar et al. 2013; Zakaria & Shabri 2013; Dodangeh et al. 2014; Feng et al. 2014; Goyal & Gupta 2014; Sarhadi & Heydarizadeh 2014).
Detection of the possible outlying data in a region is an important matter for RFA (Saf 2010). Hosking & Wallis (1993, 1997) recommended the discordancy measure based on the mean and covariance of the LMOM ratios using the Mahalanobis distance. This measure can be strongly affected by the discordant sites and is not robust for outlying data. To improve this problem, Neykov et al. (2007) proposed the robust discordancy measure. A few studies have then used this measure in hydrology. Saf (2010) analyzed the effect of discordancy detection measures on regional flood probability types and the accuracy of the estimation on regional analysis. Ilorme & Griffis (2013) utilized the robust discordancy measure to determine physically discordant sites and proposed a new methodology to identify the physical attributes that are the greatest demonstrator of extreme hydrological response. The mentioned studies showed the capability of robust discordancy measure against the classic one for outliers’ identification; however, the importance of subset size for all data sets was not investigated. In this study, the different subset sizes are considered to evaluate the accuracy of the robust discordancy measure.
The LMOM method is also utilized for parameter estimation of statistical distribution for RFA, although there are no simple expressions and explicit solutions for fulfilling this aim. Therefore, several approximations or simplifications are suggested for solving such nonlinear and complex problems (Hassanzadeh et al. 2011; Reddy & Singh 2014). Meta-heuristic algorithms (MHA) are one useful tool to solve these problems and various MHA have been extensively used to solve the optimization problems of hydrologic phenomena (e.g., Rai et al. 2009; Hassanzadeh et al. 2011; Reddy & Singh 2014; Shin et al. 2014).
Recently, an MNA known as charged system search (CSS) was introduced by Kaveh & Talatahari (2010), and applied to different engineering optimization problems. CSS utilizes the Coulomb and Gauss's laws from electrostatics, and the Newtonian laws of mechanics. There are few studies on applying CSS in the field of water resources. Kaveh et al. (2012) developed four different models of composite open channels and then optimized them by CSS. The results of this method were compared to those of ant colony optimization and genetic algorithm (GA). Talatahari et al. (2013) applied CSS and particle swarm optimization (PSO) for identifying the parameters of the linear Muskingum model. In order to evaluate these algorithms, a numerical example was utilized and the results were compared to those of other algorithms. Sheikholeslami et al. (2014) utilized CSS to optimal cost design of water distribution networks. They utilized some well-known benchmark instances as case studies. Comparison of the results of CSS with some other MHA demonstrated the high performance of this algorithm. In this study, an adjusted CSS (ACSS) algorithm is applied for estimating the regional parameters of Kappa distribution, and then the results are compared with those of the LMOM as the conventional approach.
Considering these points, the main contributions of this paper are:
- (1)
presenting the robust discordancy measure with different values of subset factor and comparing with the classic discordancy measure;
- (2)
applying the ACSS algorithm for the parameter estimation of statistical distributions; and
- (3)
proposing a new model to give regional parameters for statistical distribution.
METHODOLOGY
ACSS
The ACSS algorithm as a meta-heuristic algorithm contains a number of charged particles (CPs), where each one treats as a charged sphere and can insert an electric force to the others. Thus, the resultant electrical force affect a CP results in its acceleration. ACSS uses the governing laws of motion from the Newtonian mechanics to determine the position of CPs. Application of these laws provides a good balance between the exploration and the exploitation of the algorithm (Kaveh & Talatahari 2010). The pseudo-code for the ACSS algorithm is summarized as follows.
Step 2: Charged memory creation. A number of the best CPs and the values of their corresponding objective functions are saved in the charged memory.
Step 8: Updating of charged memory. The best new vectors are included in the CM and the worst ones are excluded from the charged memory.
Step 9: Terminating criterion control. Steps 3–9 are repeated until a terminating criterion is satisfied.
It is notable that the steps 4, 6, and 7 are only applied in the ACSS algorithm.
LMOM approach
Hosking (1986, 1990) defined LMOM as linear combinations of probability weighted moments, which can be interpreted as measures of the location, scale, and shape of probability distributions. Based on LMOM, Hosking & Wallis (1993) proposed three useful statistics in RFA, the discordancy, heterogeneity, and goodness-of-fit measures. Detailed explanations and formula can be found in Hosking & Wallis (1997).
Discordancy measure
is the regional unweighted average and S is the matrix of sums of squares and cross-products. The critical discordancy measure is equal to 3 and the site becomes discordant when Di > 3 (Hosking & Wallis 1993, 1997). Due to the outlier observations’ influence on the classical estimates of and S, the number of discordant sites for a region can be smaller than the real number. Therefore, it is necessary to utilize the robust discordancy measure in RFA (Saf 2010).There are many robust estimators in statistics, such as minimum covariance determinant (MCD), minimum volume ellipsoid, and M-estimators (Alameddine et al. 2010). The MCD estimator is a highly robust estimator of multivariate location and scatter (Rousseeuw & Van Driessen 1999; Hubert et al. 2005; Rousseeuw et al. 2006). Rousseeuw (1984, 1985) defined the MCD estimator, which can be computed efficiently with the FAST-MCD algorithm. Finding the h observations (out of N) whose classical covariance matrix has the lowest possible determinant, is the main objective of MCD. There is a relation between the subset size (h) and subset factor (α) when the value of α varies from 0.5 to 1. If α= 0.5, h equals , and when α= 1, h equals N (Pison et al. 2002; Neykov et al. 2007). For other values of α, the h can be calculated by linear interpolation.
and S) in the Mahalanobis distance formula for RFA:
In this study, the MATLAB library for robust analysis, LIBRA, introduced by Verboven & Hubert (2005), is used for the computation of the robust discordancy measure.
Heterogeneity measure
One of the main steps in the RFA contains identifying homogenous regions (Saf 2009). For this purpose, the heterogeneity measure was proposed to compute the degree of heterogeneity in a region. The steps of the heterogeneity measure are illustrated as below (Hosking & Wallis 1993).
Step 1: Calculate the regional LMOM ratios.
Step 2: Compute the V variables.
Step 3: Fit the four-parameter Kappa (KAP) distribution.
Step 4: Simulate a large number of synthetic regions.
- Step 5: Drive the H indices. Finally, heterogeneity measures are given as follows:where μvi and σvi are the mean and standard deviation of the simulated V variables.16
Hosking & Wallis (1993) suggested that a region is ‘acceptably homogeneous’ if Hi is less than 1, ‘possibly heterogeneous’ if Hi is between 1 and 2, and ‘definitely heterogeneous’ if Hi is greater than 2. Also, they stated that the H1 statistic compared to the statistics H2 and H3 has a discriminating power and is the most important heterogeneity measure (Rao & Hamed 2000).
Estimation of regional parameters
In order to determine the heterogeneity measure, it is necessary to estimate the regional parameters of KAP distribution for the sample data. In this study, three models are considered to fulfill this aim.
is the jth observed severity of the site i; is the jth computed severity of the site i; and is the average of observed severity of the site i.
), the at-site parameters are combined to give regional parameters () as:Model 2: Based on this model, after combining the sites, the regional parameters of these sites are estimated using the proposed methods. For calculating the cost function, the combined data are used instead of at-site data.
Due to the nature of the defined optimization problem, the regional parameters of KAP distribution can only be estimated by the ACSS method as a meta-heuristic algorithm.
Goodness-of-fit measure
The coefficient Ak was given by Hosking & Wallis (1997, Table A3).
Estimation of the regional distribution
Cluster analysis
The cluster analysis is used to assemble objects into a set of specific groups with a maximum similarity between the cluster members (Modarres 2010). Generally, the clustering algorithms can be classified into two types: hierarchical and partitional (Demirel et al. 2012). Among the different partitional algorithms, the commonest and most well-known method is the K-means method (Sönmez & Kömüşcü 2011; Kar et al. 2012; Rahman et al. 2013).
) is the average of the silhouette widths over all objects (Kaufman & Rousseeuw 1990).Adaption to droughts
Since droughts are regional in nature, it is logical to utilize RFA for drought assessment (Tallaksen et al. 1997). Various indices have been introduced to define drought characteristics, and commonly, the standardized precipitation index (SPI) is the most popular and widely used one for drought analysis (Raziei et al. 2009; Serinaldi et al. 2009; Santos et al. 2011; Mirabbasi et al. 2012; Chen et al. 2013). Among different time scales, the SPI values on a six-month time scale (SPI-6) is most useful for describing the shallow soil moisture available to agricultural crops (Reddy & Singh 2014; Abdi et al. 2016a, 2016b). After computing the SPI-6 values, the drought characteristics can be obtained from these values. In this study, we considered the drought severity (Sd) as one of the main characteristics of drought. The drought severity is defined based on the absolute value of cumulative SPI values within the duration of a drought event (Chen et al. 2013). The drought duration is defined as the number of consecutive intervals (months) where the SPI values are less than the threshold value of zero (Mirabbasi et al. 2012).
In order to adopt the LMOM approach to the drought severity, the details of the methodology can be summarized as follows:
Consider a number of sites within the region.
Assemble the monthly precipitation data of each site.
Calculate the SPI values and extract the drought characteristics.
Compute the discordancy and heterogeneity measures.
Specify the homogeneous sub-regions by using the cluster analysis.
Determine the best regional distribution for each sub-region.
Study area and data
Considered attributes in this study
| Attribute | Statistic | Range |
|---|---|---|
| Annual precipitation | Average (mm) | 199.9 to 508.3 |
| Maximum (mm) | 306.0 to 1,101.7 | |
| Minimum (mm) | 70.1 to 287.3 | |
| Coefficient of variation | 0.2 to 0.6 | |
| Drought severity | Average | 4.2 to 6.4 |
| Maximum | 14.5 to 84.7 | |
| Minimum | 0.0 to 0.1 | |
| Coefficient of variation | 0.7 to 2.9 | |
| L-CV | 0.4 to 0.8 | |
| L-CS | 0.1 to 0.8 | |
| L-CK | −0.1 to 0.7 | |
| Geographical coordination | Latitude (N) | 37.1 to 39.1 |
| Longitude (E) | 45.4 to 47.9 |
| Attribute | Statistic | Range |
|---|---|---|
| Annual precipitation | Average (mm) | 199.9 to 508.3 |
| Maximum (mm) | 306.0 to 1,101.7 | |
| Minimum (mm) | 70.1 to 287.3 | |
| Coefficient of variation | 0.2 to 0.6 | |
| Drought severity | Average | 4.2 to 6.4 |
| Maximum | 14.5 to 84.7 | |
| Minimum | 0.0 to 0.1 | |
| Coefficient of variation | 0.7 to 2.9 | |
| L-CV | 0.4 to 0.8 | |
| L-CS | 0.1 to 0.8 | |
| L-CK | −0.1 to 0.7 | |
| Geographical coordination | Latitude (N) | 37.1 to 39.1 |
| Longitude (E) | 45.4 to 47.9 |
Geographical location of the study area and sites.
RESULTS
Assumption of homogeneous region for all sites
First, the entire the study area is assumed as one homogeneous region. Hence, the classic discordancy (CDi) measure and robust discordancy (RDi) measure with five different values of α are computed and the discordant sites are shown in Table 2. As can be seen, sites 2, 43, 51, and 56 are the discordant sites according to the values of CDi and RDi. The difference between CDi and RDi with α = 0.75, 0.95 is only at site 13, and with α= 0.80, 0.85, 0.90, it is at sites 3 and 13.
The discordant sites for the study area
| Discordancy measure | |||||
|---|---|---|---|---|---|
| Robust | |||||
| Classic | α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 |
| 2(3.45), 43(5.05), 51(3.24), 56(3.11)a | 2(3.20), 13(3.67), 43(4.92), 51(3.88), 56(4.08) | 2(3.66), 3(3.06), 13(3.80), 43(4.83), 51(3.93), 56(3.85) | 2(3.82), 3(3.40), 13(4.10), 43(4.78), 51(4.18), 56(3.82) | 2(3.82), 3(3.40), 13(4.10), 43(4.78), 51(4.18), 56(3.82) | 2(3.53), 13(3.12), 43(4.87), 51(3.13), 56(3.73) |
| Discordancy measure | |||||
|---|---|---|---|---|---|
| Robust | |||||
| Classic | α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 |
| 2(3.45), 43(5.05), 51(3.24), 56(3.11)a | 2(3.20), 13(3.67), 43(4.92), 51(3.88), 56(4.08) | 2(3.66), 3(3.06), 13(3.80), 43(4.83), 51(3.93), 56(3.85) | 2(3.82), 3(3.40), 13(4.10), 43(4.78), 51(4.18), 56(3.82) | 2(3.82), 3(3.40), 13(4.10), 43(4.78), 51(4.18), 56(3.82) | 2(3.53), 13(3.12), 43(4.87), 51(3.13), 56(3.73) |
aSite number (Di value).
Then, the regional parameters of the KAP distribution are estimated based on three proposed models. For model 1, the results of at-site parameters and the values of two main goodness-of-fit tests, the root mean square error (RMSE) and efficiency coefficient (EC), based on the LMOM and ACSS methods, are presented in Tables 3 and 4 for the selected sites, respectively. Tables 3 and 4 indicate that the results of the ACSS method are better than the LMOM method for all sites, especially for sites 13, 51, and 60 where the EC value becomes smaller than 0.9 for LMOM, whereas the performances of ACSS are adequate for these sites. Also, the values of ξ and α parameters of the LMOM results for sites 22, 42, 52, and 58 are significantly different in comparison with these parameters for the other sites. Therefore, the regional parameters can be extensively affected by the weighted ratio () of model 1. In contrast, the values of these parameters have some small variations for the results of the ACSS method.
The at-site parameters and goodness-of-fit tests obtained by the LMOM method for the selected sites
| Parameters of KAP distribution | Goodness-of-fit tests | |||||
|---|---|---|---|---|---|---|
| Site number | ξ | α | k | h | EC | RMSE |
| 13 | −0.3371 | 0.3914 | −0.6539 | 2.0483 | 0.8579 | 1.0119 |
| 22 | −15.8244 | 19.6832 | 0.9210 | 4.3378 | 0.9809 | 0.1882 |
| 42 | −21.7919 | 23.9107 | 0.8254 | 5.4338 | 0.9853 | 0.1879 |
| 51 | −0.3343 | 0.3368 | −0.6962 | 2.2837 | 0.8293 | 1.1716 |
| 52 | −41.4305 | 95.4353 | 2.1477 | 3.5142 | 0.9515 | 0.2266 |
| 58 | −16.3784 | 25.8720 | 1.2713 | 3.6308 | 0.9894 | 0.1179 |
| 60 | −3.9450 | 3.6644 | 0.2142 | 3.4678 | 0.8728 | 0.5099 |
| Parameters of KAP distribution | Goodness-of-fit tests | |||||
|---|---|---|---|---|---|---|
| Site number | ξ | α | k | h | EC | RMSE |
| 13 | −0.3371 | 0.3914 | −0.6539 | 2.0483 | 0.8579 | 1.0119 |
| 22 | −15.8244 | 19.6832 | 0.9210 | 4.3378 | 0.9809 | 0.1882 |
| 42 | −21.7919 | 23.9107 | 0.8254 | 5.4338 | 0.9853 | 0.1879 |
| 51 | −0.3343 | 0.3368 | −0.6962 | 2.2837 | 0.8293 | 1.1716 |
| 52 | −41.4305 | 95.4353 | 2.1477 | 3.5142 | 0.9515 | 0.2266 |
| 58 | −16.3784 | 25.8720 | 1.2713 | 3.6308 | 0.9894 | 0.1179 |
| 60 | −3.9450 | 3.6644 | 0.2142 | 3.4678 | 0.8728 | 0.5099 |
The at-site parameters and goodness-of-fit tests obtained by the ACSS method for the selected sites
| Parameters of KAP distribution | Goodness-of-fit tests | |||||
|---|---|---|---|---|---|---|
| Site number | ξ | α | k | h | EC | RMSE |
| 13 | −0.3525 | 0.2005 | −1.0000 | 3.1951 | 0.9975 | 0.1334 |
| 22 | −4.7943 | 5.7087 | 0.5174 | 2.9539 | 0.9863 | 0.1595 |
| 42 | −7.4693 | 7.8411 | 0.5289 | 3.6934 | 0.9885 | 0.1658 |
| 51 | −0.7334 | 0.2095 | −1.0000 | 4.9999 | 0.9988 | 0.0996 |
| 52 | −3.9892 | 7.5960 | 1.0636 | 2.0900 | 0.9585 | 0.2096 |
| 58 | −5.9095 | 8.9804 | 0.9055 | 2.6776 | 0.9897 | 0.1163 |
| 60 | −8.9238 | 9.9841 | 0.6695 | 3.9585 | 0.9314 | 0.3744 |
| Parameters of KAP distribution | Goodness-of-fit tests | |||||
|---|---|---|---|---|---|---|
| Site number | ξ | α | k | h | EC | RMSE |
| 13 | −0.3525 | 0.2005 | −1.0000 | 3.1951 | 0.9975 | 0.1334 |
| 22 | −4.7943 | 5.7087 | 0.5174 | 2.9539 | 0.9863 | 0.1595 |
| 42 | −7.4693 | 7.8411 | 0.5289 | 3.6934 | 0.9885 | 0.1658 |
| 51 | −0.7334 | 0.2095 | −1.0000 | 4.9999 | 0.9988 | 0.0996 |
| 52 | −3.9892 | 7.5960 | 1.0636 | 2.0900 | 0.9585 | 0.2096 |
| 58 | −5.9095 | 8.9804 | 0.9055 | 2.6776 | 0.9897 | 0.1163 |
| 60 | −8.9238 | 9.9841 | 0.6695 | 3.9585 | 0.9314 | 0.3744 |
The regional parameters of KAP distribution based on model 1 are given in Table 5. The results of the goodness-of-fit tests by using the regional parameters of model 1 for each site demonstrate that LMOM performs worse than ACSS. For 97% of sites, the results of EC are less than zero for LMOM while for ACSS only 8% of sites have EC < 0 and in return 65% and 30% of sites have EC greater than 0.5 and 0.8, respectively.
The regional parameters of the KAP distribution based on the three models
| Parameters of the KAP distribution | |||||
|---|---|---|---|---|---|
| Model number | Method | ξ | α | K | h |
| 1 | LMOM | −3.3795 | 5.0803 | 0.1223 | 2.3682 |
| ACSS | −2.3047 | 2.9404 | 0.0998 | 2.338 | |
| 2 | LMOM | −1.0803 | 1.4655 | −0.0731 | 2.0612 |
| ACSS | −1.0179 | 1.4671 | −0.0613 | 1.9648 | |
| 3 | ACSS | −0.7107 | 1.4442 | 0.0256 | 1.6082 |
| CSS | −0.6209 | 1.3793 | 0.0109 | 1.5196 | |
| PSO | −0.4925 | 1.2782 | −0.0167 | 1.3974 | |
| GA | −0.5197 | 1.2441 | −0.0371 | 1.4939 | |
| Parameters of the KAP distribution | |||||
|---|---|---|---|---|---|
| Model number | Method | ξ | α | K | h |
| 1 | LMOM | −3.3795 | 5.0803 | 0.1223 | 2.3682 |
| ACSS | −2.3047 | 2.9404 | 0.0998 | 2.338 | |
| 2 | LMOM | −1.0803 | 1.4655 | −0.0731 | 2.0612 |
| ACSS | −1.0179 | 1.4671 | −0.0613 | 1.9648 | |
| 3 | ACSS | −0.7107 | 1.4442 | 0.0256 | 1.6082 |
| CSS | −0.6209 | 1.3793 | 0.0109 | 1.5196 | |
| PSO | −0.4925 | 1.2782 | −0.0167 | 1.3974 | |
| GA | −0.5197 | 1.2441 | −0.0371 | 1.4939 | |
The regional parameters of KAP distribution for model 2 are given in Table 5. The results of the goodness-of-fit tests obtained by the regional parameters show that the performance of the LMOM and ACSS methods is close. However, the value of the sum of differences of the EC (SDEC) measure is equal to 0.4631, which shows the relative superiority of ACSS. In both methods, the value of EC for site 43 is less than zero. As shown in Table 2, site 43 is a discordant site and has the maximum value of the discordancy measure among the others. Therefore, this model is less sensitive than the first one to the outliers. In addition, the results demonstrate that the percentage of sites which have EC values greater than 0.5 and 0.8, obtained by LMOM are equal to 97% and 72%, respectively. These values are 97% and 75% for ACSS, respectively.
Table 5 presents the regional parameters of the KAP distribution estimated from model 3. Based on these parameters, the value of the SDEC measure for this model is equal to 1.8301. The results of the ACSS method for this model are more accurate than those for the second model. The value of EC in model 3 for site 43 is equal to 0.3345, whereas less than zero was obtained for the previous model. Indeed, this indicates that the sensitivity of model 3 to the discordant sites is less than model 2. The values of EC obtained by the ACSS results for model 3 show that 98% and 87% of the sites have values greater than 0.5 and 0.8, respectively. Consequently, model 3 yields the global optimal parameters, while the optimum parameters of model 2 are local.
The mean values of the at-site EC and RMSE values using the regional parameters based on model 3
| Method | Best value of cost function | Mean of EC | Mean of RMSE |
|---|---|---|---|
| ACSS | 7.6924 | 0.8718 | 0.5023 |
| CSS | 7.6958 | 0.8717 | 0.5028 |
| PSO | 7.7079 | 0.8715 | 0.5029 |
| GA | 7.7391 | 0.8710 | 0.5033 |
| Method | Best value of cost function | Mean of EC | Mean of RMSE |
|---|---|---|---|
| ACSS | 7.6924 | 0.8718 | 0.5023 |
| CSS | 7.6958 | 0.8717 | 0.5028 |
| PSO | 7.7079 | 0.8715 | 0.5029 |
| GA | 7.7391 | 0.8710 | 0.5033 |
Comparison of the convergence histories for model 3 obtained by different meta-heuristic methods.
Comparison of the convergence histories for model 3 obtained by different meta-heuristic methods.
According to the results of various models and methods, model 3 with the ACSS method is selected for estimating the regional parameters of the KAP distribution. However, with regional parameters of the KAP distribution, sample sizes the same as at-site historical data and Nsim = 1,000, the heterogeneity and goodness-of-fit measures are calculated as shown in Table 7, which indicate the study area with 60 sites is ‘definitely heterogeneous’ before removal of the discordant sites according to the H1 value. Also, based on the results of goodness-of-fit measure (ZDIST), only the PE3 distribution is selected for the study area.
The heterogeneity and goodness-of-fit measures for the study area
| After removing discordant sites based on discordancy measure | ||||||||
|---|---|---|---|---|---|---|---|---|
| Measure | Before removing discordant sites | Classic | Robust | |||||
| α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 | ||||
| Heterogeneity | H1 | 3.3 | 2.13 | 1.51 | 1.29 | 1.29 | 1.29 | 1.51 |
| H2 | 5.06 | 4.12 | 3.49 | 3.33 | 3.33 | 3.33 | 3.49 | |
| H3 | 6.68 | 5.58 | 4.66 | 4.53 | 4.53 | 4.53 | 4.66 | |
| Goodness-of-fit | GPA | −3.8 | −4.71 | −4.68 | −4.62 | −4.62 | −4.62 | −4.68 |
| GEV | −7.28 | −8.54 | −8.48 | −8.36 | −8.36 | −8.36 | −8.48 | |
| GLO | −8.01 | −9.73 | −9.3 | −9.18 | −9.18 | −9.18 | −9.3 | |
| LN3 | −4.39 | −5.64 | −5.4 | −5.38 | −5.38 | −5.38 | −5.4 | |
| PE3 | 0.51 | 0.53 | −0.45 | −0.31 | −0.31 | −0.31 | −0.45 | |
| After removing discordant sites based on discordancy measure | ||||||||
|---|---|---|---|---|---|---|---|---|
| Measure | Before removing discordant sites | Classic | Robust | |||||
| α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 | ||||
| Heterogeneity | H1 | 3.3 | 2.13 | 1.51 | 1.29 | 1.29 | 1.29 | 1.51 |
| H2 | 5.06 | 4.12 | 3.49 | 3.33 | 3.33 | 3.33 | 3.49 | |
| H3 | 6.68 | 5.58 | 4.66 | 4.53 | 4.53 | 4.53 | 4.66 | |
| Goodness-of-fit | GPA | −3.8 | −4.71 | −4.68 | −4.62 | −4.62 | −4.62 | −4.68 |
| GEV | −7.28 | −8.54 | −8.48 | −8.36 | −8.36 | −8.36 | −8.48 | |
| GLO | −8.01 | −9.73 | −9.3 | −9.18 | −9.18 | −9.18 | −9.3 | |
| LN3 | −4.39 | −5.64 | −5.4 | −5.38 | −5.38 | −5.38 | −5.4 | |
| PE3 | 0.51 | 0.53 | −0.45 | −0.31 | −0.31 | −0.31 | −0.45 | |
The LMOM ratio diagram in terms of L-CV and L-CS.
The LMOM ratio diagram in terms of L-CS and L-CK.
At last, after removing the discordant sites based on the classic and robust discordancy measures, the heterogeneity and goodness-of-fit measures are computed as demonstrated in Table 7. A value of H1 = 2.13 in this table indicates that the study area is ‘definitely heterogeneous’. On the contrary, using the robust discordancy measure with different values of subset factor for recognizing and removing the discordant sites shows an improvement in heterogeneity measure and suggests that the whole region is ‘possibly heterogeneous’. Also, based on the goodness-of-fit measure, the PE3 distribution is selected as the regional distribution for this region.
Clustering and identification of homogeneous regions
In order to identify the optimum number of clusters, the average silhouette width for some number of clusters are computed and illustrated in Table 8. According to this table, the maximum value of average silhouette width corresponds to the best number of clusters (k), which is equal to 3.
The average silhouette width for a number of clusters
| Number of clusters (k) | Average silhouette width |
|---|---|
| 2 | 0.4315 |
| 3 | 0.5122 |
| 4 | 0.4907 |
| 5 | 0.4704 |
| 6 | 0.4815 |
| 7 | 0.5050 |
| Number of clusters (k) | Average silhouette width |
|---|---|
| 2 | 0.4315 |
| 3 | 0.5122 |
| 4 | 0.4907 |
| 5 | 0.4704 |
| 6 | 0.4815 |
| 7 | 0.5050 |
The values in bold indicate the optimum number of cluster and its average silhouette width.
The silhouette plot.
The location of sites for three clusters.
The classic discordancy measure and robust discordancy measure with five different values of α are computed and the discordant sites for three clusters are presented in Table 9. The heterogeneity and goodness-of-fit measures are demonstrated in Table 10, as well. The results in Table 10 show that the third cluster is ‘acceptably homogeneous’, the first cluster is ‘possibly heterogeneous’, and the second one is ‘definitely heterogeneous’ according to the H1 test.
The discordant sites for the all clusters
| Discordancy measure | ||||||
|---|---|---|---|---|---|---|
| Robust | ||||||
| Cluster | Classic | α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 |
| 1 | – | 23, 37, 43 | 37, 43 | 37, 43 | 43 | 43 |
| 2 | 56 | 13, 51, 56 | 13, 51, 53, 56 | 13, 51, 56 | 13, 51, 56 | 53, 56 |
| 3 | – | 42, 48 | 42, 48 | 42, 48 | 42, 48 | 42 |
| Discordancy measure | ||||||
|---|---|---|---|---|---|---|
| Robust | ||||||
| Cluster | Classic | α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 |
| 1 | – | 23, 37, 43 | 37, 43 | 37, 43 | 43 | 43 |
| 2 | 56 | 13, 51, 56 | 13, 51, 53, 56 | 13, 51, 56 | 13, 51, 56 | 53, 56 |
| 3 | – | 42, 48 | 42, 48 | 42, 48 | 42, 48 | 42 |
The heterogeneity and goodness-of-fit measures for all clusters
| After removing discordant sites based on discordancy measure | |||||||||
|---|---|---|---|---|---|---|---|---|---|
| Classic | Robust | ||||||||
| Measure | Cluster | Before removing discordant sites | α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 | ||
| Heterogeneity | H1 | 1 | 1.58 | 1.58 | −0.93 | −0.70 | −0.70 | 0.29 | 0.29 |
| 2 | 2.14 | 2.03 | 1.37 | 0.87 | 1.37 | 1.37 | 1.63 | ||
| 3 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01 | ||
| H2 | 1 | 2.09 | 2.09 | 0.50 | 0.59 | 0.59 | 1.24 | 1.24 | |
| 2 | 3.44 | 3.55 | 2.75 | 2.19 | 2.75 | 2.75 | 3.13 | ||
| 3 | 0.85 | 0.85 | 0.85 | 0.85 | 0.85 | 0.85 | 0.85 | ||
| H3 | 1 | 3.25 | 3.25 | 1.90 | 1.72 | 1.72 | 2.36 | 2.36 | |
| 2 | 5.10 | 5.31 | 3.93 | 3.49 | 3.93 | 3.93 | 5.09 | ||
| 3 | 1.10 | 1.10 | 1.10 | 1.10 | 1.10 | 1.10 | 1.10 | ||
| Goodness-of-fit | GPA | 1 | −2.48 | −2.48 | −2.52 | −2.58 | −2.58 | −2.52 | −2.52 |
| 2 | −2.88 | −3.25 | −3.66 | −3.71 | −3.66 | −3.66 | −3.13 | ||
| 3 | −1.93 | −1.93 | −1.93 | −1.93 | −1.93 | −1.93 | −1.93 | ||
| GEV | 1 | −5.15 | −5.15 | v4.81 | −4.86 | −4.86 | −4.81 | −4.81 | |
| 2 | −5.07 | −5.50 | −6.00 | −5.96 | −6.00 | −6.00 | −5.17 | ||
| 3 | −3.84 | −3.84 | −3.84 | −3.84 | −3.84 | −3.84 | −3.84 | ||
| GLO | 1 | −5.88 | −5.88 | −5.46 | −5.49 | −5.49 | −5.40 | −5.40 | |
| 2 | −5.44 | −5.88 | −6.44 | −6.37 | −6.44 | −6.44 | −5.49 | ||
| 3 | −4.26 | −4.26 | −4.26 | −4.26 | −4.26 | −4.26 | −4.26 | ||
| LN3 | 1 | −3.57 | −3.57 | −3.50 | −3.50 | −3.50 | −3.33 | −3.33 | |
| 2 | −2.87 | −3.23 | −3.83 | −3.79 | −3.83 | −3.83 | −3.02 | ||
| 3 | −2.33 | −2.33 | −2.33 | −2.33 | −2.33 | −2.33 | −2.33 | ||
| PE3 | 1 | −0.87 | −0.87 | −1.25 | −1.19 | −1.19 | −0.80 | −0.80 | |
| 2 | 0.83 | 0.61 | −0.15 | −0.13 | −0.15 | −0.15 | 0.60 | ||
| 3 | 0.23 | 0.23 | 0.23 | 0.23 | 0.23 | 0.23 | 0.23 | ||
| After removing discordant sites based on discordancy measure | |||||||||
|---|---|---|---|---|---|---|---|---|---|
| Classic | Robust | ||||||||
| Measure | Cluster | Before removing discordant sites | α = 0.75 | α = 0.80 | α = 0.85 | α = 0.90 | α = 0.95 | ||
| Heterogeneity | H1 | 1 | 1.58 | 1.58 | −0.93 | −0.70 | −0.70 | 0.29 | 0.29 |
| 2 | 2.14 | 2.03 | 1.37 | 0.87 | 1.37 | 1.37 | 1.63 | ||
| 3 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01 | ||
| H2 | 1 | 2.09 | 2.09 | 0.50 | 0.59 | 0.59 | 1.24 | 1.24 | |
| 2 | 3.44 | 3.55 | 2.75 | 2.19 | 2.75 | 2.75 | 3.13 | ||
| 3 | 0.85 | 0.85 | 0.85 | 0.85 | 0.85 | 0.85 | 0.85 | ||
| H3 | 1 | 3.25 | 3.25 | 1.90 | 1.72 | 1.72 | 2.36 | 2.36 | |
| 2 | 5.10 | 5.31 | 3.93 | 3.49 | 3.93 | 3.93 | 5.09 | ||
| 3 | 1.10 | 1.10 | 1.10 | 1.10 | 1.10 | 1.10 | 1.10 | ||
| Goodness-of-fit | GPA | 1 | −2.48 | −2.48 | −2.52 | −2.58 | −2.58 | −2.52 | −2.52 |
| 2 | −2.88 | −3.25 | −3.66 | −3.71 | −3.66 | −3.66 | −3.13 | ||
| 3 | −1.93 | −1.93 | −1.93 | −1.93 | −1.93 | −1.93 | −1.93 | ||
| GEV | 1 | −5.15 | −5.15 | v4.81 | −4.86 | −4.86 | −4.81 | −4.81 | |
| 2 | −5.07 | −5.50 | −6.00 | −5.96 | −6.00 | −6.00 | −5.17 | ||
| 3 | −3.84 | −3.84 | −3.84 | −3.84 | −3.84 | −3.84 | −3.84 | ||
| GLO | 1 | −5.88 | −5.88 | −5.46 | −5.49 | −5.49 | −5.40 | −5.40 | |
| 2 | −5.44 | −5.88 | −6.44 | −6.37 | −6.44 | −6.44 | −5.49 | ||
| 3 | −4.26 | −4.26 | −4.26 | −4.26 | −4.26 | −4.26 | −4.26 | ||
| LN3 | 1 | −3.57 | −3.57 | −3.50 | −3.50 | −3.50 | −3.33 | −3.33 | |
| 2 | −2.87 | −3.23 | −3.83 | −3.79 | −3.83 | −3.83 | −3.02 | ||
| 3 | −2.33 | −2.33 | −2.33 | −2.33 | −2.33 | −2.33 | −2.33 | ||
| PE3 | 1 | −0.87 | −0.87 | −1.25 | −1.19 | −1.19 | −0.80 | −0.80 | |
| 2 | 0.83 | 0.61 | −0.15 | −0.13 | −0.15 | −0.15 | 0.60 | ||
| 3 | 0.23 | 0.23 | 0.23 | 0.23 | 0.23 | 0.23 | 0.23 | ||
In order to improve the heterogeneity measure, Table 9 demonstrates that the classic discordancy measure could not detect any discordant site for cluster 1. Also, only one discordant site is detected for cluster 2, whereas the different discordant sites are identified based on the robust discordancy measure with different values of α. Finally, after removing the discordant sites for clusters 1 and 2, the heterogeneity and goodness-of-fit measures are also shown in Table 10. In this condition, the value of H1 demonstrates that cluster 2 is not ‘acceptably homogeneous’. On the contrary, after removing the discordant sites (based on robust discordancy measure with all values of α), cluster 1 is found to be ‘acceptably homogeneous’. Just α = 0.80 yields ‘acceptably homogeneous’ for cluster 2 and the other values of α yield ‘possibly heterogeneous’.
The LMOM ratio diagram in terms of L-CS and L-CK for three clusters.
After the selection of PE3 distribution as a suitable regional frequency distribution, the regional parameters of the respective distribution are estimated using model 3 and by the ACSS method, then the regional quantiles of various return periods for all clusters are illustrated in Table 11.
Estimation of regional parameters and quantiles for the PE3 distribution
| Cluster | Regional parameters | Regional quantiles | ||||||
|---|---|---|---|---|---|---|---|---|
| γ | α | β | T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | |
| 1 | 1.0000 | 2.7212 | 0.1482 | 9.92 | 15.88 | 21.52 | 28.01 | 35.26 |
| 2 | 1.0000 | 3.2201 | 0.1524 | 11.61 | 19.12 | 26.21 | 34.33 | 43.38 |
| 3 | 1.0000 | 3.1939 | 0.1461 | 12.04 | 19.86 | 27.28 | 35.84 | 45.41 |
| Cluster | Regional parameters | Regional quantiles | ||||||
|---|---|---|---|---|---|---|---|---|
| γ | α | β | T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | |
| 1 | 1.0000 | 2.7212 | 0.1482 | 9.92 | 15.88 | 21.52 | 28.01 | 35.26 |
| 2 | 1.0000 | 3.2201 | 0.1524 | 11.61 | 19.12 | 26.21 | 34.33 | 43.38 |
| 3 | 1.0000 | 3.1939 | 0.1461 | 12.04 | 19.86 | 27.28 | 35.84 | 45.41 |
The estimated regional growth curve values and the relative RMSE for each cluster are determined in Table 12. The results verify that by increasing the value of the return period, the accuracy of the regional growth curve decreases. In addition, Monte Carlo simulation is used to assess the accuracy of the estimates of the determined homogeneous regions, by generating the same sample sizes as those of the original sample. The number of random samples, Nsim, is set to 10,000 for selected quantiles corresponding to return periods of 10, 25, 50, 100, and 200 years. Table 13 shows the results of the BR(F), AR(F), and RR(F). According to Tables 12 and 13, the values of the RMSE of the estimated growth curve values are always less than the estimated quantiles for the same return periods, which represent the uncertainty in the estimated scale factor. These results are verified by Saf (2009) and Hussain (2011).
The regional growth values and the relative RMSE values
| Cluster | Regional growth curve | RMSE | ||||||||
|---|---|---|---|---|---|---|---|---|---|---|
| T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | |
| 1 | 2.01 | 3.22 | 4.37 | 5.69 | 7.16 | 0.1174 | 0.1582 | 0.1820 | 0.2021 | 0.2196 |
| 2 | 2.25 | 3.70 | 5.07 | 6.65 | 8.40 | 0.1464 | 0.1951 | 0.2252 | 0.2516 | 0.2752 |
| 3 | 2.17 | 3.57 | 4.91 | 6.45 | 8.17 | 0.1368 | 0.1677 | 0.1792 | 0.1860 | 0.1901 |
| Cluster | Regional growth curve | RMSE | ||||||||
|---|---|---|---|---|---|---|---|---|---|---|
| T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | |
| 1 | 2.01 | 3.22 | 4.37 | 5.69 | 7.16 | 0.1174 | 0.1582 | 0.1820 | 0.2021 | 0.2196 |
| 2 | 2.25 | 3.70 | 5.07 | 6.65 | 8.40 | 0.1464 | 0.1951 | 0.2252 | 0.2516 | 0.2752 |
| 3 | 2.17 | 3.57 | 4.91 | 6.45 | 8.17 | 0.1368 | 0.1677 | 0.1792 | 0.1860 | 0.1901 |
Summary of simulation results
| Measure | Cluster | Return period (years) | ||||
|---|---|---|---|---|---|---|
| T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | ||
| BR(F) | 1 | −0.1280 | −0.2256 | −0.2604 | −0.2791 | −0.2892 |
| 2 | −0.1583 | −0.2580 | −0.2906 | −0.3070 | −0.3150 | |
| 3 | −0.1453 | −0.2601 | −0.3026 | −0.3271 | −0.3420 | |
| AR(F) | 1 | 0.1435 | 0.2280 | 0.2604 | 0.2791 | 0.2892 |
| 2 | 0.1634 | 0.2629 | 0.3014 | 0.3243 | 0.3408 | |
| 3 | 0.1540 | 0.2601 | 0.3026 | 0.3271 | 0.3420 | |
| RR(F) | 1 | 0.2056 | 0.2753 | 0.3052 | 0.3237 | 0.3358 |
| 2 | 0.2306 | 0.3059 | 0.3384 | 0.3591 | 0.3730 | |
| 3 | 0.2311 | 0.3067 | 0.3405 | 0.3611 | 0.3739 | |
| Measure | Cluster | Return period (years) | ||||
|---|---|---|---|---|---|---|
| T = 10 | T = 25 | T = 50 | T = 100 | T = 200 | ||
| BR(F) | 1 | −0.1280 | −0.2256 | −0.2604 | −0.2791 | −0.2892 |
| 2 | −0.1583 | −0.2580 | −0.2906 | −0.3070 | −0.3150 | |
| 3 | −0.1453 | −0.2601 | −0.3026 | −0.3271 | −0.3420 | |
| AR(F) | 1 | 0.1435 | 0.2280 | 0.2604 | 0.2791 | 0.2892 |
| 2 | 0.1634 | 0.2629 | 0.3014 | 0.3243 | 0.3408 | |
| 3 | 0.1540 | 0.2601 | 0.3026 | 0.3271 | 0.3420 | |
| RR(F) | 1 | 0.2056 | 0.2753 | 0.3052 | 0.3237 | 0.3358 |
| 2 | 0.2306 | 0.3059 | 0.3384 | 0.3591 | 0.3730 | |
| 3 | 0.2311 | 0.3067 | 0.3405 | 0.3611 | 0.3739 | |
CONCLUSIONS
In the standard CSS algorithm, search agents or CPs affect each other according to laws from electrostatics and mechanics. Here, a new ACSS is developed to improve the performance of CSS. The ACSS and CSS algorithms use the governing laws from electrical physics and Newtonian mechanics to determine the amount and the direction of a CPs movement. The CSS and ACSS algorithms can distinguish the local search and the global search phases, and utilize suitable relationships in these phases resulting in a good balance between exploration and exploitation. However, applying some controlling limits and the constraints is the main difference between the ACSS and CSS algorithms.
On the other hand, a new model is introduced to estimate the regional parameters of the statistical distributions. Then, the ACSS method is utilized to estimate the accurate parameters and quantile estimates for a region. Also for the Monte Carlo simulation, this new model is employed to generate samples from the selected distribution. In addition, the robust discordancy measure with different subset factors is developed to detect and remove the discordant sites, which causes some improvements in the heterogeneity measure. The data used in this study are drought severity extracted for 60 sites in East-Azarbaijan, Iran. The relevant conclusions can be summarized as follows:
The ACSS has the ability to improve the search and convergence speed compared to the standard CSS.
The ACSS algorithm has stable convergence characteristics and good computational ability, and is an effective algorithm compared to other algorithms, such as the GA and PSO, for RFA.
The ACSS can estimate the global optimum (or at least near it) of parameters without using nonlinear complex equations and any approximations.
The new model is the best model among the others. This model is based on the participation of all the data of sites in a region for estimating the regional parameters and quantiles, while other models are based on the weighted at-site parameters or the combination of the sites.
Identifying discordant sites using the robust discordancy measure is more accurate than the classic discordancy measure. Among the different values for the subset factor, the value of 0.8 is found suitable for robust discordancy measure and to detect the discordant sites of a heterogeneous region.
The study area is divided into three sub-regions, south, northeast, and northwest by using the K-means clustering method and the Pearson type 3 distribution is found the best regional distribution for all sub-regions.
Finally, this method can be applied for RFA of other hydrological phenomena, such as rainfall and stream flow. Also, the multivariate RFA framework based on the multivariate LMOM approach is suggested for future application of this work. In addition, another goodness-of-fit measure, the LMOM ratios-based test, can be applied and compared with the ZDIST-based test for selecting the best regional distribution.
%20cropped.png?versionId=5948)