Geographic information system models with fuzzy logic for susceptibility maps of debris flow using multiple types of parameters: a case study in Pinggu District of Beijing, China

Debris flow is one of the main causes of loss of life and infrastructure damage in mountainous areas. This hazard should be recognized in the early stage of land development planning. According to field investigation and expert experience, a scientific and effective quantitative susceptibility assessment model was established in Pinggu District of Beijing. This model is based on geographic information system (GIS) combined with grey relational, data-driven and fuzzy logic methods. The influence factors, which are divided into two categories and consistent with the system characteristics of a debris flow gully, are selected, but also a new important factor is proposed. The results of the 17 models are verified using data published by the authority and validated by two other indexes, as well as area under curve (AUC). Through the comparison and analysis of the results, we believe that the streamlining of factors and scientific classification should attract attention from other researchers to optimize a model. We also propose a good perspective to make better use of the watershed feature parameters. These parameters fit well with the watershed units. With full use of insufficient data, scientific calculation and reliable results, the final optimal susceptibility map could potentially help decision makers in determining regional-scale land use planning and debris flow hazard mitigation. The model has advantages in economically weak areas with insufficient data in mountainous areas because of its simplicity, interpretability and engineering usefulness.

1 Introduction

Debris flows are processes of rapid transport of water and soil materials in mountain watersheds, with sudden and destructive outbreaks (Di et al., 2019). Some debris flows can often cause devastating disasters and huge losses (Zhang et al., 2021), seriously threaten the lives and properties of people in the mountains and the safety of major projects, and restrict social and economic development (Iverson, 1997; Hungr et al., 2005; Hu et al., 2011; Takahashi, 2014; Wu et al., 2019). Mass movements in Beijing range in scale from shallow slope failures and rockfalls to catastrophic rock avalanches frequently mobilize to form debris flows, threatening the ecological environment of the mountainous area (Zhong et al., 2004). Especially in recent years, due to the superposition of extreme rainstorm weather and human engineering activities, debris flow events have increased gradually (Z. Li et al., 2021). As the capital of China, Beijing also has strong influence at home and abroad where geological disasters are widely concerned (Xie et al., 2004; Li et al., 2020b). With the deepening understanding of debris flow disaster and the updating of a database, a new and more accurate evaluation is also very necessary. Therefore, it is of great significance to establish an accurate and scientific debris flow susceptibility map.

Through previous studies, it can be summarized that the current research on debris flow mainly focuses on the following aspects: study on mechanism of debris flow, study on early warning and prediction of debris flow, study on numerical simulation of debris flow, and study on debris flow hazard analysis. Especially studies on debris flow hazard analysis gain the attention of the researchers as soon as they appear (Dong et al., 2009). Communicating information about debris flow hazard analysis is a crucial component of preparedness and hazard mitigation (Chiou et al., 2015). Susceptibility assessment, an important part of a hazard assessment of geological processes, is more flexible (Y. Li et al., 2021). In the early days, the susceptibility assessment of debris flows was mainly qualitative research using geomorphological information (Guzzetti et al., 1999). In 1976, the United Nations commissioned the International Union of Engineering Geology to conduct a risk assessment of debris flows, which marked the beginning of research on the susceptibility assessment of debris flows as an important research direction for disaster prevention and prediction (Li et al., 2020b). Many methods and techniques have been proposed to evaluate debris flow susceptibility assessment based on different qualitative and quantitative approaches along with geo-environmental information (Liu and Wang, 1995), such as the analytic hierarchy process (Wu et al., 2016), logistic regression method (Regmi et al., 2013; Conoscenti et al., 2015), information value (Akbar and Ha, 2011; Melo et al., 2012), support vector machine (Pourghasemi et al., 2017), frequency ratio (FR) (Sun et al., 2018), certainty factor (CF) (Tsangaratos and Ilia, 2015), neural network (Lee et al., 2003; Liu et al., 2005) and Bayesian network algorithm (Liang et al., 2012; Tien Bui et al., 2012). These methods have corresponding advantages and limitations for research subjects with different geological conditions. Generally speaking, it is easier to get satisfactory results by combining and comparing various methods (Meyer et al., 2014; Di Napoli et al., 2020; Fang et al., 2020). In summary, with the development of mathematical theory, the susceptibility assessment of debris flows has been extensively and quantitatively studied, and the research methods have also changed from single to comprehensive.

