A fuzzy hybrid clustering method for identifying hydrologic homogeneous regions

ANN

Artificial Neural Networks

CN

Curve Number

DEM

Digital Elevation Model

FA

Factor Analysis

FCM

Fuzzy C-means Algorithm

HHR

Hydrologic Homogeneous Regions

HLR

Hydrologic Landscape Regions

HRU

Hydrological Response Units

KMO

Kaiser–Meyer–Olkin statistic

Ksat

Soil saturated hydraulic conductivity

MARD

Mean Absolute Relative Distance

NDVI

Normalized Difference Vegetation Index

P25

The 24-hour rainfall with return period of 25 years

P50

The 24-hour rainfall with return period of 50 years

P100

The 24-hour rainfall with return period of 100 years

PCA

Principal Component Analysis

REA

Representative Elementary Area

REW

Representative Elementary Watersheds

SOFM

Self-Organizing Features Maps

TPI

Topographic Position Index

VAR

Variable

WGEW

Walnut Gulch Experimental Watershed

In recent decades, water resources studies at watershed scale have been enhanced through the deployment of geographic information systems, satellite imagery data, and water-related softwares. In this regard, spatially distributed analysis and tools are becoming vital for hydrologists and watershed managers. Accordingly, classification and clustering of the watershed landscape into homogeneous regions facilitates prioritization of watershed management measures.

In recent decades, hydrological spatial discretization has also enjoyed new terminologies such as representative elementary area (REA) by Wood et al. (1988), hydrological response units (HRU) by Flügel (1995), representative elementary watersheds (REWs) by Reggiani et al. (1998, 1999, 2000), hydrologic landscape regions (HLRs) by Wolock et al. (2004), REWs for cold and snowy regions by Tian et al. (2006), a flexible methodology for discretization into REWs by Dehotin & Braud (2008) and Kuentz et al. (2017). In the meantime, clustering algorithms such as Ward algorithm (Ward 1963), self-organizing feature mapping (SOFM) method (Kohonen 1982), and fuzzy c-means (FCM) clustering algorithm (Bezdek 1981) have been used in numerous applications.

The FCM method has been widely adopted in natural sciences (e.g., Bruin & Stein 1998), discretization of natural watersheds into groups (e.g., Rao & Srinivas 2006), water resources (e.g., Golian et al. 2010; Irwin 2015), and regionalization studies (e.g., Basu & Srinivas 2015). One of the advantages of FCM is the ability to determine the cells/clusters uncertainty through confusion index. However, the main FCM weakness is the way to determine the FCM's parameters, i.e., the fuzzification parameter and the optimum number of clusters. Therefore, fuzzy clustering, similar to other clustering methods, may not be used without determination of the optimum number of clusters. Fuzzy clustering, FCM in particular, has been recommended to partition watersheds into natural landscapes or groups with the homogeneous hydrologic response by different researchers (e.g., Bruin & Stein 1998; Rao & Srinivas 2008).

The SOFM has also been used in numerous fields, including in water resources applications (e.g., Nourani & Parhizkar 2013; Nourani et al. 2016; Sharghi et al. 2018). This method is considered as a useful way to determine homogeneous regions as it enjoys high performance in visualizing and summarizing the contributing features as well as the ability to display the distribution of each component (Farsadnia et al. 2014).

Hybrid clustering methods, for example, the combination of SOFM and FCM, have been successfully used in the regionalization of sites or catchments (e.g., Srinivas et al. 2008; Farsadnia et al. 2014; Ahani & Nadoushani 2016). The SOFM-FCM combination can identify homogeneous regions more quickly and accurately while having sufficient flexibility to minimize the FCM objective function.

The main objective of this study is to derive hydrologic homogeneous regions (HHRs) through the adoption of effective rainfall-runoff characteristics. This is achieved via application of a fuzzy hybrid clustering method that will be evaluated through new criteria. In this regard, in addition to the application of well-known and conventional clustering validation approaches traditionally developed and used for FCM, two and one technique(s) are proposed for determining the optimum number of clusters and the fuzzification parameter, respectively. Furthermore, with emphasis on the outcome of this study, the advantages of integrating FCM with the SOFM are outlined.

