ABSTRACT
The clustering of small watersheds based on hydrological similarity serves as an effective technique for identifying similarities in watershed runoff generation and routing conditions. It also addresses the challenge of parameter transplantation in undocumented areas. In this research, 545 small watersheds in hilly areas within Shandong Province were studied using 22 selected indicators to represent their climate and underlying surface characteristics. The study employed a two-stage clustering method combining the self-organizing map (SOM) neural network with the K-means algorithm, facilitating the classification of these watersheds into various groups. Each group of small watersheds was then analyzed for its unique characteristics. To validate the reasonableness of the classification results, the flood peak modulus of each watershed was calculated using a hydrologic-hydraulic method, while a parameter transplantation study was carried out and generalized for the clustering results. The findings indicate that the SOM-K-means clustering method efficiently classified the watersheds into 12 similar groups, validating its effective application in small watershed classification. This classification assists in solving the problem of flood forecasting in the ungauged watersheds in Shandong Province and developing more effective flood risk management strategies.
HIGHLIGHTS
The study focused on 545 watersheds in the hilly areas of Shandong Province.
A technical system of 'small watershed division-watershed similarity classification-small watershed hydrological model construction-parameter transplantation' has been proposed.
The basic features and hydrological characteristics of each group's subwatersheds were analyzed, and flood modulus was used to verify the reasonableness of the results.
INTRODUCTION
Carrying out hydrological forecasting in hilly watersheds can effectively prevent flooding and reduce disaster losses, and is the key to solving floods in hilly areas (Xue 2016). Currently, imperfect runoff monitoring networks and delayed initiation of monitoring activities result in scarce or nonexistent runoff data in most hilly watersheds. In addition, rapid socio-economic developments have led to significant changes in underlying surface conditions, further complicating hydrological forecasting in ungauged areas (Zhu et al. 2020). Thus, hydrological forecasting in these ungauged areas remains a significant challenge.
In recent years, more and more researchers and scholars have adopted the parameter regionalization method based on basin hydrological similarity to carry out hydrological forecasting in ungauged areas, which is to derive the model parameters of ungauged basins from the model parameters of the gauged basins by some means. Watershed hydrologic similarity is usually defined as the watersheds with similar underlying surface characteristics, driving forces, hydrodynamic conditions, etc., and satisfy the transplantation of hydrologic model parameters between each other (Hrachowitz et al. 2013; Wu et al. 2023). The primary methods of parameter regionalization include the spatial proximity, physical attribute similarity, and regression method (Wu et al. 2023). The regression method establishes the regression relationship between the watershed attributes and the hydrological model parameters and applies it to ungauged watersheds, but this method is prone to the phenomenon of equifinality, which limits its practical application. The spatial proximity method is generally based on the spatial distance between the watersheds to determine the reference watershed for parameter transplantation, which includes spatial distance averaging and spatial interpolation (Guo et al. 2020). However, its applicability is reduced in areas with diverse watershed characteristics (Sun et al. 2023a, b).
In contrast, the physical attribute similarity method is based on the premise that watersheds with similar physical characteristics will exhibit similar hydrological responses (Yi et al. 2014). The method has the characteristics of wide applicability and small limitations, increasingly making it a focal point for solving hydrological forecasting problems in ungauged watersheds. It typically employs various clustering algorithms to identify and group watersheds with similar characteristics (Wu et al. 2023). These clustering algorithms (Hassan & Ping 2012; Zhou et al. 2014) are different from classification methods (Cao et al. 2020) and belong to unsupervised learning algorithms without training data. Notable examples include hierarchical clustering as demonstrated by the analytic hierarchy process (Fan & Liu 2015; Alves et al. 2023), division-based clustering using the K-means algorithm (Clare & Cohen 2001), the fuzzy clustering algorithm (Yu & Li 2014), density clustering such as DBSCAN (Ram et al. 2010), model clustering through the self-organizing map (SOM) neural network (Kohonen 1982), graph clustering represented by spectral clustering (Liu et al. 2023), etc. These algorithms find extensive utilizations in watershed classification and have achieved desirable results. For instance, Yi et al. (2014) classified watersheds using the SOM and hierarchical clustering analysis (HCA), later employing the Hydrologiska Byråns Vattenbalansavdelning (HBV) model for parameter transposition in ungauged watersheds. Similarly, Yang et al. (2022) used the K-means algorithm based on physical similarity to group 64 small watersheds in hilly areas into 11 similar groups and conducted tests for hydrological model parameter transplantation. Merz & Blöschl (2004) focused on hydrological forecasting through watershed similarity, and Mayer et al. (2014) applied cluster analysis to classify the Great Lakes basin in the United States. SOM neural networks also now have more applications in hydrogeology. Varouchakis et al. (2022) utilized the SOM neural network to process the relevant data and calculate the groundwater stress index, proving its effectiveness and reliability in large-scale applications, as well as combining SOM with geostatistics to provide a new method for estimating the spatial distribution of groundwater hydrology in complex hydrogeological systems (Varouchakis et al. 2023).
However, inherent limitations in single clustering algorithms, such as K-means, might result in suboptimal clustering effects. For instance, the selection of initial clustering centers in the K-means algorithm is typically random, significantly influencing the final results. Furthermore, the SOM neural network is difficult to explain the clustering results better. To overcome these issues, some researchers have combined the strengths of the SOM and K-means algorithms, applying this integrated approach across various fields. This synthesis has proven to enhance the effectiveness of clustering results significantly, as evidenced by studies such as those by Huang et al. (2022) and Bigdeli et al. (2022). In this paper, aiming at the difficult problem of watershed similarity division feature index construction and classification algorithm, we constructed a feature index library considering physical index and attribute index, synthesized the two-stage clustering algorithm, and proposed a small watershed similarity evaluation method for coupling two-stage clustering model.
At present, there are limited studies on large-scale hydrological similarity analysis for specific regions combined with machine learning. In this paper, a hydrological similarity analysis framework based on machine learning for watershed clustering-categorization determination is proposed and applied to Shandong Province, China. Shandong Province is a large economic province in northern China, has an important strategic position, is an important part of the Northeast Asian Economic Circle, and serves as a vital link between the Beijing–Tianjin–Hebei area and the Yangtze River Delta region. This geographical and economic significance underscores the importance of ungauged watershed flood forecasting for the region. Consequently, the focus of this study is to analyze the hydrological similarity of 545 upstream watersheds in Shandong's hilly area, based on their physical characteristics, using a combined approach of the SOM neural network and K-means method for clustering. Based on the clustering results, two similar basins were selected for parameter transplantation, and their flood processes were simulated using the HEC-HMS model to validate the applicability of this method in hydrological forecasting for ungauged watersheds. The study also used hydrologic-hydraulic methods to calculate the flood peak modulus of these small watersheds, providing a means to verify the effectiveness of the clustering. Finally, based on the watershed clustering results, a supervised machine learning classification technique (Random Forest Classification Model) (Jonathan et al. 2022) was introduced for classification determination and parameter transplantation of other small watersheds in Shandong Province, confirming that the method can be used to solve the problem of flood forecasting for any ungauged watersheds in Shandong Province. The results of this research not only clarify the hydrological characteristics of various small watersheds in Shandong Province but also simplify their study and facilitate the implementation of flood risk management strategies. It can also be used to identify reference watersheds and provide an effective tool for solving the flood forecasting problems in ungauged watersheds within Shandong Province, as well as providing some experiences and insights into flood forecasting in other ungauged areas.
The remainder of the paper is structured as follows: Section 2 describes the study area and the data sources used. Section 3 details the methodology and the research process. Section 4 presents the results and includes a discussion of these findings, and Section 5 concludes.
STUDY AREA AND DATA
Shandong Province, situated in northeastern China, is part of the North China Plain coastal region, adjacent to the east of the Taihang Mountains. The province spans approximately 155,800 km², with diverse terrain. The central region is dominated by the Taiyi Mountain Range, characterized by elevated terrain. In contrast, the southwest and northwest regions are notably flat, while the eastern section comprises gently undulating hills. Particularly, over 40% of Shandong's territory comprises mountainous and hilly areas, which are prone to frequent flooding that threatens the lives and property of tens of millions of residents across approximately 12 cities. Shandong's water systems are expansive, encompassing parts of the Yellow River, Huai River, and the Sea basins, along with the Southern Four Lakes Region. The climate exhibits a warm temperate profile with monsoonal influences, marked by synchronized seasonal rainfall and temperature, brief spring and autumn, and extended winter and summer periods.
General topography of typical hilly areas of Shandong Province, including the research area.
General topography of typical hilly areas of Shandong Province, including the research area.
METHOD
Methodological framework
SOM neural network
The SOM neural network is an unsupervised learning technique capable of transforming high-dimensional data into a more manageable low-dimensional space. This transformation is achieved through self-organizing mapping processes. The model segregates sample points into different discrete areas in accordance with the similarity of their data features, thereby achieving clustering. The SOM neural network mimics the two-dimensional spatial structure of neurons and incorporates self-organization and learning capabilities through mutual competition and interaction among the neurons. Its applications in pattern recognition and classification are well regarded, owing to its proficiency in automatic clustering, robust nonlinear mapping, and high fault tolerance. However, it is also known for its high complexity and the inability to provide precise clustering information post-clustering (Kohonen 1982).
The network structure of the SOM neural network comprises an input layer and a competitive layer. The neurons in the input layer are connected to the weights of the output neurons in the competitive layer. There exists localized connectivity between neurons within a certain vicinity in the competitive layer. The network manages the response weights of the neurons to the input patterns and the lateral inhibition among neurons, continuously adjusting these weights throughout the model training process. The neuron nodes in the model are topologically connected, and this connectivity reflects the strength of associations between the nodes (Bigdeli et al. 2022).
K-means clustering algorithm