The economy in mountainous areas is often limited; we cannot supervise and verify every basin due to limited funds. The debris flow susceptibility assessment can give decision makers a basis for rational allocation of resources and determine which gullies should be focused on. In other words, the study plays a link role for other studies. Recently, with the development of mathematical theory, computer technology and the application of 3S (remote sensing, geography information systems, global positioning systems), the susceptibility assessment of debris flows has been extensively and quantitatively studied (Li et al., 2020a). As research progresses, debris flows are increasingly seen as an open system. There are many factors influencing the system, the combination of factors is nonlinear, and the interactions are chaotic. Therefore, it is very difficult to find a unified and standard evaluation model. At present, when the information is insufficient, field investigation and the experience of experts are necessary. However, the experience is often subjective and needs a lot of professional experience accumulation. It is very important to express the experience of experts objectively and understandably to serve decision makers. The application of fuzzy set theory in geographic information system (GIS) environments is effective for similar problems (Luo and Dimitrakopoulos, 2003; Porwal et al., 2006).

The main objective of this paper is to propose a GIS-based model. The results of expert experience scoring and site surveys are used as guidance and reference in the modeling process. We have tried to apply methods that can indicate the nonlinearity of the debris flow system. Finally, the modeling process should respect the laws of geomorphological evolution and the geological basis. Otherwise, the result will tend to be simply data fitting (Porwal et al., 2006).

2 Study area

The study area is located to the northeast of Beijing, China (Fig. 1), with a total area of 948.24 km². The elevation of Pinggu is high in the northeast and low in the southwest. It is surrounded by mountains and accounts for about two-thirds of the total area on three sides in the southeast and north. The central and southern parts are alluvial plains. The area, geologically, is the west extension of the famous Jixian section, whose bedrock is mainly Middle and Late Proterozoic dolomite (Lü et al., 2017). The administrative unit of Pinggu District is used as the study area boundary, mainly considering that geological hazards frequently influence human economic activities, so political factors must be taken into account. Within the administrative region, inconsistent decision-making can be effectively avoided.

Figure 1 Study area.

3 Data and methodology

In this study, the susceptibility assessment of debris flow hazard was based on the drainage basin unit. In such a model, a hydrological response unit can fully represent the hillside hydrological process and will make the results more meaningful (Khan et al., 2013, 2016; Zou et al., 2019). First, drainage networks were extracted from the ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer) digital elevation model (DEM) by using the ArcGIS ArcHydro Toolbox, and regions without obvious watershed characteristics were directly deleted. Then for each drainage basin, 19 controlling and triggering factors divided into two types were calculated. In addition, since these factors have different characteristics, different methods were used to calculate the fuzzy membership for different type factors. Field investigation is generally required in geological hazard surveys. If these data are applied to the model, it can help with the model building and reduce the time for model training. The weights derived from the grey relational analysis method used in the following section (Sect. 3.4.1) are based on the data from the field investigation. While geology and geomorphology factors are independent of watershed characteristics, it is suitable to use statistical methods to determine the objective weight. Finally, the debris flow susceptibility index (DFSI) map was derived by overlaying the factor thematic layers with the fuzzy logic method. The workflow of debris flow susceptibility assessment is shown in Fig. 2. First, a DEM map of the Pinggu area was downloaded. Then, the basin units were generated from the DEM map using the ArcHydro tool. The derived results were analyzed, and units that did not fit the characteristics of the watershed were removed. During the analysis, the field investigation data and Google images were referenced. After that, the controlling and triggering factors for the remaining 135 catchments were counted. For the fuzzy memberships, watershed characteristic parameters were determined by grey correlation, and the geological and geomorphological factors were determined by the frequency ratio (FR) method and the cosine amplitude method. Finally, the individual layers were overlaid by fuzzy logic operations to obtain the final map. As there were different combinations of factors, 17 results were derived. Three indexes (area under curve, AUC, accuracy ratio, AR, and resolution ratio, RR) were used to evaluate advantages and disadvantages of these results.