Factor analysis (FA) is one of the early multivariable statistical approaches that are able to determine the data structure by summarizing and reducing the data size (Suhr 2006). This method is frequently applied in climatic (e.g., Bartzokas et al. 2003) and hydrologic studies (e.g., Farsadnia et al. 2014). FA method is fully described by Harman (1976) and Rencher (1995).

FA can improve the effectiveness of SOFM classification, leading to an increase in accuracy of cluster centers. In this study, clustering is performed based on the most effective variables resulting from the FA pre-processing.

Natural geographical landscapes, unlike the man-made domains, often do not have a clear border. Therefore, the use of fuzzy clustering may be more justified than hard clustering (crisp) in the discretization of terrestrial and, in particular, hydrological regions.

FCM is performed through an iterative procedure of optimizing the OF by updating the degree of membership u_ij (Equation (5)) and the cluster centers c_i (Equation (6)). The procedure has been described in detail by Bezdek (1981).

In order to ensure proper clustering, special attention needs to be paid to: (1) determine the fuzzification parameter (m); (2) select the degree of membership for the winning cluster (degree of membership threshold); (3) control clustering uncertainty; (4) validate fuzzy clustering to optimize the number of clusters; and ultimately (5) measure the similarity of cells to its cluster center as a second control.

Chiu (2005) proposed two α-cut (α_low, α_high) for the degree of membership based on the number of clusters (c). If α_high≥ 1−1/c, then that cell definitely belongs to the cluster, and if α_low <1/c, the membership of a cell to the cluster is decisively rejected. Therefore, the value of 1/c is believed to be an acceptable choice for the threshold fuzzy membership (Rao & Srinivas 2008). In this study, the cells which firmly belong to each cluster are determined by specifying a threshold fuzzy membership value equal to 1/c (c= 4, … ,11). Therefore, the cells in the cluster which have the memberships greater than the specified threshold value can produce a fuzzy cluster.

For qualitative evaluation of this parameter, the CI values in the interval [0,1] may be divided into five classes denoting insignificant (0–0.2), very low (0.2–0.4), low (0.4–0.6), high (0.6–0.8), and very high (0.8–1) uncertainty.

Validation in fuzzy clustering involves maximization of two important criteria: first, the density of the data in each cluster and second, the distance between the clusters, which indirectly leads to the optimization of the number of clusters. Various fuzzy cluster validation indices have been proposed to determine optimum c in a data set (Pal & Bezdek 1995). In the present study, five cluster validity indices, as presented in Table 1, are examined.

Nonlinearity and flexibility of artificial neural networks (ANN) are attractive features in a variety of applications in hydrology (Govindaraju & Rao 2000). A particular category of ANN, known as SOFM (Kohonen 1982), has been used as a clustering tool in order to convert the complex nonlinear statistical relationship among high-dimensional input data into a simple and geometric relationship on a low-dimensional display (Nourani et al. 2015).

SOFMs are one of the unsupervised training networks such that no output is required to classify the existing data, yet they convert complex multi-dimensional data into visible clusters. SOFMs attempt to find the topological structure in the input data by mapping the available data on an attribute mapping, which may also be considered as the output layer (Rao & Srinivas 2008).

In the absence of clear understandable patterns in existing data, the interpretation of SOFM clusters, regardless of size and dimensions, is rarely possible. However, in such cases, SOFMs may be considered as a useful tool for supervised clustering algorithms (Rao & Srinivas 2008).

In this study, SOFM method was used to provide the initial cluster centers for FCM algorithm. The combination of SOFM and FCM is denoted as SOMFCM, that was assembled using R features including ‘fclust’ (Ferraro & Giordani 2015), ‘kohonen’ (Wehrens & Buydens 2007), ‘cluster’ (Maechler et al. 2015), and ‘e1071’ (Meyer et al. 2015).

Generally speaking, determination of the number of clusters before clustering is a disadvantage of most clustering techniques (Irwin 2015). Furthermore, the main disadvantage of FCM algorithm is its sensitivity to the fuzzification parameter (m) and the number of clusters (c), both of which should be optimized simultaneously. As m increases, the degree of membership of a given cell to a cluster decreases, causing overlap of the membership of the cell to adjacent clusters in boundary regions (the overlapping of adjacent clusters leads to the gradual display of variations in nature). In contrast, as c increases, the similarity of cells within a cluster to its cluster center increases and the entropy in each cluster decreases, while the average cell values of the membership degree also decrease.

