Arts, Faculty of
Non UBC
Geography, Department of
DSpace
Remote Sensing 17 (6): 946 (2025)
Li, Yifan
Zhang, Chendi
Cui, Peng
Hassan, Marwan A.
Duan, Zhongjie
Bhattacharyya, Suman
Yao, Shunyu
Zhao, Yang
The Hengduan Mountains region (HMR) is vulnerable to flash flood disasters, which account for the largest proportion of flood-related fatalities in China. Flash flood regionalization, which divides a region into homogeneous subdivisions based on flash flood-inducing factors, provides insights for the spatial distribution patterns of flash flood risk, especially in ungauged areas. However, existing methods for flash flood regionalization have not fully reflected the spatial topology structure of the inputted geographical data. To address this issue, this study proposed a novel framework combining a state-of-the-art unsupervised Graph Neural Network (GNN) method, Dink-Net, and Shapley Additive exPlanations (SHAP) for flash flood regionalization in the HMR. A comprehensive dataset of flash flood inducing factors was first established, covering geomorphology, climate, meteorology, hydrology, and surface conditions. The performances of two classic machine learning methods (K-means and Self-organizing feature map) and three GNN methods (Deep Graph Infomax (DGI), Deep Modularity Networks (DMoN), and Dilation shrink Network (Dink-Net)) were compared for flash-flood regionalization, and the Dink-Net model outperformed the others. The SHAP model was then applied to quantify the impact of all the inducing factors on the regionalization results by Dink-Net. The newly developed framework captured the spatial interactions of the inducing factors and characterized the spatial distribution patterns of the factors. The unsupervised Dink-Net model allowed the framework to be independent from historical flash flood data, which would facilitate its application in ungauged mountainous areas. The impact analysis highlights the significant positive influence of extreme rainfall on flash floods across the entire HMR. The pronounced positive impact of soil moisture and saturated hydraulic conductivity in the areas with a concentration of historical flash flood events, together with the positive impact of topography (elevation) in the transition zone from the Qinghai–Tibet Plateau to the Sichuan Basin, have also been revealed. The results of this study provide technical support and a scientific basis for flood control and disaster reduction measures in mountain areas according to local inducing conditions.
Academic Editors: Ying Zhang,Peng Dou and Chunlin HuangReceived: 20 January 2025Revised: 19 February 2025Accepted: 5 March 2025Published: 7 March 2025Citation: Li, Y.; Zhang, C.; Cui, P.;Hassan, M.; Duan, Z.; Bhattacharyya,S.; Yao, S.; Zhao, Y. Flash FloodRegionalization for the HengduanMountains Region, China, CombiningGNN and SHAP Methods. RemoteSens. 2025, 17, 946. https://doi.org/10.3390/rs17060946Copyright: © 2025 by the authors.Licensee MDPI, Basel, Switzerland.This article is an open access articledistributed under the terms andconditions of the Creative CommonsAttribution (CC BY) license(https://creativecommons.org/licenses/by/4.0/).ArticleFlash Flood Regionalization for the Hengduan MountainsRegion, China, Combining GNN and SHAP MethodsYifan Li 1,2 , Chendi Zhang 1,* , Peng Cui 1,2,3, Marwan Hassan 4 , Zhongjie Duan 5, Suman Bhattacharyya 4,Shunyu Yao 6 and Yang Zhao 71 Key Laboratory of Land Surface Pattern and Simulation, Institute of Geographic Sciences and NaturalResources Research, Chinese Academy of Sciences, Beijing 100101, China; liyf.20b@igsnrr.ac.cn (Y.L.)2 University of Chinese Academy of Sciences, Beijing 100049, China3 China-Pakistan Joint Research Center on Earth Sciences, CAS-HEC, Islamabad 45320, Pakistan4 Department of Geography, University of British Columbia, Vancouver, BC V6T1Z2, Canada5 Department of Data Science and Engineering, East China Normal University, Shanghai 200062, China6 China Institute of Water Resources and Hydropower Research, Beijing 100038, China7 Sichuan Zipingpu Development Co., Ltd., Chengdu 610091, China* Correspondence: zhangchendi@igsnrr.ac.cnAbstract: The Hengduan Mountains region (HMR) is vulnerable to flash flood disasters,which account for the largest proportion of flood-related fatalities in China. Flash floodregionalization, which divides a region into homogeneous subdivisions based on flashflood-inducing factors, provides insights for the spatial distribution patterns of flash floodrisk, especially in ungauged areas. However, existing methods for flash flood regionaliza-tion have not fully reflected the spatial topology structure of the inputted geographical data.To address this issue, this study proposed a novel framework combining a state-of-the-artunsupervised Graph Neural Network (GNN) method, Dink-Net, and Shapley AdditiveexPlanations (SHAP) for flash flood regionalization in the HMR. A comprehensive datasetof flash flood inducing factors was first established, covering geomorphology, climate,meteorology, hydrology, and surface conditions. The performances of two classic machinelearning methods (K-means and Self-organizing feature map) and three GNN methods(Deep Graph Infomax (DGI), Deep Modularity Networks (DMoN), and Dilation shrinkNetwork (Dink-Net)) were compared for flash-flood regionalization, and the Dink-Netmodel outperformed the others. The SHAP model was then applied to quantify the impactof all the inducing factors on the regionalization results by Dink-Net. The newly developedframework captured the spatial interactions of the inducing factors and characterized thespatial distribution patterns of the factors. The unsupervised Dink-Net model allowedthe framework to be independent from historical flash flood data, which would facilitateits application in ungauged mountainous areas. The impact analysis highlights the sig-nificant positive influence of extreme rainfall on flash floods across the entire HMR. Thepronounced positive impact of soil moisture and saturated hydraulic conductivity in theareas with a concentration of historical flash flood events, together with the positive impactof topography (elevation) in the transition zone from the Qinghai–Tibet Plateau to theSichuan Basin, have also been revealed. The results of this study provide technical supportand a scientific basis for flood control and disaster reduction measures in mountain areasaccording to local inducing conditions.Keywords: flash flood; Graph Neural Network (GNN); regionalization; the HengduanMountains region (HMR); interpretation modelRemote Sens. 2025, 17, 946 https://doi.org/10.3390/rs17060946Remote Sens. 2025, 17, 946 2 of 271. IntroductionFloods are among the most severe and destructive natural disasters worldwide [1,2].Among them, flash floods, the rapid flooding processes in small mountainous watersheds(<200 km2 suggested by the Ministry of Water Resources of China) triggered by heavyrainstorms [3], stand out as particularly devastating and deadly. Between 1949 and 2015,over 60,000 flash floods resulted in approximately 280,000 deaths in China [4,5]. Thenumber of missing and dead caused by flash floods accounted for 66.4% of the totalcasualties owing to flood-related disasters in China during 2000–2021 in average [6]. TheHengduan Mountains region (HMR) in Southwest China is particularly susceptible to flashflood disasters, accounting for 42.3% of the total fatalities in China during 2011–2015 [4,7].In 2014, the deaths by flash floods in the HMR even exceeded the sum of the remainingareas of China [4]. The HMR is also the stronghold for many national key engineeringprojects, e.g., the hydro–wind–solar multi-energy complementation bases in the JinshaRiver and Yalong River [8,9] and the Sichuan–Tibet Railway [10]. The need to address thegreat threats from flash floods for the safety of local people and infrastructure in this regionis great.Understanding the spatial distribution of flash floods is crucial for risk analysis and pre-vention. Traditional methods have advanced by integrating high-resolution ground-basedhydrological data with numerical models that capture hydraulic dynamics and morpho-logical features of river networks [11]. However, their reliance on ground-based networksresults in limited coverage, especially in mountainous regions. In contrast, satellite mon-itoring provides broader and near-real-time data on flood extents, channel morphology,and post-flood topographic changes [12]. Moreover, integrating satellite data with GIS andmachine learning techniques in flash flood regionalization refines our ability to characterizethe spatial distribution patterns [13]. This is particularly valuable for mountainous areas inthe HMR, where monitoring facilities are sparse [14]. Although regionalization analysison climate, vegetation [15], and hydrology [16] in the HMR has been conducted, thesestudies apply traditional indicator thresholds or geographical analysis methods withoutconsidering the spatial relationship between geographical data dimensions and were notdesigned specifically for flash floods in this region.Flash flood regionalization utilizes the inducing factors as input parameters andpartitions the geographical space into homogeneous regions considering both the intra-regional similarity and inter-regional differences of the input parameters [17]. The inducingfactors are closely related to the generation processes of flash floods, which involve complexinteractions between various aspects of geomorphology, climate, meteorology, hydrology,vegetation, and soil conditions [11]. In each homogeneous region, the influence of inducingfactors on flash floods is similar in space, and hence, the risk management strategies andhazard mitigation measures can be consistent for the same region [18]. Additionally, flashflood regionalization facilitates the extrapolation of monitored data from gauged sites toungauged areas within the same region [14,18], which would be of great help to addressthe challenges of limited data access in mountainous areas due to measuring difficulty.Classic machine learning techniques, including K-means, tree-based ensemble mod-els [19], support vector machine (SVM, e.g., [20]), and self-organizing-map-based (e.g., [18]),have been widely employed in flash flood regionalization. These methods perform as thespatial expansion of various clusters obtained by integrating the impacts of the inputinducing factors (e.g., [18]). Zhang et al. [18] applied a two-stage hybrid self-organizingmap-based clustering algorithm to delineate 18 homogeneous flash flood regions accordingto 13 key factors in Jiangxi Province. Deep learning techniques, such as convolutionalneural network (CNN, e.g., [2]) and LSTM neural network (e.g., [1]), have also been ap-plied in flash flood regionalization and showed superior computational efficiency andRemote Sens. 2025, 17, 946 3 of 27feature extraction capabilities. Flash flood regionalization utilizes multi-dimensional ge-ographical data that relates to the formation of flash floods and shows significant spatialcorrelations [18]. The factors considered in regionalization typically include topographiccharacteristics (e.g., slope, elevation [1]), climatic conditions (e.g., precipitation [2]), hydro-logical properties (e.g., drainage density [19]), soil attributes (e.g., type, depth [21]), andland cover types (e.g., vegetation index [2]). Therefore, it is crucial to account for both thegeographical attributes at individual map grids and the spatial connections between thesegrids during the clustering process. Despite their effectiveness, classic machine learningapproaches primarily focus on geographical attributes and often fail to capture the spatialcorrelations crucial for flash flood regionalization [22].In contrast, graph neural networks (GNNs) take both attribute characteristics andspatial topology structures into account in classification tasks in geoscience [23] and hence,show great potential in flash flood regionalization. Recent advancements have enhancedthe performance of GNNs by jointly optimizing embedding learning and graph cluster-ing [24], as well as enhancing the scalability of GNNs for better classification performancein processing large-scale data. Notable examples, such as GraphSAGE [25], S3GC [26], andDink-Net [27], are capable of dealing with large-scale input data at the level of millions.These advancements in GNN would facilitate flash flood regionalization for large areas withmulti-dimensional input attributes [28]. However, since these advanced GNN methodshave only been developed very recently, their performance in flash flood regionalizationstill needs to be tested and evaluated [28].An inherent issue regarding the predictions or classifications by deep learning (DL)models (including GNNs) lies in that these models perform as “black-box” approaches.Fortunately, post-hoc interpretation models are capable of quantifying and visualizing thecomplex interactions and non-linear relationships that DL models have captured from theinput data [29]. Particularly, the Shapley Additive exPlanations (SHAP) model studiesthe sensitivity aspect of DL models by perturbing the original inputs [29,30]. This modelcan not only examine the influence of each input feature on the entire region, parts of theregion, or even each grid but also indicate whether the influence is positive or negative.Therefore, the SHAP model has the promise to quantify the contributions of the inducingfactors to flash flood regionalization results obtained by the GNN methods.The objectives of this study are: (i) to establish an automatic flash flood regionalizationframework by integrating GNNs with interpretation models for the HMR; and (ii) to evalu-ate the framework’s effectiveness in revealing spatial distribution patterns of flash floodsand elucidating the impacts of inducing factors within the HMR. To achieve these researchgoals, we first compiled and analyzed a comprehensive set of attributes related to flashflood generation for the entire HMR. We tested three GNN methods (Deep Graph Infomax(DGI), Deep Modularity Networks (DMoN), and Dilation shrink Network (Dink-Net)) andtwo classical machine learning methods (K-means and the self-organizing feature map(SOFM)) for flash flood regionalization. The SHAP model was then employed to quan-tify the impact of inducing factors on the optimal regionalization result. Last, the spatialdistribution features of flash floods and the impacts of inducing factors were analyzedand discussed.2. Study Area and Data2.1. Study AreaThe area of the Hengduan Mountains region (HMR) is about 430,000 km2, covering5 provinces of China (Figure 1). The terrain elevation in the HMR decreases from northwestto southeast, with altitudes typically exceeding 4000 m in the northwest and falling to lessthan 1500 m in the southeast. This region features seven prominent north–south orientedRemote Sens. 2025, 17, 946 4 of 27mountain ranges (i.e., the Boshulaling Mountain, Taniantaweng Mountain, MangkangMountain, Shaluli Mountain, Daxueshan Mountain, Qionglai Mountain, and Min Moun-tains) and six large rivers (i.e., the Nu (Salween) River, Lancang (Mekong) River, JinshaRiver, Yalong River, Dadu River, and Min River). This topography results in the alternationof high mountains and deep valleys in the east–west direction [31,32].Remote Sens. 2025, 17, x FOR PEER REVIEW 4 of 28 2. Study Area and Data 2.1. Study Area The area of the Hengduan Mountains region (HMR) is about 430,000 km2, covering 5 provinces of China (Figure 1). The terrain elevation in the HMR decreases from northwest to southeast, with altitudes typically exceeding 4000 m in the northwest and falling to less than 1500 m in the southeast. This region features seven prominent north–south oriented ountain ranges (i.e., the Boshulaling ountain, Tanianta eng ountain, angkang ountain, Shaluli Mountain, Daxueshan Mountain, Qionglai Mountain, and Min Mountains) and six large rivers (i.e., the Nu (Salween) River, Lancang (Mekong) River, Jinsha River, Yalong River, Dadu River, and Min River). This topography results in the alternation of high mountains and deep valleys in the east–west direction [31,32]. Figure 1. (a) The digital elevation model of the Hengduan Mountains region (HMR) in China, with the historical flash flood events from 1950 to 2015 marked as black dots; (b) the location of the HMR with the provinces involved marked in grey. The HMR spans the middle subtropical zone and the plateau temperate zone, with an average annual temperature and precipitation of 11.33 °C and 816.20 mm [33]. Influenced by the East Asian monsoon and South Asian monsoon and blocked by north–south mountain ranges, rainfall diminishes from the southeast and southwest towards the northwest hinterland (Figure S5f). The wet seasons, generally from May to October, contribute approximately 90% of the annual precipitation [34]. In summer, the occurrence rate of nocturnal precipitation generally exceeds 50% here [35]. Over the past 50 years, an increase in both the frequency and intensity of rainfall has been observed, which in turn i re 1. (a) e i ital ele ati el f t e e a tai s re i ( ) i i a, it.The R spans the middle subtropical zone and the plateau temperate zone, with anaverage annual temperature and precipitation of 11.33 ◦C and 816.20 mm [33]. Influencedby the East Asian monsoon and South Asian monsoon and blocked by north–south moun-tain ranges, rainfall diminishes from the southeast and southwest towards the northwesthinterland (Figure S5f). The wet seasons, generally from May to October, contribute approx-imately 90% of the annual precipitation [34]. In summer, the occurrence rate of nocturnalprecipitation generally exceeds 50% here [35]. Over the past 50 years, an increase in both thefrequency and intensity of rainfall has been observed, which in turn has elevated flash floodrisks in vulnerable areas [36]. Both the frequency and intensity of extreme rainfall eventsare higher in the southern and western fringe areas of the HMR but lower in the northernpart (Figure S5d,e) [37]. As for soil characteristics, brown soil, red soil, and alpine meadowand steppe soil are the dominant soil types and sandy loam soil is the main soil texturein the region [38]. The highly spatial heterogeneity of soil characteristics and landscapesresults in high diversity of plant communities, with distinct zonation in both vertical andRemote Sens. 2025, 17, 946 5 of 27horizontal gradients [32]. The most widely distributed vegetation types in the HMR areevergreen coniferous forests and shrub grass [32,39]. Because of the abundant rainfall andsteep terrains in the mountainous areas, flood-related disasters are widely common in theHMR (Figure 1a). In general, historical flash flood events appear more frequently in thesouth of the HMR than in the north [38] and are mainly concentrated along the valleys ofthe large rivers flowing in the HMR (Figure 1a).2.2. DataTo comprehensively capture the factors influencing flood formation, an input databasecomprising 36 inducing factors was established (see details in Section 3.1 and Table 1). Thisdataset was selected from an initial set of 43 factors. To reduce redundancy and ensure therobustness of the analysis, 7 factors with correlation coefficients > 0.7 and high varianceinflation were eliminated through Spearman correlation analysis, multicollinearity analysis,and flood mechanism analysis [40]: residual soil moisture, saturated soil moisture, soilbulk density, topsoil sand content, topsoil silt content, topsoil clay content, and streampower index (see details in Supplementary Section S1). The selected input data coveredfive key aspects, i.e., geomorphology, climate, meteorology, hydrology, and underlyingsurface conditions.The original dataset included both raster and vector formats. Runoff data, originally invector form, was converted to raster format, while the remaining data were already in rasterform. Among the 36 variables, geomorphic type, vegetation type, soil type, soil texture,and soil depth were classified data, whereas the others were continuous data. Additionally,factors listed in Table 1 without specified time scales or durations represented averagestates or static indicators. For factors spanning extended time periods, annual averageswere calculated first, followed by multi-year averages (e.g., annual mean rainfall).Extreme rainfall is an important direct trigger for the occurrence of flash floods [18].We have included several key extreme rainfall factors in the input data, and the calculationprocess is as follows. 25 extreme rainfall factors under different statistical frequencies wereinitially calculated based on the 3-h resolution dataset from 1979 to 2018. After testing theclustering results under different data combinations, 12 in Table 1 were selected. To quantifyextreme rainfall factors, the maximum rainfall values at the scales of 3 h, 6 h, 12 h, and24 h in each year were calculated and ranked in descending order. The empirical frequencyof the maximum values for each temporal scale was computed by fitting Pearson III-typeempirical frequency curves with the method of moments [28]. From this, the maximumrainfall values for each temporal scale with frequencies of 1% and 50% were then extracted,representing rainfall with a recurrence interval of 100 years and 2 years, respectively. Thesevalues were denoted as (P3h)1%, (P3h)50%, (P6h)1%, (P6h)50%, (P12h)1%, (P12h)50%, (P24h)1%,and (P24h)50%. Additionally, the annual mean of the maximum values for each temporalscale were also calculated and added to the input dataset as (P3h)mean, (P6h)mean, (P12h)mean,and (P24h)mean. Since flash floods due to snow- or ice-melting only occur in very limitedareas in the HMR, data on snow or ice were not considered in this work [41].The point-based records of 1680 historical flash flood events cover the period during1950–2015 in the HMR (Figure 1), as documented by the National Flash Flood Investiga-tion and Evaluation Project of China [6]. This dataset primarily consists of events causedby heavy rainfall. Flash floods caused by snowmelt or ice melt were not prominentlyreported during this period. The data were used in the validation of the regionalizationresults in Section 3.3, and the analysis was shown in Section 4.2. Each event was georefer-enced with latitude and longitude coordinates, allowing for spatial analysis of flash floodoccurrence patterns.Remote Sens. 2025, 17, 946 6 of 27Table 1. Influencing factors for flash floods in the HMR.Aspect Factors Abbreviation SpatialResolutionTemporalResolution orSpanGeographicalMeanings Data SourceGeomorphologyElevation Dem 30 × 30 m2Averagestates orstaticindicatorHeight above sea.Geospatial data cloud(www.gscloud.cn(accessed on 28 February2025))Slope Slope 1 × 1 km2The ratio of elevationincrement to horizontalincrement.Calculated fromDEM dataAspect Aspect 1 × 1 km2 The topographic slopeorientation.Elevationdifference ED 1 × 1 km2The difference betweenthe maximum andminimum of elevation.Topographicwetness index(TWI)TWI 1 × 1 km2TWI = In[sink flow perunit area/Tan(slope)]The influence oftopography on runoffdirection andaccumulation.Geomorphictype GT —Classification based onelevation andtopographic relief.Resource andEnvironment Science andData Center(https://www.resdc.cn/(accessed on 28 February2025))ClimateAnnual meantemperature Temp 0.1◦ × 0.1◦ 2000~2023(yearly data)The annual averagevalue of temperature.European Centre forMedium- Range WeatherForecasts (ECMWF),ERA5-Land Dataset (https://www.ecmwf.int/(accessed on 28 February2025))Annual meanrainfall Rainfall 1 × 1 km2 2000~2022(yearly data)The annual averagevalue of rainfall.National Tibetan PlateauData Center(http://data.tpdc.ac.cn/(accessed on 28 February2025))Annual meanpotential evapo-transpirationPET 1 × 1 km2 2000~2022(yearly data)The annual averagevalue of potentialevapotranspiration.Terrestrialactual evapo-transpirationAET 1 × 1 km2 2001~2019(yearly data)The annual averagevalue of actualevapotranspiration.National Tibetan PlateauData Center(http://data.tpdc.ac.cn/(accessed on 28 February2025)), ETMonitor GlobalActual EvapotranspirationDataset with 1-kmResolution [42]MeteorologyMaximum of3-h/6-h/12-h/24-h withfrequency of 1%and 50%rainfall, annualaveragemaximumrainfall for3 h/6-h/12-h/24-h(P3h)1%,(P3h)50%,(P3h)mean,(P6h)1%,(P6h)50%,(P6h)mean,(P12h)1%,(P12h)50%,(P12h)mean,(P24h)1%,(P24h)50%,(P24h)mean0.1◦ × 0.1◦ 1979~2018Flash floods areprimarily influenced byshort-termmeteorological indiceswithin a day. Therefore,extreme rainfall dataacross four time scaleswithin 24 h wereselected. The 1% and50% percentilescorrespond to eventfrequencies occurringonce every 100 yearsand once every twoyears, respectively.A Big Earth Data Platformfor Three Poles(https://poles.tpdc.ac.cn/(accessed on 28 February2025)), Chinameteorological forcingdataset [43].Remote Sens. 2025, 17, 946 7 of 27Table 1. Cont.Aspect Factors Abbreviation SpatialResolutionTemporalResolution orSpanGeographicalMeanings Data SourceHydrologyRunoff Runoff 1 × 1 km2 1979–2013(yearly data)Index of the annualdaily runoff for 35 yearsfor each river segment,reflecting both thestructure of the rivernetwork and the flux inthe river network,obtained from VICsimulations.National Tibetan PlateauData Center(http://data.tpdc.ac.cn/(accessed on 28 February2025)), GlobalReconstruction ofNaturalized RiverDischarge at 2.94 MillionRiver Reaches (GRADES)[44,45].Vegetationtranspiration Ec 500 × 500 m2 2000–2020(yearly data)Water lost in plants tothe atmosphere as watervapor.PML-V2(China):evapotranspiration andgross primary productiondataset. A Big Earth DataPlatform for Three Poles(https://poles.tpdc.ac.cn/(accessed on 28 February2025)) [46,47].Vaporization ofinterceptedrainfallEi 500 × 500 m2 2000–2020(yearly data)Water from vegetationcanopy interceptingrainfall evaporates intothe atmosphere.Soilevaporation Es 500 × 500 m2 2000–2020(yearly data)Water in soil thatevaporates into theatmosphere.UnderlyingsurfaceconditionVegetation type Vege 1 × 1 km2Averagestates orstaticindicatorVegetation type.Resource andEnvironment Science andData Center(https://www.resdc.cn/(accessed on 28 February2025))NDVI NDVI 1 × 1 km2 2001~2020(yearly data)Normalized differencevegetation index.Resource andEnvironment Science andData Center(https://www.resdc.cn/(accessed on 28 February2025))LAI LAI 1 × 1 km2 2001~2023(yearly data) Leaf area index.Google Earth Engine(https://developers.google.cn/earth-engine/datasets/catalog/MODIS_061_MOD15A2H(accessed on 28 February2025)), MOD15A2H LAISoil type Stype 1 × 1 km2Averagestates orstaticindicatorSoil type.HWSD v1.1(http://www.ncdc.ac.cn(accessed on 28 February2025))Soil texture Stexture 1 × 1 km2Assemblage of mineralparticles of differentsizes in soil.Soil depth Sdepth 1 × 1 km2 The depth between thesurface and the bedrock.Soil gravelcontent Sgravel 1 × 1 km2 Proportion of gravelvolume in surface soil.Soil saturatedhydraulicconductivityKs 1 × 1 km2In water-saturated soil,the amount of waterpassing through the soilper unit area in a unittime under a unit waterpotential gradient.A High-Resolution GlobalMap of Soil HydraulicProperties(https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/UI5LCE(accessed on 28 February2025)) [48].Remote Sens. 2025, 17, 946 8 of 27Table 1. Cont.Aspect Factors Abbreviation SpatialResolutionTemporalResolution orSpanGeographicalMeanings Data SourceUnderlyingsurfaceconditionSoil moisture Smoisture 1 × 1 km2 2000~2020(yearly data)The amount of watercontained in the soil.A Big Earth Data Platformfor Three Poles(https://poles.tpdc.ac.cn/(accessed on 28 February2025)), China SoilMoisture Dataset [49,50].Soil erodibility Uslek 250 × 250 m2Averagestates orstaticindicatorSoil erodibility reflectsthe sensitivity of soil toerosion and thetransport capacity ofsoil particles.National Tibetan PlateauData Center(http://data.tpdc.ac.cn/(accessed on 28 February2025)), Soil ErodibilityDataset of Pan-Third Pole20 Countries [51,52].3. MethodsThe framework for flash flood regionalization is presented in Figure 2. The frameworkinitiated with collecting and preprocessing the input data on inducing factors of flashfloods. Clustering analysis was then conducted using five methods and different clusternumbers. By postprocessing and mapping of the flash flood regionalization, the optimalmethod and regionalization map for the HMR were determined. Last, the interpretation ofthe regionalization result was achieved based on the SHAP model and the impact of eachinput factor on the regionalization result was revealed.Remote Sens. 2025, 17, x FOR PEER REVIEW 8 of 28 in a unit time under a unit water potential gradient. u/dataset.xhtml?persistentId=doi:10.7910/DVN/UI5LCE (accessed on 28 February 2025)) [48]. Soil moisture Smoisture 1 × 1 km2 2000~2020 (yearly data) The amount of water contained in the soil. A Big Earth Data Platform for Three Poles(https://poles.tpdc.ac.cn/ (accessed on 28 February 2025)), China Soil Moisture Dataset [49,50]. Soil erodibility Uslek 250 × 250 m2 Average states or static indicator Soil erodibility reflects the sensitivity of soil to erosion and the transport capa ity of soil particles. National Tibetan Plateau Data Center (http://data.tpdc.ac. n/ (accessed on 28 February 2025)), Soil Erodibility Dataset of Pan-Third Pole 20 Countries [51,52]. 3. Methods The framework flash flood region lization is pr se ted in Figure 2. The framework initiated with collecting and preprocess g the input data on indu ing factors of flash flood . Clustering analysis was then con ucted using five metho s and differen cluster numbers. By p stprocessing and mapping of the flash flood regionalization, he optimal method and regi alizati n map for the HMR w re determined. Last, the interp tation of the r gionalization result was achieved based on the SHAP model and the impact of each input factor on the regionalization result was revealed. Figure 2. The framework of flash flood regionalization.Remote Sens. 2025, 17, 946 9 of 273.1. Data PreprocessingThe raw database containing 43 inducing factors for flash floods were re-rasterized to2 × 2 km2 in ArcGIS (v10.6). The correlation between the input data will affect the subse-quent influence analysis of the factors and 7 factors were eliminated through correlationanalysis, multicollinearity analysis, and flood mechanism analysis. Consequently, the finalinput data contained 36 inducing factors. Spearman correlation analysis and multicollinear-ity analysis were performed with the pandas library in Python 3.8. The Stream Power Index,soil saturated water content, soil residual water content, soil bulk density, sand content, siltcontent, and clay content were removed and 36 attributes were used for clustering analysis(Table 1 and Section 2.2).For the GNN methods, the attributes and locations of all the grids needed to betransformed into a graph structure G = (V, E, X) first [22]. Here, V = {v1, . . . , vN} representsa set of nodes, where each node is associated with corresponding d-dimension attributes,and E ={eij}denotes a set of edges connecting nodes. The adjacency matrix of G isdenoted as A ∈ {0, 1}N×N , where Aij = 1 if (vi, yi) ∈ E, else Aij = 0. The node attributematrix of G is denoted as X = {x1, . . . xN} ∈ RN×d, where the i-th row of xi ∈ Rd representsthe d-dimensional attributes vector of node i [22].The GNN methods aims to partition the nodes in graph G into k disjoint partitionsor subgraphs {G1, G2G3, . . . , Gk} to ensure that the nodes within the same cluster exhibitsimilarity or proximity in terms of both topology structure and attribute characteristics [27].To convert the grid data into graph structure data, each grid was initially regardedas a node in the graph and the node set V (V = 108503) was constructed. The distancesbetween the nodes in V were calculated, and the edges connecting the nodes were thenestablished within the 8 neighboring nodes according to the given nearest node numberand the upper limit of the node-to-node distances. Nodes represent uniformly distributedgrids in our study, so these edges were undirected and unweighted and the connectionsbetween the nodes did not have any specific direction. This setting allowed informationto flow in both directions in each connection between nodes so that all connections wereequal in importance or influence in graph G [53].The 36 filtered factors were input as the node attributes, which formed the nodeattribute matrix X (X ∈ R108503×36). The adjacency matrix A was constructed to representthe edges connecting nodes in V, where A ∈ R108503×108503.3.2. Clustering Analysis3.2.1. Classic Machine Learning MethodsWe employed two classic machine learning methods for flash flood regionalization: K-means and self-organizing feature map (SOFM). The K-means method repeatedly updatesthe cluster centers to segment all the data points into K clusters by minimizing the sumof squared distances between each data point and its closest cluster center [54]. SOFMefficiently maps high-dimensional data onto a lower-dimensional space based on the ideaof competitive learning and self-organization [55]. The grid data of all the 36 inducingfactors was transformed into a text data format for these two classic machine learningmethods, and both methods were implemented by the scikit-learn library and minisomlibrary in Python 3.8.3.2.2. GNN MethodsGraph Neural Networks (GNNs) employ self-supervised manners to embed the inputdata into the latent space and separate the embeddings into several disjoint clusters [53].Three recently developed GNN methods with the capacity to process large-scale input data,Remote Sens. 2025, 17, 946 10 of 27i.e., Deep Graph Infomax (DGI), Deep Modularity Networks (DMoN), and Dilation shrinkNetwork (Dink-Net), were chosen and implemented by Python 3.8.DGI utilizes the principle of contrastive learning and trains the Graph ConvolutionalNetwork by optimizing local mutual information between patch node representationsand corresponding high-level summaries of graphs. This method employs a binary cross-entropy loss function, leading to an embedding representation that encapsulates the globalinformation content of the entire graph [56]. The DGI method is a two-step process:(i) generating the embedding representation of the graph data, and (ii) applying K-meansto produce the final clustering results.DMoN and Dink-Net are one-step approaches that integrate the clustering processwithin GNNs to yield the ultimate clustering results. Inspired by the modularity qualityfunction and spectral optimization, DMoN applies activation functions to derive soft clusterassignments and utilizes a fully differentiable unsupervised clustering objective to optimizesoft cluster assignments, with a null model to manage graph inhomogeneities [57].Dink-Net first learns representations in a self-supervised manner by discriminatingnodes. The K-means method was then used to initialize the clustering centers of thenode representations, and the clustering centers were assigned as learnable parametersof GNN. Motivated by the dilation and shrink of galaxies in the universe, Dink-Netoptimizes the final clustering distribution by minimizing both cluster dilation loss andcluster shrink loss in an adversarial manner [27]. Unlike DGI and DMoN, Dink-Netunifies representation learning and clustering optimization within a cohesive end-to-endframework. This integration facilitates the extraction of clustering-friendly features andimproves the overall clustering performance.3.3. Evaluation of Clustering PerformanceThe clustering performance by the classic machine learning and GNN methods wasevaluated by two clustering validity indices and the optimal regionalization was identified.First, both the clustering quality index (CQI) and Davies–Bouldin Index (DBI) werecalculated (Equations (1) and (2)) across all methods and all cluster numbers to quantifythe overall performance of each regionalization result. CQI is the sum of the average coeffi-cient of variation of the inducing factors within the clustering results and quantifies thedispersion level of these characteristics [58]. A lower CQI value indicates that the variationin inducing factors within each cluster is smaller, meaning that the flash flood conditioningfactors exhibit higher internal consistency within clusters. This suggests that the clusteringmethod effectively differentiates regions based on similar flash flood conditioning envi-ronments. DBI assesses the compactness within the clusters and separation between theclusters by the distance between cluster centers [59]. A lower DBI value signifies that theaverage distance of all grid points to their cluster centroid is smaller, while the distancebetween different cluster centroids is larger. This indicates that the clusters are spatiallycompact internally and well-separated from each other. Therefore, lower values of bothCQI and DBI indicate better clustering performance. The optimal regionalization methodand result were achieved by simultaneously considering the performance in CQI and DBIof each regionalization result.CQI =1NK∑i=1niP∑j=1CV ji (1)Remote Sens. 2025, 17, 946 11 of 27where N is the number of total grids; K is the cluster number; ni is the number of gridsfor cluster i; P is the number of inducing factors; and CV ji is the coefficient of variation offactor j in cluster i.DBI =1KK∑i,j=1maxi ̸=j{di + djdij}(2)where di is the average distance of all grids in the i-th cluster to the centroid of the cluster;and dij is the distance between the centroid of cluster i and the centroid of cluster j.The point data of historical flash flood events was used to assess the reasonability ofthe optimal regionalization result and how well the regionalization captured the spatialdistribution patterns of flash flood occurrences First, the quantity and density of historicalflash flood events in each subdivision of the optimal regionalization result were calcu-lated. Then, we applied the hotspot analysis to examine the statistically significant spatialstructure of the historical flash flood events [60]. Hotspot analysis calculates the Z-score,which represents the heterogeneity degree of geographical element (i.e., historical flashflood events in this study) distributions relative to the uniform distribution. The Z-scorecan be used to identify the spatial structure features of the geographical elements in themap [60]. A positive Z-score suggests spatial clustering of the geographical elements, whilea negative Z-score suggests a sparser distribution. A larger magnitude of Z-score meansthe distribution is increasingly concentrated or sparse for positive and negative Z-scores,respectively [18]. In this work, the distribution of Z-score was calculated for the entire HMRand compared with the regionalization result (see details in Supplementary Section S2).If some subdivisions from the regionalization overlapped with the densely or sparselydistributed areas denoted by the hotspot analysis, it means that the regionalization couldwell capture the spatial distribution structure of historical flash flood events.3.4. Data Postprocessing and MappingTo enhance the regional continuity and integrity of the results, we refined the resultsthrough a post-processing step. Polygons smaller than 400 km2 were merged with adjacentlarger polygons with the longest common boundaries using the Eliminate tool in ArcGIS(v10.6). The polished regionalization result in raster format was then converted into shapefileformat to generate the regionalization map for each regionalization result.3.5. Interpretation of Regionalization ResultsThe Shapley Additive Explanations (SHAP) interpretation model was used to quantifythe contribution of 36 inducing factors to the overall regionalization result and to each sub-division (Figure 2). The SHAP model is based on the Shapley value concept in cooperativegame theory and calculates the marginal contribution of the input factors to the modeloutput by perturbing the original inputs for the tested model [29,30].Given that the SHAP model can only be applied to supervised models, a RandomForest (RF) model with 100 decision trees was trained using the optimal regionalizationresult obtained by Dink-Net with 12 clusters and the 36 input inducing factors, and theSHAP model was then applied to the trained RF model. The accuracy of the RF modelreached 0.9, indicating that its results well represented the regionalization result by theoriginal GNN model. The SHAP values of all the input factors were calculated for eachgrid in the regionalization result as the local interpretation result to quantify the influenceof each inducing factor at grid scale. The average SHAP absolute values of all the gridsin a subdivision and the entire HMR for each input factor were calculated as the globalinterpretation results. The SHAP and RF models were implemented by the scikit-learnlibrary and shap library in Python 3.8, respectively.Remote Sens. 2025, 17, 946 12 of 274. Results4.1. Clustering Performance of the Five MethodsFigure 3 presents the variations in the DBI and CQI values with the cluster number(K) from 2 to 20 for all the five tested clustering methods. DBI values generally increasewith K for each tested clustering method (Figure 3a). The DBI values of Dink-Net and DGImethods increase for K = 2–10. Then, they stabilize around 2.8 at K = 11–20. The K-meansmethod shows a continuous increase in DBI values in the tested range of K. For the DMoNmethod, the DBI values fluctuate around 4.5 at K = 3–14. And the DBI value remains around4 at K > 14. The DBI values for the SOFM method keep increasing at K ≤ 12 and fluctuatearound 7 at higher cluster numbers. The SOFM has DBI values significantly higher than theother four methods after K reaches 7. Overall, the Dink-Net and DGI methods outperformthe other methods as lower DBI values indicate higher clustering quality.Remote Sens. 2025, 17, x FOR PEER REVIEW 13 of 28 Figure 3. Variations of (a) DBI and (b) CQI with the cluster number (K) for K-means, SOFM, DGI, DMoN, and Dink-Net. Figure 4 presents the regionalization maps for the HMR with all five tested methods with 12 clusters. The regionalization maps by both K-means and SOFM methods generally show much more fragmented subdivisions than those generated by GNN methods. In the regionalization map by K-means, the fragmented subdivisions are mainly located in the western and northeastern parts of the HMR. These areas correspond to the Three Parallel Rivers region (Nu River, Lancang River, and Jinsha River, Figures 4 and 1a) and the transition area from the western Sichuan Plateau to Chengdu Plain, respectively, both characterized by the alternation of deep valleys and high mountains. Similarly, in the regionalization map by SOFM, the subdivisions are seriously broken along the areas with intense river incision, especially the Jinsha River and Yalong River in the central HMR. This indicates that both K-means and SOFM methods are unsuitable for flash flood regionalization in areas with deep incised rivers. The DMoN method generates flawed clustering results due to a significant imbalance in the number of grids across classes in the initial clustering output. In contrast, the regionalization maps generated by the Dink-Net and DGI methods feature fewer subdivisions, well-defined and continuous boundaries, and a better reflection of the spatial characteristics of geographical elements. For instance, subdivisions in the northwestern part of the regionalization maps by these two methods can reflect the drainage structure of the Nu, Lancang, and Jinsha Rivers (Figure 4c,d). Additionally, the Dink-Net method requires less computation time compared to the other methods. Overall, Dink-Net achieves the lowest DBI and CQI values while demonstrating the highest Figure 3. Variations of (a) DBI and (b) CQI with the cluster number (K) for K-means, SOFM, DGI,DMoN, and Dink-Net.In terms of CQI values, Figure 3b shows an overall decreasing trend as the clusternumber inc ases for all five methods. The CQI values of all five m thods show an veralldecreasing trend with the increase of the cluster number (Figure 3b). The CQI of the fivemethods show no clear variation trends at K < 6. The CQI values of the Dink-Net methodRemote Sens. 2025, 17, 946 13 of 27are the lowest of the five methods at K ≥ 6. Although the CQI values of DGI are alsokept at a low level, the fluctuation is much more significant than that for Dink-Net. Insome cases, the CQI values of DGI even become the highest among the five methods at thepeaks (e.g., K = 13). After K reaches 6, the CQI values of DMoN generally exceed those ofK-means, which is higher than those of SOFM. Overall, the values of both CQI and DBI forthe Dink-Net method remain the lowest among all the tested methods, and this methodshows higher potential as the optimal method for regionalization.Figure 4 presents the regionalization maps for the HMR with all five tested methodswith 12 clusters. The regionalization maps by both K-means and SOFM methods generallyshow much more fragmented subdivisions than those generated by GNN methods. Inthe regionalization map by K-means, the fragmented subdivisions are mainly locatedin the western and northeastern parts of the HMR. These areas correspond to the ThreeParallel Rivers region (Nu River, Lancang River, and Jinsha River, Figures 4 and 1a) andthe transition area from the western Sichuan Plateau to Chengdu Plain, respectively, bothcharacterized by the alternation of deep valleys and high mountains. Similarly, in theregionalization map by SOFM, the subdivisions are seriously broken along the areaswith intense river incision, especially the Jinsha River and Yalong River in the centralHMR. This indicates that both K-means and SOFM methods are unsuitable for flash floodregionalization in areas with deep incised rivers.The DMoN method generates flawed clustering results due to a significant imbalancein the number of grids across classes in the initial clustering output. In contrast, the region-alization maps generated by the Dink-Net and DGI methods feature fewer subdivisions,well-defined and continuous boundaries, and a better reflection of the spatial characteris-tics of geographical elements. For instance, subdivisions in the northwestern part of theregionalization maps by these two methods can reflect the drainage structure of the Nu,Lancang, and Jinsha Rivers (Figure 4c,d). Additionally, the Dink-Net method requires lesscomputation time compared to the other methods. Overall, Dink-Net achieves the lowestDBI and CQI values while demonstrating the highest computational efficiency, making itthe optimal clustering method among those tested in this study.We then narrowed the range of cluster numbers to 10–16 for the regionalizationresults by Dink-Net based on the two following considerations. First, lower DBI and CQIvalues are preferred, as they refer to higher clustering quality. Second, a smaller clusternumber facilitates the visualization and interpretation of the regionalization results. Theseregionalization results at K = 10–16 were evaluated against the geographic characteristics ofthe HMR and consequently, the cluster number of 12 was selected to generate the optimalregionalization map (see Supplementary S3 for details). The optimal regionalization resultreflected the spatial distribution patterns of inducing factors for flash floods, and thesubdivisions matched with the heterogeneous spatial distributions of these factors (seeSupplementary S4 for details). Through the above analysis, the optimal clustering method,number of clusters and regionalization result are determined.Remote Sens. 2025, 17, 946 14 of 27Remote Sens. 2025, 17, x FOR PEER REVIEW 14 of 28 computational efficiency, making it the optimal clustering method among those tested in this study. Figure 4. The regionalization maps by (a) K-means, (b) SOFM, (c) DGI, (d) DmoN, and (e) Dink-Netwith 12 clusters. In the clustering result obtained by DMoN method with 12 clusters, the number ofgrids among the clusters is extremely unbalanced. Therefore, there are only 7 clusters in (d).Remote Sens. 2025, 17, 946 15 of 274.2. Regionalization Result vs. Historic Flash Flood EventsThe regionalization result was compared with the spatial distribution of historicalflash flood events (see details in Supplementary Section S2), with the quantity, density,and Z-score of the historical flash flood events for the subdivisions shown in Figure 5.Regions SW-1, SW-7, and SE-8 show relatively high quantity and density of the historicalflash flood events (Figure 5a,b). The larger positive mean values of the Z-score in thesethree subdivisions indicates that these areas are spatially clustered areas for historical flashflood events (Figure 5c). In contrast, the regions NW-3, NW-4, M-6, and M-12 show arelatively low quantity and density of the historical events. Correspondingly, the meanvalues of the Z-score in these subdivisions are negative, indicating that these subdivisionsare sparse areas of historical flash flood events. On balance, these findings suggest thatthe regionalization subdivisions are capable of characterizing the densely and sparselydistributed areas of the historical flash flood events.Remote Sens. 2025, 17, x FOR PEER REVIEW 15 of 28 Figure 4. The regionalization maps by (a) K-means, (b) SOFM, (c) DGI, (d) DmoN, and (e) Dink-Net with 12 clusters. In the clustering result obtained by DMoN method with 12 clusters, the number of grids among the clusters is extremely unbalanced. Therefore, there are only 7 clusters in (d). We then narrowed the range of cluster numbers to 10–16 for the regionalization results by Dink-Net based on the two following considerations. First, lower DBI and CQI values are preferred, as they refer to higher clustering quality. Second, a smaller cluster number facilitates the visualization and interpretation of the regionalization results. These regionalization results at K = 10–16 were evaluated against the geographic characteristics of the HMR and consequently, the cluster number of 12 was selected to generate the optimal regionalization map (see Supplementary S3 for details). The optimal regionalization result reflected the spatial distribution patterns of inducing factors for flash floods, and the subdivisions matched with the heterogeneous spatial distributions of these factors (see Supplementary S4 for details). Through the above analysis, the optimal clustering method, number of clusters and regionalization result are determined. 4.2. Regionalization Result vs. Historic Flash Flood Events The regionalization result was compared with the spatial distribution of historical flash flood events (see details in Supplementary Section S2), with the quantity, density, and Z-score of the historical flash flood events for the subdivisions shown in Figure 5. Regions SW-1, SW-7, and SE-8 show relatively high quantity and density of the historical flash flood events (Figure 5 ,b). The larg r positive mean values of the Z-score in these three subdivisions indicates that these areas are spatially clustered areas for historical flash flood events (Figure 5c). In contrast, the regions NW-3, NW-4, M-6, and M-12 show a relatively low quantity and density of the historical events. Correspondingly, the mean values of the Z-score in these subdivisions are negative, indicating that these subdivisions are sp rse areas of historical flash flood events. On bal nce, thes finding suggest that the regionalization subdivisions are capable of characterizing the densely and sparsely distributed areas of the historical flash flood events. Figure 5. (a) The regionalization map by Dink-Net with 12 subdivisions; (b) quantity and density, and (c) Z-score of the historical flash flood events in each subdivision of the regionalization map. The locations of historical flash flood events are marked with black dots in (a). The event Figure 5. (a) The regionalization map by Dink-Net with 12 subdivisions; (b) quantity and density,and (c) Z-score of the historical flash flood events in each subdivision of the regionalization map. Thelocations of historical flash flood events are marked with black dots in (a). The event number/eventdensity has been marked for several subdivisions in panel (b), and the mean Z-score values havebeen presented in panel (c).To further analyze the characteristics of these regions, we examined their geographicaland environmental attributes. Regions SW-1 and SW-7 are located in the southwest of theHMR, encompassing the lower basins of the Lancang River and Jinsha River (Figure 6a,b).Flash flood events are relatively evenly distributed in Region SW-1 but concentrated in thelower basins of the Lancang River in Region SW-7. Region SE-8, situated in the Pan-xi area(southwest of Sichuan Province), experiences widespread flash flood events throughout theregion (Figure 6a,d). Region SW-1 is relatively flat, while Regions SW-7 and SE-8 exhibitslight undulations. Additionally, the soil moisture conditions vary among these regions:Region SW-1 has a relatively low soil moisture, whereas Region SE-8 has a higher soilmoisture content (Figure S5g). Despite these differences, all three regions share a commoncharacteristic of low soil saturated hydraulic conductivity (Figure S5h).Remote Sens. 2025, 17, 946 16 of 27Remote Sens. 2025, 17, x FOR PEER REVIEW 17 of 28 Figure 6. (a) The locations of the Regions SW-1, SW-7, SE-8, NW-3, NW-4, M-6, and M-12 in the HMR; (b,c): topography and vegetation for Regions SW-1 and SW-7, with low relief and sub-high latitude mountain and mixed coniferous broad-leaved forest; (d,e): topography and vegetation for Region SE-8, with low relief and sub-high latitude mountain and the vegetation of shrub or farmland; (f,g) show the deep dry, hot valley canyon with sparse shrubs in Regions NW-3 and NW-4; (h,i) refer to the wide river valley on Western Sichuan plateau and the vegetation of alpine meadows and shrubland in Regions M-6 and M-12. All photos were taken by Yifan Li. 4.3. Influence of Inducing Factors on Flash Flood Regionalization The impact of inducing factors for the entire HMR based on the global interpretation of the SHAP model is presented in Figure 7. Notably, the five most influential factors on flash flood regionalization are temperature, (P12h)mean, (P24h)mean, soil moisture, and (P24h)50%, and their impacts occupy proportions of 8.29%, 7.41%, 6.49%, 5.71%, and 5.69% in the sum of 36 factors’ absolute SHAP values. Importantly, six of the top ten factors are related to the extreme rainfall and account for 55.03% of the impact from all the inducing factors, highlighting the significant influence of extreme rainfall in driving flash floods within the HMR. Furthermore, Figure 7 shows the influence of the inducing factors and their variability across the 12 subdivisions. For the flash flood-prone areas (e.g., Regions SW-1, SW-7, and SE-8), soil moisture, temperature, and extreme rainfall-related factors (e.g., (P12h)mean, (P24h)mean, and (P24h)50%) emerge as the most significant effects in general. Additionally, soil saturated hydraulic conductivity (Ks) plays an important role in Region SE-8, ranking as the second most important factor for this subdivision. In contrast, flash flood sparse areas, such as Regions NW-3, NW-4, M-6, and M-12, exhibit a different set of dominant factors. Here, temperature, potential evapotranspiration, and (P12h)mean are the Figure 6. (a) The locati s of the Regions SW-1, SW-7, SE-8, NW-3, -4, M-6, and -12 in the HMR;(b,c): topography and vegetation for Regi ns SW-1 and SW-7, with low relief and sub-high latitudemountain and mixed coniferous broad-leaved forest; (d,e): topography and vegetation for RegionSE-8, with low relief and sub-high latitude mountain and the vegetation of shrub or farmland;(f,g) show the deep dry, hot valley canyon with sparse shrubs in Regions NW-3 and NW-4;(h,i) refer to the wide river valley on Western Sichuan plateau and the vegetation of alpine meadowsand shrubland in Regions M-6 and M-12. All photos were taken by Yifan Li.As the sparsely distributed areas of the historical flash flood events, Regions NW-3and NW-4 are located in the northwest of the HMR and the southeast of Tibet, characterizedby high elevation and significant topographic relief (Figure 6a,f), and Regions M-6 andM-12 are located in the plateau mountains of western Sichuan Province (Figure 6a,h).Regions NW-3 and NW-4 differ significantly in vegetation and geomorphology types fromRegions M-6 and M-12. The former r gions are characterized by dry-hot valleys with strongdowncutting and lateral constraints, teep topography, and exposed surface soils. Sparsshrub is the dominating veget tion. The rare flash flood events in these regions may berelated to the low rainfall (Figures S5d–f). In contrast, Regions M-6 and M-12 are located inthe broad, flat river valleys on the plateau, where the river valleys are generally ‘U’ shapedwith gentler slopes (Figure 6h,i). Unlike the steep terrain of Regions NW-3 and NW-4, theseregions experience weaker incision dynamics and benefit from better vegetation coverage,primarily consisting of alpine meadows and shrublands. Consequently, the improvedvegetation conditions and reduced geomorphic potential energy conditions might lead torare occurrence of flash flood events even though the rainfall in Regions M-6 and M-12 ishigher than that in Regions NW-3 and NW-4.Remote Sens. 2025, 17, 946 17 of 274.3. Influence of Inducing Factors on Flash Flood RegionalizationThe impact of inducing factors for the entire HMR based on the global interpretationof the SHAP model is presented in Figure 7. Notably, the five most influential factorson flash flood regionalization are temperature, (P12h)mean, (P24h)mean, soil moisture, and(P24h)50%, and their impacts occupy proportions of 8.29%, 7.41%, 6.49%, 5.71%, and 5.69%in the sum of 36 factors’ absolute SHAP values. Importantly, six of the top ten factors arerelated to the extreme rainfall and account for 55.03% of the impact from all the inducingfactors, highlighting the significant influence of extreme rainfall in driving flash floodswithin the HMR.Remote Sens. 2025, 17, x FOR PEER REVIEW 18 of 28 three most influential factors. Moreover, in Region NW-3, the geomorphic type is found to be particularly critical, ranking the second among all inducing factors. Figure 7. The average of SHAP absolute values (average impact on model output magnitude) of inducing factors on the regionalization result from the SHAP model for the entire Hengduan Mountains region (HMR). To further illustrate the local impact of these factors, taking the three clustered areas for flash floods (i.e., Regions SW-1, SW-7, and SE-8) as examples, Figure 8 illustrates the local impact of the inducing factors and the way of their influence. In both Regions SW-1 and SW-7, soil moisture and temperature are the top two influential factors, with reversed rankings in these two regions. These are followed by extreme rainfall-related factors. Notably, the way in which soil moisture and temperature influence flash floods are the opposite (Figure 8a,b): higher soil moisture negatively impacts flash flood occurrence, reducing the likelihood of an event, while higher temperatures show a positive impact, demonstrating that elevated temperatures promote flash flood occurrence. In Region SE-8, soil moisture is the most significantly influencing factor, followed by soil saturated hydraulic conductivity (Figure 8c). The higher the soil moisture and lower soil saturated hydraulic conductivity (Ks) would lead to the stronger positive impact effect on flash floods. As for the extreme rainfall-related factors, higher values generally have a greater positive influence on the occurrence of flash floods in all three regions. Specifically, the annual mean values of the maximum rainfall at 6 h-, 12 h-, and 24 h-scales showed more significantly positive influence than the remaining rainfall-related factors in the three regions except that (P24h)50% had a higher positive effect in Regions SW-7 and SE-8 (Figure 8b,c). Figure 7. The average of SHAP absolute values (average impact on model output magnitude)of inducing factors on the regionalization result from the SHAP model for the entire HengduanMountains region (HMR).Furthermore, Figure 7 shows the influence of the inducing factors and their variabilityacross the 12 subdivisions. For the flash flood-prone areas (e.g., Regions SW-1, SW-7,and SE-8), soil moisture, temperature, and extreme rainfall-related factors (e.g., (P12h)mean,(P24h)mean, and (P24h)50%) emerge as the most significant effects in general. Additionally,soil saturated hydraulic conductivity (Ks) plays an important role in Region SE-8, rankingas th second ost import nt factor for this subdivision. In c ntrast, flash flood sparseareas, such as Regions NW-3, NW-4, M-6, and M-12, exhibit a different set of dominantfactors. Here, temperature, potential evapotranspiration, and (P12h)mean are the threemost influential factors. Moreover, in Region NW-3, the geomorphic type is found to beparticularly critical, ranking the second among all inducing factors.Remote Sens. 2025, 17, 946 18 of 27To further illustrate the local impact of these factors, taking the three clustered areasfor flash floods (i.e., Regions SW-1, SW-7, and SE-8) as examples, Figure 8 illustrates thelocal impact of the inducing factors and the way of their influence. In both Regions SW-1and SW-7, soil moisture and temperature are the top two influential factors, with reversedrankings in these two regions. These are followed by extreme rainfall-related factors.Notably, the way in which soil moisture and temperature influence flash floods are theopposite (Figure 8a,b): higher soil moisture negatively impacts flash flood occurrence,reducing the likelihood of an event, while higher temperatures show a positive impact,demonstrating that elevated temperatures promote flash flood occurrence. In RegionSE-8, soil moisture is the most significantly influencing factor, followed by soil saturatedhydraulic conductivity (Figure 8c). The higher the soil moisture and lower soil saturatedhydraulic conductivity (Ks) would lead to the stronger positive impact effect on flash floods.As for the extreme rainfall-related factors, higher values generally have a greater positiveinfluence on the occurrence of flash floods in all three regions. Specifically, the annual meanvalues of the maximum rainfall at 6 h-, 12 h-, and 24 h-scales showed more significantlypositive influence than the remaining rainfall-related factors in the three regions except that(P24h)50% had a higher positive effect in Regions SW-7 and SE-8 (Figure 8b,c).Remote Sens. 2025, 17, x FOR PEER REVIEW 19 of 28 Figure 8. The SHAP value of the top 20 important inducing factors for each 2 × 2 km2 grid in (a) Region SW-1, (b) Region SW-7, and (c) Region SE-8. Each dot in a panel represents the data for a grid. The inducing factors in each panel are ranked in descending order according to the factors’ local impact in each subdivision. The dot color indicates the attribute values of the corresponding factors, with red referring to the high attribute values of a factor. The x-axis value represents the SHAP values that quantify the factors’ impact on the classification tendency for that grid. Figure 9 illustrates the spatial distributions of the subdivision-scale global impacts from several representative inducing factors. Overall, the impacts of all of these inducing factors on flash flood regionalization exhibit significant spatial heterogeneity at the subdivision scale in general. Temperature shows a pronounced positive impact along the river valleys in the HMR, particularly in Region W-2, which encompasses the Three Parallel Rivers region (Figure 9a). The two key extreme rainfall factors, (P12h)mean and (P24h)mean, exhibit significantly positive contributions in the south and east of the HMR (i.e., Regions E-11, SW-1, and SE-8), with the strongest impact emerging along the eastern edge of the HMR (Figure 9b,c). Elevation exhibits a highly positive contribution in the central and northeast parts of the HMR, generally along the Yalong, Dadu, and Min Rivers in the HMR, which locate in the transition zone between the Qinghai-Tibet Plateau and the Sichuan Basin (Figure 9d). Soil moisture contributes positively in the southern part of HMR, where the historical flash flood events are densely distributed (Figures 6a and 9e). Lastly, soil saturated hydraulic conductivity exhibits a pronounced influence mainly in the eastern and western parts of HMR, while showing less impacts in the northwestern areas (Figure 9f). Figure 8. The SHAP value of the top 20 important inducing factors for each 2 × 2 km2 grid in(a) Region SW-1, (b) Region SW-7, and (c) Region SE-8. Each dot in a panel represents the data fora grid. The inducing factors in each panel are ranked in descending order according to the factors’local impact in each subdivision. The dot color indicates the attribute values of the correspondingfactors, with red referring to the high attribute values of a factor. The x-axis value represents theSHAP values that quantify the factors’ impact on the classification tendency for that grid.Figure 9 illustrates the spatial distributions of the subdivision-scale global impactsfrom several representative inducing factors. Overall, the impacts of all of these inducingfactors on flash flood regionalization exhibit significant spatial heterogeneity at the subdi-vision scale in general. Temperature shows a pronounced positive impact along the rivervalleys in the HMR, particularly in Region W-2, which encompasses the Three ParallelRivers region (Figure 9a). The two key extreme rainfall factors, (P12h)mean and (P24h)mean,exhibit significantly positive contributions in the south and east of the HMR (i.e., RegionsE-11, SW-1, and SE-8), with the strongest impact emerging along the eastern edge of thRemote Sens. 2025, 17, 946 19 of 27HMR (Figure 9b,c). Elevation exhibits a highly positive contribution in the central andnortheast parts of the HMR, generally along the Yalong, Dadu, and Min Rivers in the HMR,which locate in the transition zone between the Qinghai-Tibet Plateau and the SichuanBasin (Figure 9d). Soil moisture contributes positively in the southern part of HMR, wherethe historical flash flood events are densely distributed (Figures 6a and 9e). Lastly, soilsaturated hydraulic conductivity exhibits a pronounced influence mainly in the eastern andwestern parts of HMR, while showing less impacts in the northwestern areas (Figure 9f).Remote Sens. 2025, 17, x FOR PEER REVIEW 20 of 28 Figure 9. The SHAP value spatial distribution of inducing factors for each subdivision: (a) temperature; (b) mean of maximum 12 h rainfall ((P12h)mean); (c) mean of maximum 24 h rainfall ((P24h)mean); (d) elevation (Dem); (e) soil moisture (Smoisture); and (f) soil saturated hydraulic conductivity (Ks). 5. Discussion 5.1. Strengths of the New Regionalization Framework The newly developed framework that combines GNN and SHAP models for flash flood regionalization effectively applied the spatial distribution information within the input data and provided a quantitative assessment of the influence of inducing factors. The framework has shown the following advantages over the existing methods for flash flood regionalization through the comparisons among the tested methods and with historical flash flood events. First, the GNN methods accounted for both the attribute characteristics and spatial topology structure of the input geographic data in regionalization, which resulted in better performance over the classic methods. This had been reflected in the lower values of both DBI and CQI for the Dink-Net and DGI methods than the K-means and SOFM methods across the tested cluster numbers (see details in Section 3.3 and Figure 3). In contrast, the other GNN method, the DMoN, showed higher or similar DBI values compared to the K-means (Figure 3a) and the highest CQI values among all tested methods at K > 8 (Figure 3b). The DMoN method produced clustering results in which some clusters had very few grids, leading to significant discrepancies with other clusters. The subdivisions corresponding to these smaller clusters were further merged into adjacent larger subdivisions during postprocessing, resulting in an imbalanced Figure 9. The SHAP value spatial distribution of inducing factors for each subdivision: (a) tempera-ture; (b) mean of maximum 12 h rainfall ((P12h)mean); (c) mean of maximum 24 h rainfall ((P24h)mean);(d) elevation (Dem); (e) soil moisture (Smoisture); and (f) soil saturated hydraulic conductivity (Ks).5. Discussion5.1. Strengths of the New Regionalization FrameworkThe newly developed framework that combines GNN and SHAP models for flashflood regionalization effectively applied the spatial distribution information within theinput data and provided a quantitative assessment of the influence of inducing factors. Theframework has shown the following advantages over the existing methods for flash floodregionalization through the comparisons among the tested methods and with historicalflash flood events.Remote Sens. 2025, 17, 946 20 of 27First, the GNN methods accounted for both the attribute characteristics and spatialtopology structure of the input geographic data in regionalization, which resulted in betterperformance over the classic methods. This had been reflected in the lower values of bothDBI and CQI for the Dink-Net and DGI methods than the K-means and SOFM methodsacross the tested cluster numbers (see details in Section 3.3 and Figure 3). In contrast, theother GNN method, the DMoN, showed higher or similar DBI values compared to the K-means (Figure 3a) and the highest CQI values among all tested methods at K > 8 (Figure 3b).The DMoN method produced clustering results in which some clusters had very few grids,leading to significant discrepancies with other clusters. The subdivisions correspondingto these smaller clusters were further merged into adjacent larger subdivisions duringpostprocessing, resulting in an imbalanced regionalization result (Figure 4e). These resultsindicate that not all GNN methods are suitable for flash flood regionalization, and theirsuitability should be assessed before application.Compared to previous studies that applied K-means, hierarchical clustering [14], tree-based ensemble models [20], or self-organizing-map-based models [18], the proposed GNNapproach effectively captures spatial dependencies and preserves geographical coherencein clustered regions. For example, in the regionalization results obtained by both K-means and SOFM methods, highly fragmented subdivisions were observed in the westernregion of the HMR (Figure 4a,b). This region corresponded to the Three Parallel Riversregion, where alternation of high mountains and deep valleys lead to abrupt changes intopography, landforms, and related meteorological and climatic conditions [38]. In contrast,the regionalization results obtained by Dink-Net and DGI methods well characterized thedrainage structure and topographic differences in this region (Figure 4c,e). These resultssuggest that the two traditional clustering methods faced challenges in handling the regionswith significant and successive topographic relief.Second, the Dink-Net model in the framework was unsupervised, meaning thatthe historical flash flood event records were not needed in the regionalization process,in contrast to supervised approaches [2,21,61]. This makes it easier to implement ourframework in flash flood-prone regions, where access to the historical data is normallynecessary to train the supervised models but is also limited, especially in mountainousregions [7]. Even though the historical flash flood event records exist sometimes, thequality and coverage are difficult to guarantee, and hence, the sparse data points in certainareas does not necessarily indicate a lower flash flood risk. Such historical data wouldintroduce biases into the supervised model through the training process. Nevertheless,without applying the historical records in training, the regionalization results obtainedby the unsupervised Dink-Net model will not be affected by the biases in the historicaldata. For example, in Region SW-7, the eastern part had fewer historical flash flood eventsthan the west (Figure 5a). If historical records were used together with a supervised neuralnetwork, the subdivision might be divided into two or more subdivisions, although theinducing factors in this area showed the similar features.Third, our framework quantified the impact of all the input data on regionalization atvarious spatial scales. It not only revealed the general impact of all the input factors on theentire study area, but also quantified the impact on each subdivision and even each grid(Figures 7–9). Owing to this multi-scale quantification, much more details of the spatialheterogeneity of the impact from the input inducing factors had been captured. For instance,our results visualized the different impacts of soil moisture across various subdivisions(Figures 8 and 9). In Regions SW-1 and SE-8, soil moisture showed the most significantinfluence at the subdivision scale but its influence operated in opposite directions at thegrid scale. In Region SW-1, lower soil moisture levels corresponded to a stronger positiveeffect, whereas in Region SE-8, higher soil moisture levels exerted a greater positive impact.Remote Sens. 2025, 17, 946 21 of 27In comparison, the probabilistic model (e.g., [62,63] or the information gain ratio method(e.g., [1]) in previous studies were only able to illustrate the effects of the inducing factorsfor the whole study area, without further information at subdivision or even grid scale.5.2. Impact of Inducing Factors on Flash Floods in the HMRBased on the interpretation of the regionalization results through the SHAP model,the influence of the inducing factors was quantified and their spatial distribution wasilluminated. These quantitative results provided insights into the formation mechanismsof flash floods in the HMR. The four most influential factors for flash floods in the HMRwere identified and discussed as follows, i.e., temperature, extreme rainfall, soil moisture,and elevation.Our results showed that the extreme rainfall-related factors had the most significantinfluence on the regionalization (Figure 7). The average cumulative influence of extremerainfall-related factors across all the subdivisions reached 46.64% (Figure 7). In all thethree clustered flash flood regions (i.e., Regions SW-1, SW-7, and SE-8), the cumulativeimpact of these factors exceeded 50% of the entire impact from all factors (Figure 8). Thesignificant influence of extreme rainfall-related factors in this study is consistent with thefinding in Tan et al. [64] that extreme rainfall-related factors contributed 42.24–58.42% offlash flood occurrences. The intense and frequent extreme rainfalls in the HMR mainlyresulted from the strong impacts of the monsoons, e.g., East Asian monsoon and SouthAsian monsoon [38]. These monsoonal systems bring substantial moisture, resultingin extreme rainfall predominantly concentrated in the southern, eastern, and westernfringe regions of the HMR. This spatial pattern exhibits high heterogeneity, as shownin Figure S5d–f. The pronounced spatial variability of extreme rainfall-related factorsconsequently amplifies their influence in the clustering process, increasing their relativecontribution to the classification.The role of topography, such as elevations and slopes, was also identified in thisstudy. The local elevations may affect the formation of orographic precipitation, andsteep slopes might facilitate the rapid rise in streamflow [65]. The significant influenceof topographic factors was particularly pronounced in the regions with strong valleydowncutting, steep bank slopes, and active geological structures, such as Regions M-5,NE-10, and M-12 (Figures 7 and 9). These regions are mainly located in the transition zonebetween the Qinghai–Tibet Plateau and the Sichuan Basin [38]. However, in areas withlower elevations and relief (e.g., Regions SW-1, SW-7, and SE-8), topography played asecondary role compared to soil characteristics (Figures 6–8). These results suggest thata large topographic relief affects flash floods in a more significantly way than elevation.Previous studies have also reported the important effect of topography on the formationof flash floods, contributing up to 35% to the general flash flood risk [13,66]. Our resultsshowed consistency with previous research but further quantified the different roles ofvarious topographic factors.Our results revealed the significant effects of temperature and evapotranspirationon flash floods, which had not been reported frequently in previous studies (e.g., [19,21]).Temperature had the most significant impact on flash floods across the entire HMR, andit was also the top influential factor for Region SW-7 with densely distributed historicalflash flood events (Figures 7 and 8b). For flash flood sparse Regions NW-3 and M-12, thepotential evapotranspiration showed the greatest impact (Figure 7) because the attributevalue was significantly lower than that of the surrounding regions (Figure S5i). Bothtemperature and potential evapotranspiration serve as the indicators of heat and influencerunoff generation through complex interactions with soil moisture and erosion, vegetationgrowth, and bedrock weathering [67,68]. Additionally, temperature is also closely relatedRemote Sens. 2025, 17, 946 22 of 27to local elevation, rainfall, and potential evapotranspiration. The broad interactions withother inducing factors might help explain the overwhelming effect of temperature tosome degree.The impact of soil moisture ranked fourth among all the input factors to the en-tire HMR and showed significant spatial heterogeneity, in terms of both magnitude andinfluence direction, in the clustered flash flood areas in the southern part of the HMR(Figure 7, Figure 8c, Figure 9e and Figure S5g). For example, in term of influence magni-tude, the influence of soil saturated hydraulic conductivity is prominent in Region SE-8 butmuch less prominent in Regions SW-1 and SW-7. The influence direction of soil moisture isthe same in Regions SW-1 and SW-7 but opposite to that in Region SE-8 (Figures 8 and 9e,f).In Region SE-8, the influence of soil moisture and saturated hydraulic conductivity indicatesthat soil conditions with high water content and poor water conductivity are favorableto flash flood formation. Under such humid and poorly permeable soil conditions, trans-forming to infiltration would be limited during rainfall [69], making it more likely togenerate surface runoff. Thus, we infer that excess infiltration runoff is the dominant runoffgeneration mechanism in this area. This aligns with findings by Liu et al. [70] that thesurface runoff was primarily controlled by saturated excess infiltration mechanism undersaturated soil conditions in the purple soil in this region. In contrast, in Regions SW-1 andSW-7, although dry soil has a high initial infiltration capacity under certain conditions,long-term low water content will gradually make the soil hydrophobic, thereby reducingthe actual infiltration rate under extreme rainfall conditions, resulting in more rainfall beingconverted into surface runoff [11]. In addition, dry soils tend to be associated with lowervegetation cover and organic matter content, which also weakens the soil’s water retentioncapacity [38], resulting in the rapid surface runoff and increasing the risk of flash floods.This explains the favorable conditions for higher occurrence of flash flood under lower soilmoisture conditions in Regions SW-1 and SW-7 (Figure 8a,b) and highlights the crucial roleof water transport characteristics of surface soil in flash flood generation.5.3. Limitations and Future WorkAlthough the presented framework shows satisfying performance in flash flood re-gionalization in the HMR and effectively quantifies the impact of the inducing factors, thisstudy also has several limitations, as discussed below.First, the input dataset might not cover all the factors that affect the physical formationprocesses of flash floods (e.g., soil porosity, maximum consecutive days of rainfall, andaccumulated temperature in different seasons), despite our efforts to make the dataset ascomprehensive as possible. Because of the limited sources, low resolution, incompletecoverage, and low accuracy of the remote sensing data in the mountainous areas [12], somedata dimensions cannot fit the spatial resolution requirement of this work. In future, theexpansion of the input data for the flash flood regionalization can be expected with theimprovement of the hydrometeorological monitoring system in the HMR.Second, the selection of the inducing factors is inevitably subjective when buildingthe input dataset. When the dimension number of the input data for a given aspect islarge, the weight of this aspect within the entire input dataset will increase and mightaffect its influence on flash flood regionalization. For instance, there are 12 factors relatedto extreme rainfall in the input data, and this might intensify the contribution of extremerainfall in the regionalization result, although the proportion of extreme rainfall-relateddata in our dataset has already been lower than in previous research (e.g., [18]). A morecomprehensive input dataset as mentioned in the first limitation might address this issuein the future. Additionally, an automatic method to determine the optimal combination ofRemote Sens. 2025, 17, 946 23 of 27input factors [21] will help to reduce the subjectivity in data selection, which will be theobjective of our next research work.Third, a trained Random Forest model was built specially for the SHAP model tointerpret the regionalization results as SHAP model was not designed for unsupervisedneural networks [30]. Despite this limitation of the SHAP model, we still chose SHAP for itsoutstanding ability to quantify the impact of input factors at various spatial scales [29,71].Should new interpretation models be invented specifically for unsupervised GNNs in thefuture, we will consider replacing the SHAP model.Last, we applied the framework only in the HMR to address the significant threatof flash floods in this area. Our framework may also work in other mountainous areasprone to flash floods as long as the dataset for the inducing factors is available. Furtherapplications of our framework in other areas are required to assess its generality and tomake necessary modifications in flash flood-prone areas in the following research.6. ConclusionsThis study presented a newly developed framework for flash flood regionalizationthat combined Dink-Net Graph Neural Network (GNN) and SHAP model for a moreaccurate and interpretable analysis of flash floods in the Hengduan Mountains Region.The comparisons between the traditional and GNN clustering methods on the regional-ization performance, together with the visualization of impacts from various inducingfactors, showed great promise in our new framework. First, our framework considersboth the attribute characteristics and spatial topology structures of the input data and theregionalization result effectively captured the spatial distribution patterns of the inducingfactors compared to previous methods. Second, since the GNN model is unsupervised,the regionalization would not be limited by the access to historical flash flood event data,which is normally scarce in mountainous areas. Last, our framework quantifies the influ-ence of inducing factors on flash floods and visualizes the spatial distributions of theseimpacts in the HMR. The interpretation results of the new framework reveals the dominantimpact of extreme rainfall-related factors and their integrated impact, which accounted for46.64% of all the tested inducing factors in the HMR. The influence of elevation on flashfloods becomes significant mainly in the transition zone between the Qinghai–Tibet Plateauand the Sichuan Basin. In the clustered flash flood areas, the impact of soil moisture andconductivity water conduction were particularly prominent and spatially heterogeneous.In conclusion, this study advances flash flood regionalization by integrating cutting-edgedeep learning techniques with explainable AI, thereby enhancing both the accuracy andinterpretability of the regionalization results. The proposed framework not only provides apowerful tool to reveal the spatial distribution properties of flash flood risk but also offerscritical insights to the occurrence mechanism of flash floods. These findings provide ascientific basis for improving flood control and disaster reduction strategies in mountainousregions. Our framework has shown great potential to be applied in early warning estab-lishment, guiding land-use plan and infrastructure development. The distributions andinterpretations of key inducing factors can also provide policy-makers with key informationto efficiently launch mitigation efforts.Supplementary Materials: The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/rs17060946/s1. The Supplementary Materials mainly containsinformation on the results of spearman correlation analysis and multicollinearity analysis (Table S1,Figures S1 and S2), the hotspot analysis (Figure S3), determination of the optimal regionalization(Figure S4), and the comparison between optimal regionalization result and key inducing factors(Figure S5).Remote Sens. 2025, 17, 946 24 of 27Author Contributions: Conceptualization, Y.L. and C.Z.; methodology, Y.L. and C.Z.; software, Y.L.,Z.D. and S.Y.; validation, C.Z., P.C., M.H. and S.B.; formal analysis, Y.L. and C.Z.; data curation, Y.L.and C.Z.; writing—original draft preparation, Y.L. and C.Z.; writing—review and editing, all authors;visualization, Y.L.; funding acquisition, C.Z. and P.C; supervision, P.C. and M.H. All authors haveread and agreed to the published version of the manuscript.Funding: This research was funded by the National Key R & D Program of China (2023YFC3006701,2024YFC3012704), the National Natural Science Foundation of China (42471086), the internationalpartnership program of the Chinese Academy of Sciences (177GJHZ2022064FN), and SichuanProvince Zipingpu Development Co., Ltd. (ZPPC2024-05).Data Availability Statement: The sources of the data used in this study can be found in Table 1.Acknowledgments: The authors gratefully thank Yifan Liu from Peking University for his assistancein model construction. Guotao Zhang and Yaqiao Wu from the Institute of Geographic Sciencesand Natural Resources Research, Chinese Academy of Sciences, are kindly acknowledged for theirsuggestions on manuscript preparation. We sincerely thank the four anonymous reviewers for theirvaluable comments, which have greatly improved this paper.Conflicts of Interest: Author Yang Zhao was employed by the company Sichuan Zipingpu Develop-ment Co., Ltd. The remaining authors declare that the research was conducted in the absence of anycommercial or financial relationships that could be construed as a potential conflict of interest.AbbreviationsThe following abbreviations are used in this manuscript:HMR Hengduan Mountains regionGNN Graph Neural NetworkSHAP Shapley Additive exPlanationsSOFM self-organizing feature mapDGI Deep Graph InfomaxDMoN Deep Modularity NetworksDink-Net Dilation shrink NetworkCQI clustering quality indexDBI Davies–Bouldin IndexK cluster numberReferences1. Fang, Z.; Wang, Y.; Peng, L.; Hong, H. Predicting flood susceptibility using LSTM neural networks. J. Hydrol. 2021, 594, 125734.[CrossRef]2. Liu, J.; Wang, J.; Xiong, J.; Cheng, W.; Sun, H.; Yong, Z.; Wang, N. Hybrid Models Incorporating Bivariate Statistics and MachineLearning Methods for Flash Flood Susceptibility Assessment Based on Remote Sensing Datasets. Remote Sens. 2021, 13, 4945.[CrossRef]3. Arabameri, A.; Saha, S.; Chen, W.; Roy, J.; Pradhan, B.; Bui, D.T. Flash flood susceptibility modelling using functional tree andhybrid ensemble techniques. J. Hydrol. 2020, 587, 125007. [CrossRef]4. He, B.; Huang, X.; Ma, M.; Chang, Q.; Tu, Y.; Li, Q.; Zhang, K.; Hong, Y. Analysis of flash flood disaster characteristics in Chinafrom 2011 to 2015. Nat. Hazards 2018, 90, 407–420. [CrossRef]5. Liu, Y.; Yang, Z.; Huang, Y.; Liu, C. Spatiotemporal evolution and driving factors of China’s flash flood disasters since 1949. Sci.China Earth Sci. 2019, 49, 408–420. [CrossRef]6. Ministry of Water Resources of China. Bulletin of Flood and Drought Disasters in China; Ministry of Water Resources of China,Beijing: Beijing, China, 2021.7. Ma, M.; He, B.; Wan, J.; Jia, P.; Guo, X.; Gao, L.; Maguire, L.W.; Hong, Y. Characterizing the Flash Flooding Risks from 2011 to 2016over China. Water 2018, 10, 704. [CrossRef]8. Tian, L.; Ming, B.; Zhang, W.; Huang, K. Multi-time scale complementarity analysis of hydropower, wind power and photoelectricresources in lower reaches of Jinsha River. J. Hydroelectr. Eng. 2023, 42, 40–49. (In Chinese with English abstract)9. Wen, X.; Sun, Y.; Tan, Q.; Lei, X.; Ding, Z.; Liu, Z.; Wang, H. Risk and Benefit Analysis of Hydro-wind-solar Multi-energy SystemConsidering the One-day Ahead Output Forecast Uncertainty. Adv. Eng. Sci. 2020, 52, 32–41. (In Chinese with English abstract)Remote Sens. 2025, 17, 946 25 of 2710. Cui, P.; Ge, Y.-G.; Li, S.; Li, Z.; Xu, X.-W.; Zhou, G.G.D.; Chen, H.-Y.; Wang, H.; Lei, Y.; Zhou, L.; et al. Scientific challenges indisaster risk reduction for the Sichuan–Tibet Railway. Eng. Geol. 2022, 309, 106837. [CrossRef]11. Vichta, T.; Deutscher, J.; Hemr, O.; Tomášová, G.; Žižlavská, N.; Brychtová, M.; Bajer, A.; Shukla, M.K. Combined effects ofrainfall-runoff events and antecedent soil moisture on runoff generation processes in an upland forested headwater area. Hydrol.Process. 2024, 38, e15216. [CrossRef]12. Ding, L.; Ma, L.; Li, L.; Liu, C.; Li, N.; Yang, Z.; Yao, Y.; Lu, H. A Survey of Remote Sensing and Geographic Information SystemApplications for Flash Floods. Remote Sens. 2021, 13, 1818. [CrossRef]13. Costache, R.; Pham, Q.B.; Sharifi, E.; Linh, N.T.T.; Abba, S.I.; Vojtek, M.; Vojteková, J.; Nhi, P.T.T.; Khoi, D.N. Flash-FloodSusceptibility Assessment Using Multi-Criteria Decision Making and Machine Learning Supported by Remote Sensing and GISTechniques. Remote Sens. 2020, 12, 106. [CrossRef]14. Rami, O.; Hasnaoui, M.D.; Ouazar, D.; Bouziane, A. A mixed clustering-based approach for a territorial hydrological regionaliza-tion. Arab. J. Geosci. 2021, 15, 75. [CrossRef]15. Zhang, X.Q.; Xu, X.M.; Li, X.; Cui, P.; Zheng, D. A new scheme of climate-vegetation regionalization in the Hengduan MountainsRegion. Sci. China-Earth Sci. 2024, 67, 751–768. [CrossRef]16. Liu, Y.; Zou, Q.; Lu, Y.; Li, J.; Xiao, P. Eco-hydrological division of the watershed in the rapid topographic change basin on theeastern Tibetan Plateau. Shuili Xuebao 2022, 53, 243–252. (In Chinese with English abstract)17. Ahani, A.; Mousavi Nadoushani, S.S.; Moridi, A. Regionalization of watersheds by finite mixture models. J. Hydrol. 2020, 583,124620. [CrossRef]18. Zhang, R.; Chen, Y.; Zhang, X.; Ma, Q.; Ren, L. Mapping homogeneous regions for flash floods using machine learning: A casestudy in Jiangxi province, China. Int. J. Appl. Earth Obs. Geoinf. 2022, 108, 102717. [CrossRef]19. Ekmekciog˘lu, Ö.; Koc, K.; Özger, M.; Is¸ık, Z. Exploring the additional value of class imbalance distributions on interpretable flashflood susceptibility prediction in the Black Warrior River basin, Alabama, United States. J. Hydrol. 2022, 610, 127877. [CrossRef]20. Costache, R.; Hong, H.; Pham, Q.B. Comparative assessment of the flash-flood potential within small mountain catchmentsusing bivariate statistics and their novel hybrid integration with machine learning models. Sci. Total Environ. 2020, 711, 134514.[CrossRef]21. Hosseini, F.S.; Choubin, B.; Mosavi, A.; Nabipour, N.; Shamshirband, S.; Darabi, H.; Haghighi, A.T. Flash-flood hazard assessmentusing ensembles and Bayesian-based machine learning models: Application of the simulated annealing feature selection method.Sci. Total Environ. 2020, 711, 135161. [CrossRef]22. Wang, S.; Yang, J.; Yao, J.; Bai, Y.; Zhu, W. An Overview of Advanced Deep Graph Node Clustering. IEEE Trans. Comput. Soc. Syst.2023, 11, 1302–1314. [CrossRef]23. Liu, G.; Ouyang, S.; Qin, H.; Liu, S.; Shen, Q.; Qu, Y.; Zheng, Z.; Sun, H.; Zhou, J. Assessing spatial connectivity effects on dailystreamflow forecasting using Bayesian-based graph neural network. Sci. Total Environ. 2023, 855, 158968. [CrossRef] [PubMed]24. Wang, C.; Pan, S.; Hu, R.; Long, G.; Jiang, J.; Zhang, C. Attributed Graph Clustering: A Deep Attentional Embedding Approach.In Proceedings of the International Joint Conference on Artificial Intelligence, Macao, China, 10–16 August 2019. [CrossRef]25. Hamilton, W.; Ying, Z.; Leskovec, J. Inductive representation learning on large graphs. Adv. Neural Inf. Process. Syst. 2017, 30,1025–1035.26. Devvrit, F.; Sinha, A.; Dhillon, I.; Jain, P. S3GC: Scalable self-supervised graph clustering. Adv. Neural Inf. Process. Syst. 2022, 35,3248–3261.27. Liu, Y.; Liang, K.; Xia, J.; Zhou, S.; Yang, X.; Liu, X.; Li, S.Z. Dink-net: Neural clustering on large graphs. In Proceedings of theInternational Conference on Machine Learning, Honolulu, HI, USA, 23–29 July 2023; pp. 21794–21812. [CrossRef]28. Shi, K.; Chen, Y.; Zhang, X.; Ma, Q.; Ren, L. Flash Flood Hazard Regionalization Based on Graph Clustering Neural Network inJiangxi Province,China. Geogr. Geo-Inf. Sci. 2023, 39, 7–15. (In Chinese with English abstract)29. Jiang, S.; Sweet, L.-b.; Blougouras, G.; Brenning, A.; Li, W.; Reichstein, M.; Denzler, J.; Shangguan, W.; Yu, G.; Huang, F.; et al. HowInterpretable Machine Learning Can Benefit Process Understanding in the Geosciences. Earth’s Future 2024, 12, e2024EF004540.[CrossRef]30. Lundberg, S.M.; Lee, S.-I. A unified approach to interpreting model predictions. Adv. Neural Inf. Process. Syst. 2017, 30, 4765–4774.31. Shi, Z.; Deng, W.; Zhang, S. Spatio-temporal pattern changes of land space in Hengduan Mountains during 1990–2015. J. Geogr.Sci. 2018, 28, 529–542. [CrossRef]32. Yao, Y.; Zhang, B.; Han, F.; Pang, Y. Diversity and geographical pattern of altitudinal belts in the Hengduan Mountains in China.J. Mt. Sci. 2010, 7, 123–132. [CrossRef]33. Yang, Z.-I.; Zhang, T.-B.; Yi, G.-H.; Li, J.-J.; Qin, Y.-B.; Chen, Y. Spatio-temporal variation of Fraction of Photosynthetically ActiveRadiation absorbed by vegetation in the Hengduan Mountains, China. J. Mt. Sci. 2021, 18, 891–906. [CrossRef]34. He, Y.; Xiong, Q.; Yu, L.; Yan, W.; Qu, X. Impact of Climate Change on Potential Distribution Patterns of Alpine Vegetation in theHengduan Mountains Region, China. Mt. Res. Dev. 2020, 40, 48–54. [CrossRef]Remote Sens. 2025, 17, 946 26 of 2735. Wang, J.; Jiang, X.; Zhang, X. Characteristics of temporal and spatial variations of nocturnal precipitation in China. J. NanjingUniv. Nat. Sci. 2022, 58, 750–765. (In Chinese with English abstract)36. Chen, Z.; Wang, L.; Li, X.; Xue, Y.; Jia, H. Spatiotemporal Change Characteristics of Extreme Precipitation in South-western Chinaand its Relationship with Intense ENSO Events. Plateau Meteorol. 2022, 41, 604–616. (In Chinese with English abstract) [CrossRef]37. Bian, Y.; Sun, P.; Zhang, Q.; Liu, R.; Ma, Z.; Zou, Y.; Lyu, Y. Spatial distribution characteristics of extreme climatic events in thehengduan mountains Region. Water Resour. Hydropower Eng. 2021, 52, 1–15. [CrossRef]38. Li, Y.; Zhang, C.; Zhang, G. The development characteristics and formation modes of rainstorm-triggered flash flood disasters inthe Hengduan Mountains. Acta Geogr. Sin. 2024, 79, 600–616. (In Chinese with English abstract) [CrossRef]39. Yang, J.; Dai, J.; Yao, H.; Tao, Z.; Zhu, M. Vegetation distribution and vegetation activity changes in the Hengduan Mountainsfrom 1992 to 2020. Acta Geogr. Sin. 2022, 77, 2787–2802. (In Chinese with English abstract) [CrossRef]40. Mia, M.U.; Rahman, M.; Elbeltagi, A.; Abdullah-Al-Mahbub, M.; Sharma, G.; Islam, H.M.T.; Pal, S.C.; Costache, R.; Islam, A.R.M.T.;Islam, M.M.; et al. Sustainable flood risk assessment using deep learning-based algorithms with a blockchain technology. GeocartoInt. 2022, 38, 1–29. [CrossRef]41. Cui, M.; Zhou, G.; Zhang, D.; Zhang, S. Global snowmelt flood disasters and their impact from 1900 to 2020. J. Glaciol. Geocryol.2022, 44, 1898–1911.42. Zheng, C.; Jia, L.; Hu, G. Global land surface evapotranspiration monitoring by ETMonitor model driven by multi-source satelliteearth observations. J. Hydrol. 2022, 613, 128444. [CrossRef]43. He, J.; Yang, K.; Tang, W.; Lu, H.; Qin, J.; Chen, Y.; Li, X. The first high-resolution meteorological forcing dataset for land processstudies over China. Sci. Data 2020, 7, 25. [CrossRef]44. Lin, P.; Pan, M.; Beck, H.E.; Yang, Y.; Yamazaki, D.; Frasson, R.; David, C.H.; Durand, M.; Pavelsky, T.M.; Allen, G.H.; et al. GlobalReconstruction of Naturalized River Flows at 2.94 Million Reaches. Water Resour. Res. 2019, 55, 6499–6516. [CrossRef]45. Lin, P.; Pan, M.; Yang, Y. Global Reconstruction of Naturalized River Discharge at 2.94 Million River Reaches (GRADES); NationalTibetan Plateau Data Center/Third Pole Environment Data Center: Beijing, China, 2022. [CrossRef]46. He, S.; Zhang, Y.; Ma, N.; Tian, J.; Kong, D.; Liu, C. A daily and 500 m coupled evapotranspiration and gross primary productionproduct across China during 2000–2020. Earth Syst. Sci. Data 2022, 14, 5463–5488. [CrossRef]47. He, S.; Zhang, Y. PML-V2 (China): Evapotranspiration and Gross Primary Production Dataset (2000.02.26–2020.12.31); National TibetanPlateau Data Center/Third Pole Environment Data Center: Beijing, China, 2022. [CrossRef]48. Zhang, Y.; Schaap, M.G.; Zha, Y. A High-Resolution Global Map of Soil Hydraulic Properties Produced by a HierarchicalParameterization of a Physically Based Water Retention Model. Water Resour. Res. 2018, 54, 9774–9790. [CrossRef]49. Li, Q.; Shi, G.; Shangguan, W.; Nourani, V.; Li, J.; Li, L.; Huang, F.; Zhang, Y.; Wang, C.; Wang, D.; et al. A 1km daily soil moisturedataset over China using in situ measurement and machine learning. Earth Syst. Sci. Data 2022, 14, 5267–5286. [CrossRef]50. Shangguan, W.; Li, Q.; Shi, G. China Soil Moisture Dataset (2000–2020), A Big Earth Data Platform for Three Poles; National TibetanPlateau Data Center/Third Pole Environment Data Center: Beijing, China, 2022. [CrossRef]51. Yang, M.; Yang, Q.; Zhang, K.; Li, Y.; Wang, C.; Pang, G. Effects of Content of Soil Rock Fragments on Calculating of SoilErodibility. Acta Pedol. Sin. 2021, 58, 1157–1168. (In Chinese with English abstract)52. Yang, Q. Soil Erodibility Dataset of Pan-Third Pole 20 Countries (2020, with a Resolution of 7.5 arc Second); National Tibetan PlateauData Center/Third Pole Environment Data Center: Beijing, China, 2021. [CrossRef]53. Zhou, S.; Xu, H.; Zheng, Z.; Chen, J.; Li, Z.; Bu, J.; Wu, J.; Wang, X.; Zhu, W.; Ester, M. A Comprehensive Survey on DeepClustering: Taxonomy, Challenges, and Future Directions. ACM Comput. Surv. 2024, 57, 69. [CrossRef]54. Hartigan, J.A.; Wong, M.A. Algorithm AS 136: A K-Means Clustering Algorithm. J. R. Stat. Society. Ser. C Appl. Stat. 1979, 28,100–108. [CrossRef]55. Kohonen, T. The self-organizing map. Proc. IEEE 1990, 78, 1464–1480. [CrossRef]56. Velicˇkovic´, P.; Fedus, W.; Hamilton, W.L.; Liò, P.; Bengio, Y.; Hjelm, R.D. Deep Graph Infomax. Int. Conf. Learn. Represent.2019, 2, 4.57. Tsitsulin, A.; Palowitch, J.; Perozzi, B.; Müller, E. Graph clustering with graph neural networks. J. Mach. Learn. Res. 2023, 24, 1–21.[CrossRef]58. Mao, Q.; Peng, J.; Liu, Y.; Wu, W.; Zhao, M.; Wang, Y. An ecological function zoning approach coupling SOFM and SVM: A casestudy in Ordos. Acta Geogr. Sin. 2019, 74, 460–474. (In Chinese with English abstract) [CrossRef]59. Davies, D.L.; Bouldin, D.W. A Cluster Separation Measure. IEEE Trans. Pattern Anal. Mach. Intell. 1979, PAMI-1, 224–227.[CrossRef]60. Hazaymeh, K.; Almagbile, A.; Alomari, A.H. Spatiotemporal Analysis of Traffic Accidents Hotspots Based on GeospatialTechniques. ISPRS Int. J. Geo-Inf. 2022, 11, 260. [CrossRef]61. Yin, Y.; Zhang, X.; Guan, Z.; Chen, Y.; Liu, C.; Yang, T. Flash flood susceptibility mapping based on catchments using an improvedBlending machine learning approach. Hydrol. Res. 2023, 54, 557–579. [CrossRef]Remote Sens. 2025, 17, 946 27 of 2762. Cao, Y.; Jia, H.; Xiong, J.; Cheng, W.; Li, K.; Pang, Q.; Yong, Z. Flash Flood Susceptibility Assessment Based on Geodetector,Certainty Factor, and Logistic Regression Analyses in Fujian Province, China. ISPRS Int. J. Geo-Inf. 2020, 9, 748. [CrossRef]63. Chen, W.; Li, Y.; Xue, W.; Shahabi, H.; Li, S.; Hong, H.; Wang, X.; Bian, H.; Zhang, S.; Pradhan, B.; et al. Modeling floodsusceptibility using data-driven approaches of naïve Bayes tree, alternating decision tree, and random forest methods. Sci. TotalEnviron. 2020, 701, 134979. [CrossRef]64. Tan, X.Z.; Li, Y.; Wu, X.X.; Dai, C.; Zhang, X.L.; Cai, Y.P. Identification of the key driving factors of flash flood based on differentfeature selection techniques coupled with random forest method. J. Hydrol. Reg. Stud. 2024, 51, 101624. [CrossRef]65. Marchi, L.; Borga, M.; Preciso, E.; Gaume, E. Characterisation of selected extreme flash floods in Europe and implications forflood risk management. J. Hydrol. 2010, 394, 118–133. [CrossRef]66. Tariq, A.; Yan, J.; Ghaffar, B.; Qin, S.; Mousa, B.G.; Sharifi, A.; Huq, M.E.; Aslam, M. Flash Flood Susceptibility Assessment andZonation by Integrating Analytic Hierarchy Process and Frequency Ratio Model with Diverse Spatial Data. Water 2022, 14, 3069.[CrossRef]67. Qiu, A.-N.; Zhang, Y.-J.; Wang, G.-X.; Cui, J.; Song, Y.-X.; Sun, X.-Y.; Chen, L. A modified TOPMODEL introducing the bedrocksurface topographic index in Huangbengliu watershed, China. J. Mt. Sci. 2022, 19, 3517–3532. [CrossRef]68. Venegas-Cordero, N.; Cherrat, C.; Kundzewicz, Z.W.; Singh, J.; Piniewski, M. Model-based assessment of flood generationmechanisms over Poland: The roles of precipitation, snowmelt, and soil moisture excess. Sci. Total Environ. 2023, 891, 164626.[CrossRef] [PubMed]69. Ke, L.; Junfang, C.; Ruoxuan, L.; Yaling, Z.; Minxi, L.; Li, G.; Euihua, N. Mechanisms of soil moisture response to rainfallinfiltration in dry-hot valley of southwest China. South-North Water Transf. Water Sci. Technol. 2024, 22, 736–746. (In Chinese withEnglish abstract)70. Liu, G.; Tian, G.; Shu, D.; Lin, S.; Liu, S. Characteristics of surface runoff and throughflow in a purple soil of Southwestern Chinaunder various rainfall events. Hydrol. Process. 2005, 19, 1883–1891. [CrossRef]71. Wang, N.; Zhang, H.; Dahal, A.; Cheng, W.; Zhao, M.; Lombardo, L. On the use of explainable AI for susceptibility modeling:Examining the spatial pattern of SHAP values. Geosci. Front. 2024, 15, 101800. [CrossRef]Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individualauthor(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury topeople or property resulting from any ideas, methods, instructions or products referred to in the content.
Article
eng
Reviewed
Vancouver : University of British Columbia Library
Multidisciplinary Digital Publishing Institute
CC BY 4.0
Faculty
Researcher
flash flood
Graph Neural Network (GNN)
regionalization
the Hengduan Mountains region (HMR)
interpretation model
Flash Flood Regionalization for the Hengduan Mountains Region, China, Combining GNN and SHAP Methods
Text