Figure 2 Workflow of debris flow susceptibility assessment.

3.1 Debris flow basin division and inventory

There are many geological hazard points in mountainous areas, so it is not realistic to monitor them completely by professional teams. According to the monitoring and prevention staff and the villagers, the detailed field investigation (Fig. 3) for the evidence collection of debris flows will be carried out at the reported disaster point, aiming at recording the loose material, delineating the basin and exploring other important information on the debris flow gullies. Moreover, field investigation is also very important for model modification. Then based on the hydrology module in ArcGIS 10.2, the research object can be determined. Compared with grid unit and slope unit, the hydrological response unit for susceptibility of debris flow has greater advantages (Z. Li et al., 2021b; Zou et al., 2019). Finally, referring to the result of the field investigation and the remote sensing image, 135 basins are divided after removing the flat and irregular areas (Fig. 4), and of them 48 basins were investigated in the field, accounting for 36 %.

Figure 3 Field investigation photos. (a) Loose material, (b) Middle and Late Proterozoic dolomite, (c) colluvium deposit, (d) slope fracture and (e) channel erosion phenomenon.

Figure 4 Debris flow basin division and inventory. Note: the data of debris flow points come from the Beijing Municipal Commission of Planning and Natural Resources websites (Beijing Municipal Commission of Planning and Natural Resources, 2022).

3.2 Debris flow controlling and triggering factors

The basic requirement for the assessment of debris flows is that some factors included are easily obtainable, meaningful for susceptibility assessment, and can be used for evaluating the need for passive or active debris flow mitigation. According to previous studies, 19 factors are selected in this study. The factors are divided into two types (Table 1) because of their different characteristics. Watershed characteristic factors (Type A) can be directly quantified once the basin is determined (Fig. 5). The influence of these parameters is bounded by the watershed; geology and geomorphology factors (Type B) need to be further processed even if the watershed is determined. The scope of these parameters is independent of the watershed boundary.

Table 1 Factors for susceptibility assessment.