Accordingly, the following SOMFCM validation techniques are adopted in order to determine m and c.

• Technique I (optimum m): As fuzzification parameter (m) increases, the V_PC (Equation (8))/V_PE (Equation (9)) decreases/increases, so that the intersection point of simultaneous changes of the two for different values of m and c may be identified as the optimum m. Rao & Srinivas (2008) also suggested the indirect use of V_PC and V_PE in order to determine the range of m and its optimum value in a fuzzy clustering context.

• Technique II (optimum c): If the number of clusters (c) increases, the distance between the scaled values of the cells and that of the cluster center (expressed as the mean Euclidean distance, Euc_mean) decreases, whereas the mean confusion index (CI_mean) of the cells increases. Therefore, the interaction of these two (namely, CI_mean and Euc_mean) for different values of c and m may yield the optimum c.

• Technique III (optimum c): The correlation coefficient R², defined as the ratio of the variance of the data corresponding to the cluster centers to the variance of the actual data, may be considered as one criterion in determining the optimum c:

  

  (13)

The Walnut Gulch Experimental Watershed (WGEW) was selected as a research area by the United States Department of Agriculture (USDA) in the mid-1950s. The WGEW with an area of 150 square kilometers is located in southeastern Arizona, USA, at 31°42′ N latitude and 110°03′ W longitude.

The elevation of the watershed ranges between 1,221 m and 1,887 m. The climate is classified as semi-arid, with a mean annual temperature of 17.7 °C and mean annual precipitation of 350 mm. The precipitation regime is dominated by the North American Monsoon with slightly more than 60% of the annual total occurring during July–September (USDA 2007).

The available literature reports on various aspects of the watershed including Goodrich et al. (2008) on precipitation, Stone et al. (2008) on runoff, Skirvin et al. (2008) on vegetation, Keefer et al. (2008) on meteorology and soil, and Heilman et al. (2008) on GIS maps. The quality of watershed data has been assured and that was the motivation for choosing this watershed in the current study.

Hydrological, meteorological, soil and land use data are available online at http://www.tucson.ars.ag.gov/dap. Figure 1 shows the digital elevation model (DEM), along with the distribution of 26 selected rain gauges within or close to WGEW.

The selection of input data required for discretization of HHRs depends on the studied hydrological process. In this study, rainfall-runoff processes are the main focus, thus effective layers must deal with meteorological attributes as well as watershed physiographic/surface features. The scale of the maps was 1:24,000 and the selected cell size was 150 meters for all watershed input layers. The characteristics of the input layers entered in factor analysis are shown in Table 2.

The flow chart of the study procedure is presented in Figure 2.

Results are presented in three sections: input data preprocessing, factor analysis, and clustering.

In the present study, the selected cell size was 150 (m) for all input layers. This selection was based on the recommendation of Bloschl & Sivapalan (1995) in order to consider compatibility between the scale of observed data with scale and purpose of the simulation of hydrological processes. On the other hand, as hilly lands cover most parts of the WGEW and the river network is relatively dense, use of cell sizes larger than 150 m may introduce uncertainty in the spatial aggregation. Some selected input features are shown in Figure 3.

Since rainfall is the most important factor in rainfall-runoff events, it must be taken into account in the prioritization of the regions in terms of potential runoff generation (Saghafian & Khosroshahi 2005). In this study, rainfall layers corresponding to 25-, 50- and 100-year return periods were generated using rainfall frequency analysis of available 50-year rainfall time series (1967–2016) for 26 selected stations inside and outside the watershed (as shown in Figure 1).

According to the results of rainfall frequency analysis, the fitted distribution functions for 26 rain gauges were Log Pearson 3 (14 stations), Lognormal (5 stations), Gamma (5 stations), and Gumbel (2 stations). The 24-hour rainfall with return periods of 25-, 50- and 100-year was estimated for each rain gauge. Figure 4, for example, shows the 50-year rainfall map.