One critical aspect of the K-means algorithm is the need to predefine the cluster centers and the number of clusters. This algorithm is highly sensitive to the initial choice of cluster centers. Different initial selections can lead to varied clustering outcomes, impacting the stability and accuracy of the results. In addition, the K-means algorithm is susceptible to falling into local optimization, which can further affect the reliability of the clustering (Clare & Cohen 2001).
SOM+K-means two-stage clustering approach
To address the challenges associated with the K-means algorithm and the SOM neural network as previously mentioned, this study integrates these two algorithms to form a two-stage SOM-K-means clustering model. The initial classification of clusters is determined by the SOM algorithm. The centroids generated by SOM clustering serve as the initial values for the K-means algorithm. The final clustering results are derived from the K-means algorithm, utilizing internal indicators to optimize the number of clusters (Han & Tang 2021; Huang et al. 2022).
Flood peak modulus


Selection of similarity indicators
In this study, various characteristics were selected as similarity indicators, including multi-year average flood-season rainfall, watershed area, watershed length, average watershed slope, average elevation, topographic relief, watershed shape factor, river network density, main stream gradient, main stream length, stream tortuosity, average topographic index (Jin et al. 2017), normalized difference vegetation index (NDVI), land-use types (percentage area in grasslands, buildings, farmland, woodland, and waters), and soil texture. The soil texture was categorized into four types (Wang et al. 2016). Table 1 shows the eigenvalues for each similarity indicator across all small watersheds. In this paper, each indicator is considered to have the same impact, i.e., each indicator is equally weighted.
Eigenvalues of characterization indicators for small watersheds in the study area
Feature . | Indicator . | Max . | Min . | Mean . |
---|---|---|---|---|
Underlay | Area/km2 | 359.55 | 42.14 | 98.65 |
Average slope/° | 21.43 | 1.24 | 6.48 | |
Average elevation/m | 554.55 | 3.65 | 156.55 | |
Topographic relief/m | 1,404.00 | 11.00 | 375.50 | |
River network density/(km/km2) | 0.96 | 0.08 | 0.39 | |
Main stream gradient/% | 3.27 | 0.01 | 0.36 | |
Main stream length/km | 66.24 | 4.69 | 19.60 | |
Watershed length/km | 49.78 | 8.55 | 17.54 | |
Watershed shape factor | 0.91 | 0.06 | 0.35 | |
Average topographic index | 7.28 | 4.36 | 5.96 | |
Stream tortuosity | 2.53 | 1.00 | 1.31 | |
Land use | Grassland/% | 62.26 | 0.00 | 10.29 |
Buildings/% | 89.29 | 0.79 | 15.72 | |
Farmland/% | 93.76 | 1.05 | 60.02 | |
Woodland/% | 62.96 | 0.00 | 11.46 | |
Waters/% | 31.60 | 0.00 | 1.98 | |
NDVI | 0.844 | 0.376 | 0.713 | |
Soil texture type | Sandy or Loamy sandy soils, Sandy loam soils/% | 98.20 | 0.00 | 33.26 |
Clay loam soils, Clay soils, Silt clay loam soils, Sandy clay soils/% | 100.00 | 0.00 | 23.44 | |
Sandy clay loam soils/% | 100.00 | 0.00 | 35.29 | |
Silt loam soils/% | 100.00 | 0.00 | 7.40 | |
Climate | Multi-year average flood-season rainfall/mm | 225.82 | 133.43 | 184.79 |
Feature . | Indicator . | Max . | Min . | Mean . |
---|---|---|---|---|
Underlay | Area/km2 | 359.55 | 42.14 | 98.65 |
Average slope/° | 21.43 | 1.24 | 6.48 | |
Average elevation/m | 554.55 | 3.65 | 156.55 | |
Topographic relief/m | 1,404.00 | 11.00 | 375.50 | |
River network density/(km/km2) | 0.96 | 0.08 | 0.39 | |
Main stream gradient/% | 3.27 | 0.01 | 0.36 | |
Main stream length/km | 66.24 | 4.69 | 19.60 | |
Watershed length/km | 49.78 | 8.55 | 17.54 | |
Watershed shape factor | 0.91 | 0.06 | 0.35 | |
Average topographic index | 7.28 | 4.36 | 5.96 | |
Stream tortuosity | 2.53 | 1.00 | 1.31 | |
Land use | Grassland/% | 62.26 | 0.00 | 10.29 |
Buildings/% | 89.29 | 0.79 | 15.72 | |
Farmland/% | 93.76 | 1.05 | 60.02 | |
Woodland/% | 62.96 | 0.00 | 11.46 | |
Waters/% | 31.60 | 0.00 | 1.98 | |
NDVI | 0.844 | 0.376 | 0.713 | |
Soil texture type | Sandy or Loamy sandy soils, Sandy loam soils/% | 98.20 | 0.00 | 33.26 |
Clay loam soils, Clay soils, Silt clay loam soils, Sandy clay soils/% | 100.00 | 0.00 | 23.44 | |
Sandy clay loam soils/% | 100.00 | 0.00 | 35.29 | |
Silt loam soils/% | 100.00 | 0.00 | 7.40 | |
Climate | Multi-year average flood-season rainfall/mm | 225.82 | 133.43 | 184.79 |
The sources of the data used in this paper are digital elevated model (DEM) data from the ALOS (Advanced Land Observing Satellite, launched in 2006) satellite phased-array L-band synthetic aperture radar; NDVI from the Resources and Environmental Science Data Center (http://www.resdc.cn/); Monthly rainfall data from the National Tibetan Plateau Data Center (https://data.tpdc.ac.cn/) (Ding & Peng 2020; Peng 2020; Peng et al. 2017, 2018, 2019); Land use and soil underlying surface data from the results of the Shandong Province Flash Flood Hazard Analysis and Assessment Project.
Parameter transplantation methods
Based on the watershed clustering results and the existing hydrological data, the HEC-HMS hydrological model was used to perform the parameter transplantation test, and the specific methods and processes are as follows:
(1) Combining the clustering results and the available hydrological information to determine the reference watersheds.
(2) Parameter transplantation tests are performed using the HEC-HMS hydrological model, in which the SCS-CN is selected for the Loss model, the Snyder Unit Hydrograph for the Transform model, and the Recession for the Baseflow model, without considering channel routing (Zhao et al. 2017; Cheng et al. 2021). The parameters include curve number, initial abstraction, standard lag, peaking coefficient, initial discharge, recession constant, and ratio.
(3) Select typical flood events to determine parameter calibration for the reference watershed. After passing the calibration is qualified, the parameters are transplanted to the ungauged watersheds. Watersheds in the same category can share a set of model parameters. All the parameters are performed by using the direct transplantation method (Wu et al. 2023).
(4) Compare the simulation effect of many floods and evaluate the parameter transplantation effect.
RESULTS AND DISCUSSION
Results
Initial clustering
The index data were first standardized and then fed into the SOM neural network for the initial phase of clustering. Since this stage was followed by a secondary clustering process, complete convergence of the SOM neural network was not necessary; rather, an approximate clustering result was sufficient. The iteration count was set at 300, and the number of neurons was determined using the empirical formula , where N represents the number of samples (Astel et al. 2007). Given that the study involved 545 samples, the neural network was configured with a neuron count of
. The interclass distances and the hit markers are depicted in the Supplementary material, Figure S2. The preliminary clustering results from the SOM neural network were then input into the K-means algorithm to derive the final clustering results.
Determination of the optimal number of clusters
Determination of optimal cluster number (a) Dunn index curve (b) DB index curve.
Determination of optimal cluster number (a) Dunn index curve (b) DB index curve.
Results of final clustering
Discussion
Characterization of different groups of small watersheds
(a) Spatial distribution of topography, rainfall and other indicators for each group of small watersheds; (b) Percentage of area for each land-use type and soil texture type in each group of small watersheds. (Note: (a): 1,2, … , | 12 means group I, group II, … , group XII. (b): ‘high’, ‘relatively high’, ‘medium’, ‘relatively low’, and ‘low’ in the figure represent the distribution of the indicator values of the watershed group among all the small watersheds).
(a) Spatial distribution of topography, rainfall and other indicators for each group of small watersheds; (b) Percentage of area for each land-use type and soil texture type in each group of small watersheds. (Note: (a): 1,2, … , | 12 means group I, group II, … , group XII. (b): ‘high’, ‘relatively high’, ‘medium’, ‘relatively low’, and ‘low’ in the figure represent the distribution of the indicator values of the watershed group among all the small watersheds).
Groups I, VI, and X, located in mountainous areas, are characterized by steep slopes, high elevations, sparse populations, scattered villages, and significant undulation. Group I, in particular, has a smaller number of watersheds. Its land use is primarily forested land, farmland, and grassland, with the largest area being forested. These watersheds have a relatively high shape coefficient, dense river networks, high vegetation cover, more tortuous rivers, and predominantly sandy loam and clay loam soils. These soil types suggest average water permeability, indicating a tendency for these areas to accumulate rainfall and pose flooding risks. Group VI shares similarities with Group I but differs in soil content and the percentage area of land-use types. Group X's small watersheds have a lower river network density, are predominantly used for cropland, and their runoff processes are significantly influenced by irrigation.
Groups II, VII, VIII, and XII, which exhibit slight differences in characteristics, are situated in hilly areas characterized by large slopes, high elevations, and relatively high density of river networks. In terms of land use, these areas are predominantly comprised of arable land and zones used for housing and construction, indicating frequent human activities. The soil types in these groups are primarily sandy clay loam and sandy loam, which have weak water permeability. In the event of a disaster, these sub-basins are particularly vulnerable to causing significant casualties and economic losses, and therefore, they should be prioritized in flash flood prevention and control efforts.
Groups III, V, and XI are located in regions with gentle slopes and low elevations, which are more densely populated and have a lower likelihood of experiencing flash floods compared to other groups. However, if prolonged and persistent rainfall occurs, the weak water permeability of the soils in these small watersheds could trigger siltation flash floods or waterlogging, potentially leading to economic losses and casualties.
Groups IV and IX are situated near lakes and at the mouths of seas, featuring a small number of small watersheds with topographic characteristics similar to those in Group V. However, the soil types in Groups IV and IX are predominantly powdery loam, which offers good water permeability. The relatively high water area in Group IX suggests a higher number of small reservoirs and other water conservancy facilities within the watersheds of this group, which may influence the hydrological dynamics and management strategies in these areas.
Discussion on reasonableness
The flood peak modulus of the small watersheds can be used to further validate the rationality of the clustering approach. In this study, design peak flows for the 545 small watersheds were calculated using the hydrologic-hydraulic method (Xue 2016).
Flood-modulus results: (a) Flood peak modulus distribution in sub-catchments; (b) Flood peak modulus distribution by sub-basin within six large communities.
Flood-modulus results: (a) Flood peak modulus distribution in sub-catchments; (b) Flood peak modulus distribution by sub-basin within six large communities.
Parameter transplantation test
This paper conducted a parameter transplantation experiment using the HEC-HMS model with two watersheds from Group X (Qinshui River watershed and Dongcun River watershed) in the eastern part of Shandong Province as examples. The Qinshui River watershed was designated as the reference watershed, and the Dongcun River watershed was considered an ungauged watershed. Initially, parameters were determined by selecting 10 floods from the Qinshui River watershed, with the results detailed in Table S2. Subsequently, these parameters were transplanted to the Dongcun River watershed, and the simulation results for five flood events are presented in Table S3 and Figure S5. According to the data in Table S3 and Figure S5, only one of the five floods was deemed unqualified, resulting in a qualification rate of 80%. This indicates that the methodology of grouping basins based on similar physical attributes is a viable approach for hydrological forecasting in ungauged watersheds.
Extension and application of results
This paper further generalizes the outcomes of watershed clustering. Given that the study area encompasses virtually all types of watersheds in Shandong Province, the clustering results have been employed to train a random forest classification model (Jonathan et al. 2022). The trained model can be used to identify the category of any ungauged watershed in Shandong Province.
Taking the Gaocun River Basin as an example, the physical characteristics of the Gaocun River Basin were input into the trained random forest classification model, and the category to which the basin belonged was predicted to be Group I. The physical characteristics of the Gaocun River Basin were then used as the basis for the classification model. Combined with the existing hydrological data, the Baisha River Basin in this group is taken as the reference basin. The parameters determined for the Baisha River Basin were then transplanted to the Gaocun River Basin. The results of the flood simulation are documented in Table S5 and Figure S7. The simulation data indicate that the errors for the four simulated floods are within a reasonable range, and the average NSE is 0.685. This confirms the efficacy of this method for parameter transplantation in any ungauged watersheds within Shandong Province, demonstrating the practical applicability of the model in practical application scenarios.
CONCLUSION
In this research, a comprehensive analysis was conducted on 545 small watersheds in Shandong Province, selecting 22 pertinent indicators related to their runoff production and concentration processes. The clustering of these watersheds was achieved through the application of the SOM neural network and K-means methods. The validity of the classification was confirmed by analyzing the flood peak modules. Meanwhile, parameter transplantation using the HEC-HMS hydrological model was carried out in the ungauged watershed, which proved the feasibility of parameter transplantation based on the similarity method of physical attributes and further verified the rationality of the classification. In addition, the results of the clustering were leveraged to extend the application of these findings within the province. The key findings of this study include:
(a) The 545 watersheds were categorized into 12 similar groups. The diverse and intricate topographical features of Shandong Province, along with the wide selection of similarity indicators, led to a predominantly dispersed spatial distribution within each group, with some degree of clustering observed in a minority of cases. Noticeable disparities were found among different groups, and the clustering results aligned with the geographical distribution of the flood peak modulus. This outcome demonstrates the ability of the method to identify similarities in runoff generation and routing conditions of the watersheds and also helps to clarify the hydrologic characteristics of each watershed.
(b) The watershed clustering method based on physical attribute similarity used in this study proves to be a viable solution for flood forecasting in ungauged watersheds. Looking forward, the applicability and stability of this method can be further analyzed by increasing the number of validated watersheds, provided that sufficient data are available. In addition, assessing the influence of diverse indicators on streamflow production and catchment processes across various watersheds could inform the assignment of weights to these indicators.
ACKNOWLEDGEMENTS
The authors wish to gratefully acknowledge the financial assistance from the Natural Science Foundation of Shandong Province (ZR2020ME249), National Natural Science Foundation of China (42301046), Major Science and Technology Program of the Ministry of Water Resources of the People's Republic of China (SKS-2022013), and the other anonymous reviewer whose comments significantly improved the quality of this paper. The datasets were provided by the National Tibetan Plateau/Third Pole Environment Data Center (http://data.tpdc.ac.cn).
DATA AVAILABILITY STATEMENT
Data cannot be made publicly available; readers should contact the corresponding author for details.
CONFLICT OF INTEREST
The authors declare there is no conflict.