Note: the geological maps are provided by the Beijing Institute of Geological and Prospecting Engineering, and the digital elevation model (DEM) of the study area is from the Shuttle Radar Topography Mission (SRTM) DEM with a resolution of 30 m (http://www.gscloud.cn/sources/accessdata/310?pid=302, last access: 25 June 2022).

Figure 5 Graphical illustration of some Type A factors. A₁ is the planimetric (projected) area of the catchment, A₂ is the curved surface area of the catchment, A₄ is the perimeter of the catchment, A₆ is the curve length of the main channel, A₇ is the straight length of the main channel, and A₁₃ is basin volume.

3.3 Fuzzy logic in susceptibility modeling

Fuzzy set theory is proposed by Zadeh (1965). It is an efficient way of expressing the concept of partial set membership degree. This concept differs from classical binary (0–1 value) logic. More words with transitional fuzzy descriptions (such as low, medium and high) are used (Kritikos and Davies, 2015). This fuzzy expression is particularly applicable to geological hazard classification. In the theory of fuzzy sets, elements have different degrees of membership in the interval [0, 1]: 1 represents complete membership, and 0 represents non membership. Ross (1995) showed that fuzzy systems are useful in two general situations (Kritikos and Davies, 2015). The method is very consistent with the characteristics of debris flow system, whose predisposing factors are fuzzy in nature and whose mechanism is complex and not fully understood. In the application of the fuzzy logic method, the critical step is to find the suitable fuzzy membership of factors. Fuzzy membership degree is equivalent to the weight in the expert scoring method, which is calculated by the objective method rather than given subjectively.

3.4 Fuzzy memberships

3.4.1 Grey relational analysis (GRA) in susceptibility modeling

GRA is proposed by Deng (1982), and it is an important part of grey system theory (Wang et al., 2014). Comparing with mathematical statistics methods which need lots of sample data, typical probability distribution and large calculations, GRA is applicable to small sample size with the data whether regular or not. There will be no inconsistency between qualitative analysis and quantitative analysis (Deng, 1988). Moreover, it is to excogitate the leading and potential factors that affect the development of the system and quantitatively describe the development and change trend of the system by studying whether the relative change trend of the grey factor variables with complex relationship is consistent in the process of system development and evolution (Liu et al., 2004). Thus, grey correlation analysis is introduced to quantify the correlation between each factor and the evaluation results according to field investigation and expert experience. First, the procedure of GRA is to translate the performance of every alternative into a comparability sequence (Lin and Lin, 2002; Kuo et al., 2008; Wei et al., 2017). Therefore, according to the technical standard “Specification of geological investigation for debris flow stabilization (DZ/T0220-2006)” (Ministry of Natural Resources of the People’s Republic of China, 2006), published by the China Ministry of Lands and Resources, the preliminary assessment results of debris flow susceptibility are obtained, which are used as the reference sequence of the grey relation method (Table 2). Second, the grey correlation coefficient of all A factors is calculated by Eq. (1). Finally, the average grey relational coefficient (the correlation degree) is calculated by Eq. (2) as the fuzzy memberships (Table 3).

$$\begin{array}{}\text{(1)}& {\mathit{\xi }}_i\left(k\right)={\displaystyle \frac{\begin{array}{c}min\\ i\end{array}\begin{array}{c}min\\ k\end{array}\left|x_{\mathrm{0}}\left(k\right)-x_i\left(k\right)\right|+\mathrm{0.5}\begin{array}{c}max\\ i\end{array}\begin{array}{c}max\\ k\end{array}\left|x_{\mathrm{0}}\left(k\right)-x_i\left(k\right)\right|}{\left|x_{\mathrm{0}}\left(k\right)-x_i\left(k\right)\right|+\mathrm{0.5}\begin{array}{c}min\\ i\end{array}\begin{array}{c}min\\ k\end{array}\left|x_{\mathrm{0}}\left(k\right)-x_i\left(k\right)\right|}},\end{array}$$

where ξ_i(k) is the grey relational coefficient, i= 1, 2, …, n are the number i of type A factors, k= 1, 2, …, n are the number of basins, x₀(k) is the reference sequence (ideal target sequence), and x_i(k) is the number i of the Type A factor sequence.

$$\begin{array}{}\text{(2)}& r_i={\displaystyle \frac{\mathrm{1}}{N}}{\sum }_{i=\mathrm{1}}^n{\mathit{\xi }}_i\left(k\right),\end{array}$$

where r_i is the correlation degree in the range (0, 1). N is the total number of basins in Table 2

Table 2 Quantitative evaluation grade standard table for debris flow susceptibility.

Note: 130 ≥ score ≥ 116, VH; 115 ≥ score ≥ 87, M; 86 ≥ score ≥ 44, L; and 43 ≥ score ≥ 15, N. VH is very high susceptibility, M is moderate susceptibility, L is low susceptibility, and N is non-debris flow.

Table 3 The fuzzy memberships of Type A factors.

3.4.2 Data-driven method in susceptibility modeling

Landslide is one of the main fixed sources of debris flow in mountainous area. Shallow landslides are one of the most common categories of landslides. They frequently involve large areas and different soils in various climatic zones (Benda and Dunne, 1987; Selby, 1982; Borrelli et al., 2014). Great debris flows may result from numerous, small slope failures that subsequently coalesce (Fairchild, 1987; Roeloffs, 1996), from flow enlargement due to incorporation of bed and bank debris (Pierson et al., 1990; Bovis and Dagg, 1992), or from large, individual landslides that mobilize partially or almost totally (Vallance and Scott, 1997; Iverson et al., 1997). Debris flows may also scour steep channels to bedrock and accelerate sediment delivery to downstream, lower-gradient channels. The spatial and temporal distribution of shallow landslides are important controls on landscape evolution and a major component of both natural and management-related disturbance regimes in mountain drainage basins (Tsukamoto et al., 1982; Dietrich et al., 1986; Benda and Dunne, 1987; Crozier et al., 1990). Therefore, the landslide susceptibility assessment methods can be used for reference to debris flow susceptibility assessment.

For Type B factors which cannot be characterized by a specific number, the frequency ratio (FR) method and the cosine amplitude method can be used to derive their fuzzy memberships. The FR ratio defined as Eq. (3). Considering the fuzzy membership must be in the interval [0, 1], the FR values of the different categories are normalized by the largest FR value (Lee, 2006; Pradhan, 2010, 2011a, b) within the same type factor (Table 4) in order to derive the function.

$$\begin{array}{}\text{(3)}& \mathrm{FR}={\displaystyle \frac{N_{\left(\mathrm{D}i\right)}/N_{\left(\mathrm{C}i\right)}}{N_{\left(\mathrm{D}\right)}/N_{\left(\mathrm{A}\right)}}},\end{array}$$

where N_(Di) is the number of debris flow pixels in the category i, N(Ci) is the total number of pixels in the category i, N_(D) is total number of debris flow pixels in the study area, and N_(A) is the total number of pixels in the study area.

Table 4 Factor categories and their fuzzy membership degrees.

The cosine amplitude method (Ross, 1995) is also widely used (Ercanoglu and Gokceoglu, 2004; Kanungo et al., 2006, 2009; Ercanoglu and Temiz, 2011) to establish relationships among elements of two or more datasets (Kritikos and Davies, 2015). Assuming that n is the number of data samples (categories of a factor used in the analysis) represented as an array X= {x₁,x₂, …, x_n} and that each of its elements, x_i, is a vector of length m (i.e. the size of the raster image) and can be expressed as $X=\mathit{\{}{x}_{i_{\mathrm{1}}},x_{i_{\mathrm{2}}}$, …, x_im}, then each element of a relation r_ij results from a pairwise comparison of a factor category x_i with a category of the debris flow distribution layer x_j (debris flow or non-debris flow). The memberships can be calculated by Eq. (4):

$$\begin{array}{}\text{(4)}& r_{ij}={\displaystyle \frac{\left|{\sum }_{k=\mathrm{1}}^m{x}_{ik}{x}_{jk}\right|}{\sqrt{\left({\sum }_{k=\mathrm{1}}^m{x}_{ik}^{\mathrm{2}}\right)\left({\sum }_{k=\mathrm{1}}^m{x}_{jk}^{\mathrm{2}}\right)}}}.\end{array}$$

As an analogy with the study of Kanungo et al. (2006), we defined the r_ij value for any given factor category as the ratio of the total number of debris flow pixels in the category to the square root of the product of the total number of pixels in that category and the total number of debris flow pixels in the area. Values of r_ij close to 1 indicate similarity, whereas values close to 0 indicate dissimilarity between the two datasets (Kritikos and Davies, 2015). What is more, every thematic layer must use the same pixel size to use the method properly.

3.5 DFSI map

To derive the debris flow susceptibility index (DFSI) map by overlaying the factor thematic layers using the fuzzy logic method, the “fuzzified” factors represented by information layers in raster format with values ranging from 0 to 1 need to be combined. Compared with the other four fuzzy operators, fuzzy gamma (Eq. 5) is more suitable for the research (Kritikos and Davies, 2015). To determine the appropriate γ value, the results of different gamma values were compared by the greatest distance (Kritikos and Davies, 2015) between the average DFSI curves of the debris flow locations and non-debris flow locations (for example, flat pixels) (Fig. 6). Finally, 0.9 is determined for the γ value because there is the greatest difference between debris flow and non-debris flow location areas. In order to illustrate the superiority of our model through a comparison, 17 results are calculated in ArcGIS.

$$\begin{array}{}\text{(5)}& {\mathit{\mu }}_{\left(x\right)}={\left(\mathrm{1}-{\prod }_{i=\mathrm{1}}^n\left(\mathrm{1}-{\mathit{\mu }}_i\right)\right)}^{\mathit{\gamma }}\cdot {\left({\prod }_{i=\mathrm{1}}^n{\mathit{\mu }}_i\right)}^{\mathrm{1}-\mathit{\gamma }},\end{array}$$

where μ_(x) is the combined membership value, μ_i is the fuzzy membership function for the ith map, i= 1,2, …, n are the numbers of thematic layers to be combined, and γ is a parameter in the range (0, 1).

Figure 6 Effect of γ value on debris flow susceptibility index (DFSI). Curves d, e and f correspond to debris flow pixels, and curves a, b and c correspond to non-debris flow area where a debris flow is unlikely. According to curve i, the maximum difference between the average DFSI values is observed for γ≈0.9.

To find the optimal model, 17 results were compared (Table 5). According to the distribution map of potential geological hazard points and susceptibility map in Pinggu District published by the Beijing Municipal Commission of Planning and Natural Resources (Beijing Municipal Commission of Planning and Natural Resources, 2022), three indexes are used to verify the validity and accuracy of the model.

Table 5 Predictive performance of different models.

Note: selected A factors with fuzzy membership more than 0.6; FRR represents the product of FR and r_ij; performance index is normalized by the largest FR value.

The results of the model are independent of the model itself, so the predictive performance of the final map is not just “the goodness of fit” of the data (Chung et al., 1995; Remondo et al., 2003). A relatively reliable technique for quantitatively assessing how good a model is is the construction of validation or success rate curves (Chung and Fabbri, 1999; van Westen et al., 2003; Remondo et al., 2003; Frattini et al., 2010) based on a comparison between the spatial distribution of debris flows and modeled debris flow susceptibility. The curves illustrate the debris flow recorded in the area with respect to susceptibility values also expressed as cumulative percentages of the total area. The area under the curve (AUC) defines the success rate (Marjanović et al., 2011). Generally, AUC values above 0.7 indicate model performance can be acceptable, while below 0.7, the performance is considered poor (Kritikos and Davies, 2015).

Although AUC is an effective evaluation method, the results are not comprehensive as mathematical features for selecting the best measurement model because of insufficient data for validation. In order to ensure the objectivity of the results, we can only effectively use the recorded debris flow gully as positive and the others as negative. Thus, a two-category test is proposed to verify the model in this paper. First, the DFSI map of each model is divided into two categories by the natural breaks (Jenks) method (Fig. 7). Then the accuracy ratio (AR) is defined as the frequency of the number of debris flows both classified by model and simultaneously recorded on site to the number of debris flows recorded on site. The resolution ratio (RR) is defined as the number of debris flows classified by model and simultaneously recorded on site to the total number of debris flows classified by the model (in red color). Take R₄, for example; there are in total 135 basins in the research area but only 46 records of debris flows (Fig. 3). In the results of two categories by the natural breaks (Jenks) method, 20 basins are divided into debris flow, while there are only 14 debris flows among them. Then 14 divided by 46 is AR, and 14 divided by 20 is RR.

Figure 7 Results of two categories by natural breaks (Jenks) method.

Figure 8 Debris flow susceptibility maps. Note: AUC results of R₁–R₄ below 0.7 are not shown. VH represents very high susceptibility area, H represents high susceptibility area, M represents moderate susceptibility area, L represents low susceptibility area, and VL represents very low susceptibility area.

The higher the two values are, the better the susceptibility map is. Finally, the performance of models (P value) can be obtained by Eq. (6). AUC values less than 0.6 are directly eliminated. Comparing the results of the rest of the models, the result of R₁₆ is optimal, and the results of the DFSI map are in good agreement with those of field investigation (Fig. 8).

$$\begin{array}{}\text{(6)}& P=\mathrm{AUC}+\sqrt{(\mathrm{AR}\cdot \mathrm{RR})}\end{array}$$

4 Results and discussion

Through the modeling process, relatively satisfactory results are obtained in this paper. The predictive performance of the output debris flow susceptibility maps, obtained from 17 different models, is verified by comparing them with maps published by relevant authorities. By comparing the results, the following results are discussed.

Firstly, comparing R₁, R₂, R₃, R₄ and R₅, it can be concluded that the model based on field investigation and expert experience is more effective than data driven directly when the information is insufficient. This is mainly because when the basin area reaches a certain size, it is no longer controlled by one or several factors but becomes a complex system. It is not only the factors that affect the system, but also the system will react to each factor. Geomorphic evolution is basically the result of the interaction of the endogenic and exogenic geological processes. A geological period can be regarded as the beginning of one endogenic geological process to the next one. In the early stage of the geological period, endogenic geological processes play a major role, and in the later relatively stable period, exogenic geological processes will take on more important parts. In this large cycle, a small cycle of energy accumulation and release occurs in the basin continuously, which leads to extremely complex system changes. In addition, there is a contradiction between the scale of geological evolution and the scale of engineering activities. So limited information can be obtained under these conditions, which leads to the unreliability of data-driven evaluation. Therefore, in the current period, field investigation and expert experience are fundamental.

Secondly, by comparing R₄ and R₅, R₆ and R₉, R₇ and R₁₀, R₈ and R₁₁, R₁₂ and R₁₅, R₁₃ and R₁₆, and R₁₄ and R₁₇, it can be concluded that the accuracy and resolution of the model can be improved by simplifying the factors, which will eliminate the ones with weak correlation and independence. In practical applications, even if the susceptibility map is obtained, the classification of the susceptibility degree is still a very difficult problem because everyone's subjective definition of “susceptibility degree” is different. By simplifying the factors, the main ones can be selected, which magnifies the differences between basins, so the boundaries between different susceptibility degrees are more obvious.

Thirdly, by comparing R₆ and R₁₂, R₇ and R₁₃, R₈ and R₁₄, R₉ and R₁₅, R₁₀ and R₁₆, and R₁₁ and R₁₇, it can be concluded that the appropriate classification of factors is helpful to optimize the susceptibility assessment model because the properties of the factors divided into one category are relatively consistent, as well as the impact on the debris flow system. We can also infer that the nonlinear combination characteristics between different types are stronger and scientific classification can improve the performance of the model.

Fourthly, comparing R₁₂ and R₁₃, as well as R₁₅ and R₁₆, it can be concluded that the frequency ratio method is better than the cosine amplitude method in the study. Different from the study of Kritikos and Davies (2015), the watershed unit rather than the grid unit is used, which indicates that the former has a wide range of applications, while the latter has a disadvantage of strict conditions.

Based on the results of the above four analyses, the most optimal model should have the features of being based on expert experience, using selected factors, classifying factors before using them and using the frequency ratio method. Then the model R₁₆ is selected according to the features, which is well in accordance with theoretical method performance score, and gets fine mutual verification.

There is also the selection of factors to discuss, which is still a very complex dilemma. Although 19 factors selected cannot fully evaluate the character of a basin, it is necessary to consider that they are easily and relatively accurately obtainable for each basin. This will facilitate a wide range of applications. Moreover, rainfall and total amount of loose material source are also very important influencing factors. But according to the Beijing hydrological manual, the rainfall change in the study area is not obvious, so it is excluded in the model. The total amount of loose material source cannot be obtained for the watershed without on-site investigation, so calculations are impossible. In fact, we indirectly consider the influence of natural loose material source by evaluating geological conditions but cannot consider the impact of human activities. As for the factors describing debris flow magnitude, usually several channels have the recorded data.

The scientific and systematic principle of model building is another challenge. To correctly classify the factors, it is necessary to grasp the characteristics of the formation, movement and accumulation of debris flow. Therefore, the classification should comprehensively consider the development background (geology, geomorphology, climate, hydrology, soil, vegetation, human activities and other factors). What the practical principle refers to is that the study should not only fully obtain scientific and accurate results but also make the professional results be understood by decision makers. Although the susceptibility grade and susceptibility value of each watershed are obtained, the results are relatively effective in this study area. In addition, with the development of technology and theory, we should replace some traditional factors which are not easy to quantify with more precise quantitative factors to improve the efficiency and accuracy of evaluation, such as surface roughness instead of drainage density.

We would like to further discuss the results derived from Table 3. It can be seen from the results that the occurrence of debris flow is highly correlated with basin volume, basin area and main gully bending coefficient with fuzzy membership above 0.7 in the Beijing area. Rainfall in the study area is abundant to induce the debris flow. Loose sources and sinks in the total volume of the catchment become more important. The watershed area determines the total volume of the catchment. For the same rainfall, generally, the larger the area is, the larger the catchment is. The bending coefficient reflects the replenishment sources along the channel. The greater the coefficient is, the slower the flow is. Then loose source along the channel has more time to replenish. Basin volume characterizes the maximum amount of loose material that can be supplied. These three features reflect the development characteristics of debris flow in the study area. It also provides ideas for disaster prevention and mitigation.

Finally, we should consider decision making under uncertainty because the debris flow phenomenon is extremely complex. The classification of geologists (high, moderate and low) is ambiguous for decision makers. It is more beneficial for them to use mathematically rigorous definitions. Considering that geological conditions tend to vary greatly from region to region, it is not appropriate to define a fixed limit. The Jenks method (chosen in this paper) can be used to classify sensitivity maps according to the characteristics of the data. We can also further process the data according to the needs of decision makers, such as identifying 10 % of the watersheds in the entire region as high risk. However, the applicability of the model to extreme rainfall and seismic conditions is not considered.

5 Conclusion

In this study, a new combination model for debris flow susceptibility based on GIS was developed in Pinggu. The objective and motivation of this study are to demonstrate a simple, extensible and convenient analytical model for the debris flow prediction. Three methods are selected in the model, each with their own advantages. GRA has great advantages in the case of fewer samples, the data-driven method is mainly used to reduce subjectivity, and fuzzy logic is fitted to solve nonlinear problems with fuzzy classification. The output optimal debris flow susceptibility maps demonstrated satisfactory performance with the relative higher susceptibility values corresponding to AUC = 0.768. The predictive performance of the susceptibility maps and the spatial correlation of debris flow gully with H and VH susceptibility with recorded debris flows illustrate that the assessment at regional scale using the proposed method is feasible. Compared with the previous results (Li et al., 2020b) based on grid units, the evaluation results are basically the same, but the model is more targeted at debris flow disasters for decision makers. Moreover, considering that the meaning of the factors used is clear and the data are easy to obtain, these conditions mentioned enable the model to be widely applied. In addition, a new factor (basin) is proposed in our study, which contributes greater weight – up to 0.79. From our 17 results by comparing the control variables, we suggest that other scholars should pay more attention to the classification and streamlining of factors, whose potential value to improve model accuracy has been indicated. It was also found that the watershed characteristic parameters can better reflect the advantages of the watershed unit, but further development is needed.

In short, an effort has been made to develop a cost- and time-efficient debris flow susceptibility assessment model. The model has an acceptable degree of accuracy for regional-scale planning and helps to make susceptibility and risk maps more accessible to individuals and local authorities. The GIS-based methods and modern data availability especially through online databases are significantly beneficial to this aim. However, a challenge remains in producing results with practical accuracy for the scale of planning using available resources. Previous studies highlight that the effectiveness of the final map depends on the quality of input data. Updating and improving existing debris flow catalogues and inventories are crucial for the development of reliable susceptibility and risk assessment methods.