The total number of cells entering factor analysis was equal to 6,225 for each component/variable. Then, each of the ten attributes was rescaled by Equation (1), thus forming a 6,225 × 10 matrix. To determine whether factor analysis was an appropriate tool for data reduction, Kaiser–Meyer–Olkin (KMO) statistic of sampling adequacy was computed. The KMO statistic equaled 0.922, that confirmed the applicability of FA in this particular case study. The value of KMO also points to a significant correlation among components.

The characteristics of principal components are presented in Table 3 including eigenvalues, variance proportion, and cumulative variance proportion. Each component whose variance is at least equal to 1 (representing at least 10% of the total variance) is selected as the principal component (Wolock et al. 2004). This threshold of the variance is a qualitative criterion and is flexible depending on the subject or user experience. According to Table 3, it is clear that the first three components account for 98.49% of total variance of input variables.

Each variable with the highest loadings for each component was selected as the most effective variable/feature. According to the rotated component matrix (Table 4), TPI, aspect, curvature-plan, curvature-profile, and slope are the most effective on the first component that accounts for more than 82.52% of input variance. One of these five variables should be selected as the representative of the first component. Although all five variables are DEM derivatives, TPI indicates the tendency of each cell to saturate. Selection of TPI is reasonable in terms of runoff generation process. Furthermore, Distance and CN have the most effect on the second and third components that account for 11.08% and 4.87% of input variance, respectively. In comparison, Dehotin & Braud (2008) considered DEM, TPI, and CN layers to determine the HHR without preprocessing.

Thus, TPI, distance, and CN were identified as the most effective surface features of the watershed along with the distributed rainfall layer. The effect of using the third layer (CN) on clustering results will be further examined.

Before performing fuzzy clustering (FCM algorithm), it is necessary to initialize cluster centers. For this purpose, a 6,225 × 4 matrix, involving three surface features/variables and a rainfall variable, were subjected to SOFM and the results were evaluated based on the silhouette width (Rousseeuw 1987).

A hexagonal lattice was recommended by Kohonen (2001) since it does not favor horizontal or vertical directions. In this study, a matrix of nodes with a rectangular or hexagonal lattice neighborhood structure was used. Therefore, the best configuration of SOFM was evaluated through several runs. The best architecture including the number of clusters (c) and lattice neighborhood structure based on the silhouette width are presented in Table 5.

The best network arrangement for each number of clusters (as presented in Table 5) was entered into the fuzzy clustering. It should be noted that results of using Ward algorithm, though slightly weaker, did not differ significantly with those of the SOFM.

Based on the best type of lattice neighborhood structure (Table 5) corresponding to the different number of clusters (4 ≤ c ≤ 11), the sensitivity to variation of the fuzzification parameter (1.1 ≤ m ≤ 3) with an increment of 0.1 was examined. Pal & Bezdek (1995) reported that the FCM provides better performance for m in 1.5–2.5 range. However, in this study, the SOMFCM hybrid fuzzy clustering was performed in 160 possible runs (i.e., 20 different states for 1.1 ≤ m ≤ 3 and eight different states for 4 ≤ c ≤ 11) for each proposed scenario. To determine the optimum c, the results of the five selected validation indices introduced in available studies (e.g., Rao & Srinivas 2006, 2008; Odgers et al. 2011; Farsadnia et al. 2014; Prasad & Arora 2014) are presented in Table 6. It should be noted that for each clustering scenario, the upper bound of m was reduced due to divergent non-meaningful results. In Table 6, the case of c= 6 and m= 1.6 is marked in bold type in order to facilitate comparison.

In Table 6, Euc_mean is the mean Euclidean distance between the actual feature values of the cells and the ones corresponding to cluster center, R² represents the ratio of the variance of the data corresponding to the cluster centers to the variance of the actual data, and OF is the objective function. Since for each clustering scenario, a threshold fuzzification parameter (m) corresponding to N_Umax = 0 is defined and the number of cells with a membership degree greater than the specified value in the present study (i.e., 1−1/c) will be equal to zero (N_Umax = 0), the results of the validation indices are presented to the threshold m. It is also worth noting that, following Chiu (2005), the cells with a degree of membership less than 1/c are not assigned to any cluster.

The optimum c corresponds to the minimum values of V_EXB, V_K, V_FS, and V_PE as well as the maximum V_PC. The values of the selected clustering validation indices in this study (as evident by Table 6) indicate that three-component clustering scenario outperforms the four-component clustering scenario. Accordingly, the fourth component or third surface feature (i.e., Distance in this study) did not improve the results. Figure 5 shows the variations of the selected validation indices for different values of m and c corresponding to the three-component clustering (TPI, CN, and P50).

As shown in Figure 5(a)–5(c), it is not possible to determine optimum c and m based on V_PC, V_PE, and V_FS. The V_FS follows a slow decreasing trend as the number of clusters increases over 8 while the value of V_PC has a monotonic decreasing trend. However, the V_PE has a monotonic increasing trend as the number of clusters increases. Overall, the indices (V_PC and V_PE) suggest c= 4 as the optimum c in this study. This is because neither V_PC and V_PE have a direct connection to any property of the input data (Zhang et al. 2008) and cannot solely serve as a reasonable basis for choosing optimum c. Consequently, these indices seem to be more appropriate when well-separated clusters with crisp boundaries are formed.

On the other hand, as shown in Figure 5(d) and 5(e), the trends of V_k and V_EXB are similar so that the minimum values obtained from these indices may be the basis for choosing c= 7 clusters while m falls in the 1.5–1.7 range. However, V_FS, V_EXB, and V_k variations show monotonic decreasing or increasing the tendency for majority m and c values, even though there are a few exceptions. Hence, no decisive judgment is expected on the optimum c on the basis of these indices. Our results are similar to those of Rao & Srinivas's (2006) study in that the selected indices may lead to the poor discretization of clusters. Therefore, it is necessary to determine optimum c and m through other means. Based on Table 6, one may infer the following:

• The range of R² varies from 0.16 to 0.76 for all four-component clustering scenarios and from 0.41 to 0.84 for the three-component clustering scenario. This indicates that clustering with three input components is better than the four input components.

• Comparing the mean CI between the three- and four-component clustering scenarios for rainfall shows that the three-component scenario involving 50-year rainfall is more reliable (mean CI values in the three-component scenario are less than the ones corresponding to all four-component scenarios for different return periods).

• For all clustering scenarios, the number of cells (N_Umax) belonging to the clusters with a certain degree of membership (as α_high ≥1−1/c) has a decreasing trend with increasing m and c.

• N_Umax for clustering scenario 4 was greater than the ones in other clustering scenarios for all values of m and c, which suggests that the clusters (regions) may be obtained without the distance layer. Therefore, distance layer has no significant effect on the results and is removed in further analysis.

Due to the inadequacy of conventional validation indices (Table 6 and Figure 5), three validation techniques are proposed for determining the optimum c and m. In all three techniques, the standardization of the inputs in [0,1] interval was first performed for both vertical axes.

Technique I of determining optimum m: The results of the application of technique I corresponding to c= 6 are shown in Figure 6.

Similar results, as shown in Table 7, were obtained for the different clustering scenarios and different number of clusters (4 ≤ c ≤ 11). In most cases of c, the point of V_PC and V_PE intersection occurred at m= 1.6 (except Figure 6). Therefore, technique I suggests the optimum m= 1.6 for most scenarios.

Technique II of determining optimum c:Figure 7, as an example, shows the interaction between CI_mean and Euc_mean for optimum m= 1.7.

As shown in Figure 7, there are three intersections near c equal to 5, 8, and 10. The criterion for determination of optimum c is the highest frequency of intersections. Therefore, by drawing the interactions similar to Figure 7 for different m and c and for different clustering scenarios, the optimum c may be determined. The results are presented in Table 7.

Technique III of determining optimum c: As an example, Figure 8 shows the interaction for m= 1.7 (as optimum m).

As shown in Figure 8, in this case, the intersection occurred at only one point, close to c= 6. In technique III, there is one intersection point for each interaction while the criterion for optimum is the highest frequency of intersections for different interactions. Therefore, by drawing the interactions similar to Figure 8 for different m and c, and for different clustering scenarios, optimum c may be obtained as presented in Table 7.

Both optimum c and m were selected based on the highest frequency for intersection points. There were two reasons for limiting upper bound of m for each clustering scenario: (1) a threshold for m corresponding to N_Umax = 0 and (2) the mean threshold for CI which is equal to 0.6. The results were subjected to these thresholds.

Based on Table 7, the highest priority corresponds to c= 6. In other words, the WGEW watershed may be divided into six HHR. Therefore, out of 160 different selections in SOMFCM, the optimum values of m= 1.6 and c= 6 predominated in most clustering scenarios. Now, the best clustering scenario should be selected.

Other interpretations derived from Tables 6 and 7 are as follows:

• Techniques I (optimum m) and III (optimum c) yield similar results in different clustering scenarios. The results of technique II are somewhat different, because of variations in mean CI and mean Euclidean distance for different values of m.

• Technique III in all clustering scenarios presents similar results, which indicates the accuracy of the selected clustering features.

• Rainfall return periods only affect technique II, especially CI_mean and Euc_mean parameters. However, the results of techniques I and III are similar in all return periods.

• It may be concluded that in the present study, the distance from the nearest channel (distance layer) has no significant effect on optimum m and c. In other words, clustering can be performed with two more effective land layers, namely, TPI and CN. As such, clustering requires less data, which is desirable in the delineation of HHR.

To better compare two superior clustering scenarios (i.e., scenarios 1 and 4) as far as the distance layer is concerned, the changes in optimal OF value are presented in Figure 9.

Generally speaking, the OF value reduces when c increases for a specified value of m. According to Figure 9, after the removal of the distance layer as the third surface feature, the OF value is significantly reduced. Therefore, in order to achieve minimum OF in FCM algorithm, two effective surface features are adequate for clustering.

Figure 9 and Tables 6 and 7 indicate that clustering scenario 4 involving TPI, CN, and P50 is the best. Clustering maps corresponding to the best three-component and four-component scenarios are shown in Figures 10 and 11. Although m = 1.7 was optimum in the three-component scenario based on technique I, other indices, such as mean CI and mean Umax, were better produced for m = 1.6 (second rank in technique I). Statistical comparison between cluster cells and cluster centers, as well as MARD, are presented in Tables 8 and 9 for four-component and three-component clustering scenarios, respectively.

As presented in Tables 8 and 9, comparison of features cell values with those of cluster centers indicates that there is a significant similarity between cells and cluster centers in each cluster. MARD values imply that four-component clustering error is less than that of the three-component clustering. This is because in four-component clustering, the addition of one input component to the clustering procedure improves the accuracy of assigning cells to the clusters. On the other hand, adding one more component increases the computation time and a relative error reduction of 10.22% (3.13–2.81/3.13). This error is negligible, so the third surface feature (distance layer) has no significant effect on the results.

As shown in Figure 11, in addition to its remarkable similarity with the CN map shown in Figure 2, the discretized regions (clusters) also identify the river network, plains, and man-made areas.

Furthermore, Tables 8 and 9 indicate that the area of the smallest cluster is 4.44% (6.66 km²) and 4% (6 km²) corresponding to the four-component and three-component clustering scenarios, respectively. Such values, as reported by Boulaine (1980), should be less than the smallest distinguishable area on a map (expressed in this study as the threshold area) that is equal to 0.5 cm × 0.5 cm = 0.25 cm² on the map. With scale 1:24,000 used in this study, the threshold area is 120 m × 120 m = 14,400 (m²), almost 1.44% of the total basin area (2.16 km²). Therefore, the aforementioned smallest areas for both clustering scenarios are larger than the threshold area.

The CI_mean values were obtained as 0.477 and 0.326, indicating that the uncertainty of clusters falls in low and very low categories for the four-component and three-component clustering scenarios, respectively. The average degree of membership is equal to 0.684 for the four-component clustering and 0.788 for the three-component clustering scenarios, which are slightly less than α_High = (1–1/6) = 0.833. These results indicate that the three-component clustering map is more reliable and accurate than the four-component map.

The results of the best clustering scenarios with different inputs derived from the conventional fuzzy clustering validation indices along with the proposed techniques are presented in Table 10.

According to Table 10, the three-component clustering scenario (TPI-CN-P50), shown in bold, relies on the validation techniques proposed in the study and presents the most reliable clustering map (Figure 11). Some physical justifications are as follows:

• By comparing the final cluster map (Figure 11) with the 50-year rainfall map (Figure 4), it is evident that clusters 1 and 2 (together covering 43% of the watershed area), and clusters 4 and 5 (together covering 32% of the watershed area) correspond to areas which receive maximum and minimum 50-year rainfall, respectively. It is worth noting that the dominant rainfall occurs in the eastern half of the WGEW (within clusters 1, 2, and 3).

• By comparing the final cluster map (Figure 11) with the CN map (Figure 3(b)) it is indicated that cluster 5 is approximately representative of maximum CN areas. Clusters 2, 3, and 4 (61% of the watershed area) cover areas with minimum CN.

• No clear pattern may be observed between the final cluster map (Figure 11) with the TPI and Distance maps (Figure 3(a) and 3(c)).

The objective of this study was to provide an efficient methodology for identification/delineation of hydrologic homogeneous regions (HHR) in the context of the rainfall-runoff framework. For this purpose, a hybrid fuzzy clustering method, namely SOMFCM, was adopted. Specifying initial cluster centers in the FCM algorithm affects the final clusters. In this paper, a hybrid approach was proposed to improve the performance of FCM algorithm for partitioning of watersheds. SOMFCM is a combination of SOFM (for determining the initial centers and the location of fuzzy clusters) and the robust and the well-known FCM algorithm (for watershed discretization). At first, a number of attributes/components associated with rainfall-runoff processes were identified and mapped. The components consisted of ten physiographic/surface features and one of three rainfall values (with return periods of 25, 50 and 100 years) representing a meteorological factor. All maps were prepared in 150-meter cell size. All surface features were entered into the factor analysis (FA) in order to reduce the dimensions of the input layers. The most effective features were identified as TPI, CN, and Distance that encompass physiographical and soil/vegetation attributes within the clustering. Therefore, these three layers along with one of the three rainfall layers were forwarded into the SOMFCM clustering procedure.

SOMFCM cluster maps were produced for 1.1 ≤ m ≤ 3 and 4 ≤ c ≤ 11 in two clustering scenarios involving four-component (with distance layer) and three-component cases (without distance layer). Then, to validate the SOMFCM clustering results, in addition to conventional and well-known fuzzy validation indices, a number of techniques for determining the optimum number of clusters (c) and fuzzification parameter (m) were introduced. The main conclusions of this study are as follows:

• Although using the maximum degree of membership is a common way to determine the winning cluster, an acceptable level of α-cut is required. In this study, the criterion of (α-cut= 1−1/c= 0.833, for c = 6) was applied. Although α-cut= 0.833 value was not met for all cells, best results were obtained for the three-component clustering scenario (TPI-CN-P50). The number of clusters was determined in such a way that the mean degree of membership of all cells was 0.788 (slightly less than 0.833), while 3,261 cells (52% of all cells) had a degree of membership more than 0.833. The choice of this criterion for assigning the cells into the clusters in the context of fuzzy clustering methods appears logical. The CI_mean is another validation index that turned out to be 0.326, which is the minimum value among the results of other clustering scenarios as well as within the range of very low uncertainty (CI∈(0.2–0.4)).

• Application of the selected conventional validation indices (V_PC, V_PE, V_FS, V_EXB, and V_K) and the proposed validation techniques showed that among 160 different clustering runs for each scenario, optimum clustering involves six clusters (c= 6) with fuzzification parameter (m) of 1.6, derived from the TPI, CN, and P50 features.

• It became clear that use of conventional validation indices, despite having a strong mathematical and rational theory, will not necessarily lead to an optimum solution with acceptable performance in the delineation of HHRs. On the other hand, the techniques introduced in the study may be more effective in determining the optimum number of hydrologically homogeneous clusters.

• The performance of FCM algorithm was improved when coupled with SOFM. In general, the quality of clustering is improved via coupled SOFM-FCM. Therefore, use of the SOMFCM with a fewer number of input layers and properly prioritized via factor analysis is recommended in clustering the HHRs.

Comparison of the HHR map resulted from clustering for WGEW with a distributed map showing the basin runoff to provide a proper judgment on the effectiveness of the clustering practice is suggested for future research. Also, the application and evaluation of the proposed methodology on other watersheds using the more effective layers introduced in this study can also be suggested. In addition, the assessment of other effective watershed features and the optimal combination of various features in preparation of a clustering map along with sensitivity analysis of SOMFCM with respect to DEM cell size may be conducted in future research.