Landslide susceptibility mapping based on Inter.iamb-Tabu algorithm considering non-landslide sampling

Abstract

Accurate landslide susceptibility mapping (LSM) remains challenging because of uncertainty in both model selection and non-landslide sample definition. To address these issues, this paper developed an LSM framework based on the Inter.iamb-Tabu algorithm while explicitly considering non-landslide sampling, using Gangu County, China, as the research area. 124 landslides were established through remote sensing interpretation and field survey. From 14 initial conditioning factors, nine were retained after correlation analysis, single-factor logistic regression and variance inflation factor analysis. Using 10 m×10 m grid cells as the mapping unit, a spatially random strategy was designed to select 26,468 non-landslide samples, and the datasets were divided into training and validation subsets at a ratio of 7: 3. Three improved Bayesian network algorithms, MMPC-Tabu, Fast.iamb-Tabu and Inter.iamb-Tabu, were then compared using Accuracy, Precision, Recall, F1-score and AUC. The results showed that Inter.iamb-Tabu achieved the best predictive performance. In addition, the proposed non-landslide sampling strategy outperformed two conventional methods, confirming its superiority in preserving spatial representativeness and reducing sample contamination. The final susceptibility map indicated that 87.1% of the documented landslides were in high and extreme susceptible areas. These areas were mainly distributed in the loess hilly region, the steep front margins of thick loess layers on the northern and southern mountains, and the secondary and tertiary gullies of the Weihe River valley. The learned DAG and CPTs revealed that slope gradient, elevation, land use and distance from fault were parent nodes of landslide and exerted direct triggering effects. Profile curvature and NDVI served as mutual feedback nodes with landslides, whereas slope aspect, distance from road and SPI mainly exerted indirect effects by interacting with other hazard factors. The proposed framework provides a robust, interpretable and transferable approach for LSM, and can support geohazard mitigation and land-use planning in regions with similar geological and geomorphological settings.

1 Introduction

Landslides occur when the shear stress on a penetrating structural surface exceeds its shear strength (Mao et al., 2022). According to the China Geological Environment Monitoring Institute, a total of 34,218 geological disasters occurred in China between 2020 and 2024, of which landslides accounted for more than 50%. On average, landslides cause over 200 fatalities and billions of yuan in economic losses annually in China (Huang et al., 2025). Landslide susceptibility mapping (LSM) estimates the probability of landslide occurrence using mathematical models and analyzes its spatial differentiation, providing a theoretical basis for landslide mitigation policies and land use planning (Shojaeezadeh et al., 2020; Yin et al., 2020).

Mathematical models for LSM generally fall into two categories (Sun et al., 2020): statistical models (e.g., logistic regression, analytic hierarchy process, and multivariate statistical methods) (Ansari et al., 2023) and machine learning methods (e.g., artificial neural networks, support vector machines, and convolutional neural networks (CNN)) (Aziz et al., 2024). Wang et al. (2024) classified landslides into rock and soil types based on geological conditions and historical inventories, proposed a CNN-based LSM method, and benchmarked its performance against classification and regression trees and multilayer perceptrons. Chen et al. (2025a) developed an LSM framework for Lin’an District, Zhejiang Province, China, integrating GeoDetector with ensemble learning algorithms, specifically random forest (RF) and extreme gradient boosting (XGBoost). Yang et al. (2026) integrated multi-source historical landslide data with 15 predictive factors and used several machine learning models (RF, gradient boosting regression trees, XGBoost, and categorical boosting) to generate susceptibility maps. To overcome the limitations of single models, some researchers have proposed hybrid algorithms. Owing to their high generalization ability and robustness, hybrid algorithms often yield significantly better LSM results than single algorithms (Jiang et al., 2026). Chen et al. (2025b) proposed a Buffer-SMOTE-Transformer optimization framework that integrates geospatial buffer sampling to refine negative sample selection, employs SMOTE to address class imbalance, and incorporates a weighted hybrid Transformer network to enhance modeling of complex geographic features. Gurung et al. (2026) integrated machine learning-based landslide susceptibility with numerical runout modeling to provide an LSM method for the Bhotekoshi watershed, overcoming the limitations of traditional models that focus solely on statistical susceptibility. Bayesian networks describe probabilistic dependencies among variables using a joint probability distribution (Sun et al., 2025). The max-min hill-climbing (MMHC) algorithm is a widely used Bayesian network method (Costache and Bui, 2019). MMPC-Tabu, Fast.iamb-Tabu, and Inter.iamb-Tabu are hybrid algorithms that improve MMHC with Tabu search (An et al., 2018). However, few studies have examined the effectiveness of Bayesian networks for LSM, making this an important direction for future research (Wu RZ. et al., 2025).

A key challenge in LSM is the uncertainty associated with non-landslide sample selection. Traditional selection methods often produce inaccurate susceptibility maps due to limitations in spatial representativeness and sample reliability (Wang et al., 2026). Chen et al. (2025a) introduced a novel framework combining multiple sampling, ensemble learning, and statistical analysis, generating N sets (N = 1, 10, 100, 500, 1,000) of random non-landslide samples to comprehensively assess the impact of sampling uncertainty on LSM. Ke et al. (2025) constructed 12 scenarios (three landslide sample types combined with four machine learning models) to assess the sensitivity of LSM to different non-landslide sampling strategies. Xu et al. (2025) creatively proposed an optimized non-landslide sampling method for the Wanzhou section of the Three Gorges Reservoir, China, selecting non-landslide samples from points with ground deformation rates between -5 mm/yr and +5 mm/yr in very-low susceptibility areas. Zhai et al. (2025) proposed an improved information quantity method for assessing landslide susceptibility in Yongfeng, South China, where non-landslide samples were randomly generated in very-low and low susceptibility areas at least 1,200 m away from any landslide sample. Although these non-landslide sample selection methods have significantly improved the accuracy of LSM results, precisely defining non-landslide samples remains a major challenge (Lusiana et al., 2026). On the one hand, if only a single attribute is controlled, such as sampling from areas with a slope gradient <5°, the resulting samples have homogeneous characteristics, which may bias the assessment of factor importance (Ismailov et al., 2025). On the other hand, the stability of some selected non-landslide samples varies; they have simply not experienced sliding during a specific period and may still be susceptible to failure under extreme conditions in the future. Absolute labeling fails to reflect this continuity (Gu et al., 2024). Continuing to conduct research on non-landslide sample selection based on existing studies is of great significance for further improving LSM accuracy.

Geological environmental characteristics and landslide mechanisms vary considerably across regions, and the logical structures of mathematical models also differ. Consequently, the optimal LSM model for a given area cannot be known in advance and must be determined through comparative studies (Erzagian et al., 2025). Bayesian networks reveal direct and indirect causal relationships between landslides and hazard factors through directed acyclic graphs (DAGs) and conditional probability tables (CPTs), identifying parent nodes and mutual feedback nodes of landslides (Shang et al., 2024). The sampling strategy for non-landslide samples has an even greater impact on the accuracy of LSM than the choice of machine learning model itself. This paper applied three improved Bayesian network-based algorithms (MMPC-Tabu, Fast.iamb-Tabu, and Inter.iamb-Tabu) for LSM in the study area, selected the optimal algorithm based on multiple evaluation metrics, and revealed the influence of each hazard factor on landslide occurrence. In addition, a novel non-landslide sample selection method was proposed, and its effectiveness in improving LSM accuracy was validated. The research flowchart is shown in Figure 1.

FIGURE 1

2 Research area and data

2.1 Research area overview

Gangu County is in southeastern Gansu Province, China, within the transitional zone between the Longxi Loess Plateau and the Qinling Mountains. It lies between 104°58′-105°31′E and 34°31′-35°03′N, covering a total area of 1,572.6 km², as shown in Figure 2.

FIGURE 2

The region experiences a continental monsoon climate, with an average annual precipitation of 437.3 mm. Rainfall occurs predominantly from July to September, accounting for 51% of the annual total, and exhibits significant interannual variability. The terrain is generally higher in the south and lower in the north, with an average elevation of 1,972 m, ranging from 1,228 m to 2,716 m, giving a relative elevation difference of 1,488 m (Wu YY. et al., 2025). The southern part, formed by the western extension of the Qinling Mountains, is a low-to-medium mountain area underlain by bedrock; the northern part consists of remnant ranges of the Liupan Mountains, characterized by loess hilly and gully terrain. The total river length is 131.1 km, with the Weihe River flowing across the entire region from west to east as the main watercourse. Prolonged fluvial incision has created an undulating landscape of ridges, hillocks, gullies, and valleys. The alluvial plains along the Weihe River are relatively flat, forming a geomorphic pattern of “extensive mountains and limited valleys”. Groundwater in the area consists of unconsolidated rock pore water and clastic rock pore-fissure water, with valley phreatic water serving as the main water supply source (Liu et al., 2023).

The geological structure of Gangu County is complex and exhibits strong tectonic activity. The regional structural lineation is generally oriented northwest-west (NWW), consistent with the course of the Weihe River. Folds and faults are well developed, particularly the nearly east-west trending fault zones distributed along both banks of the Weihe River (Li et al., 2024).

2.2 Research data

The data used in this paper are shown in Table 1.

2.3 Landslides of Gangu County

To provide training and validation datasets for LSM, landslides were identified through remote sensing interpretation and field surveys in Gangu County (Wu et al., 2026). The data used were Landsat 8 OLI Lower Truncation images, acquired on 2025-10-09 (02:42:12).

2.3.1 Pre-processing of the remote sensing image

FIGURE 3

2.3.2 Establishment of interpretation signs

2.3.3 Landslide remote sensing interpretation results

Areas matching the above interpretation signs were extracted using e-Cognition 9.0 software. Combined with field surveys, a total of 124 landslides were identified. The results show that Gupo Town has the highest number of landslides (33). The total volume of landslides in Gangu County is 12,923,641 m³, and the total area is 2,645,940 m². The largest landslide is in Gupo Town, with a volume of 241,346 m³ and an area of 37,418 m². The DEM of the research area was resampled into 10 m×10 m grids using ArcGIS 10.2, resulting in a total of 15,726,358 grids, of which 26,468 correspond to landslide grids. The landslides in Gangu County are shown in Figure 2 and Table 2.

Landslide No. 1 in Table 2 is provided as an example. This landslide is in Jiacang Village, Gupo Town, with a width of approximately 112.3 m, a height of about 87.2 m, and a slope gradient of 76°. The depth of the sliding surface ranges from 2.8 m to 5.1 m. The landslide mass consists of typical Quaternary loess with medium collapsibility, covering an area of 5,124 m² with a volume of 21,110 m³. A large-scale sliding occurred in August 2022, followed by small-scale collapses during each subsequent rainy season. The surface of the landslide mass is exposed, with almost no vegetation cover. Groundwater seeps from the lower part of the landslide mass year-round, and the moisture content is significantly higher than in the upper part. Currently, cracks are developing at the rear scarp and on both sides of the landslide mass. Although the cracks on the sides and rear scarp have not yet connected, they are wide and have become infiltration pathways for surface water and groundwater. Under the influence of heavy rainfall, a large-scale sliding is highly likely to occur again, as shown in Figure 4.

FIGURE 4

2.4 Screening of hazard factors

Landslide hazard factors include slope gradient, elevation, slope aspect, land use, lithology, distance from fault, distance from river, distance from road, plane curvature, profile curvature, normalized difference vegetation index (NDVI), topographic wetness index (TWI), sediment transport index (STI) and stream power index (SPI) (Khalaj et al., 2020).

Plane curvature represents the bending of a point on the ground surface along its contour line, i.e., the horizontal curvature component. Profile curvature is the rate of elevation change along the direction of maximum descent, i.e., the vertical curvature component (Yilmaz, 2009). NDVI assesses vegetation condition by calculating the reflectance difference between visible and near-infrared bands (Aditian et al., 2018). TWI reflects the influence of regional topography on runoff flow and accumulation, quantifying topographic control over basic hydrological processes. SPI represents the erosive capacity of surface water flow and can be used to identify strong flow paths formed by water convergence and locations where gully erosion may occur. STI is a comprehensive topographic variable indicating the extent of surface sand and other materials transported by water flow (Wang et al., 2019). The calculation methods for NDVI, TWI, SPI and STI are shown in Equations 1–4.where: NIR is the near-infrared band reflectance; Red is the red band reflectance; A_s is the upstream area of surface water flowing per unit contour length, calculated from the cumulative confluence area and upstream flow length; β is the surface gradient.

2.4.1 Correlation analysis

In Bayesian networks, an excessive number of nodes increases the complexity of structure learning, and landslide hazard factors are not completely independent (Michalowski and Park, 2020). To eliminate correlation among hazard factors and improve modeling efficiency, Pearson correlation coefficients (δ_xy) were calculated for 12 quantitative hazard factors using SPSS 23.0, as shown in Figure 5.

FIGURE 5

Four pairs of hazard factors exhibited δ_xy > 0.4: STI with slope gradient; TWI with distance from fault; TWI with distance from river; and distance from river with plane curvature. Therefore, STI, TWI and distance from river were excluded.

2.4.2 Single-factor logistic regression

Single-factor logistic regression is a generalized linear model that assesses the individual association between each factor and the outcome variable, thereby screening for statistically significant factors for subsequent analysis (Conforti and Ietto, 2021). This method was further used to screen landslide hazard factors, with calculations performed using SPSS 23.0. The results are shown in Table 3.

Where: β₁ is the regression coefficient: positive values indicate that the probability of the event increases with the independent variable, while negative values indicate the opposite; OR (Odds Ratio) represents the multiplicative change in the odds of the event for each one-unit increase in the independent variable; p is the significance value, with p < 0.05 indicating a statistically significant individual association between the independent variable and the outcome variable (Kim et al., 2021). As the p-values for lithology and plane curvature were 0.127 and 0.096, respectively, these two hazard factors were excluded. It should be particularly noted that although lithology is crucial to landslide formation, the assessment model cannot reflect regional differences in the effect of lithology on landslide occurrence because the land surface in the study area is entirely covered by Quaternary loess. This is the main reason for the exclusion of lithology.

2.4.3 Variance inflation factor analysis

Variance inflation factor (VIF) is a measure of the severity of complex collinearity in multiple linear regression models. It represents the ratio of the variance of the regression coefficient estimator to the variance when the independent variables are assumed not to be linearly correlated. When the correlation between independent variables is very low, VIF is close to 1. The larger the VIF, the larger the collinearity. VIF > 10 means there is a serious collinearity problem in the variables, and VIF is shown in Equation 5.where: R_i² is the correlation coefficient of x_i for regression analysis of the other variables.

SPSS23.0 was used to test the collinearities of the 9 remaining hazard factors, and the results are shown in Table 4.

The VIF of each hazard factor is <10, and the collinearity is low. Consequently, the finally selected landslide hazard factors include: slope gradient, elevation, slope aspect, land use, distance from fault, distance from road, profile curvature, NDVI and SPI.

2.5 Classification of hazard factors

In this paper, c_n was initially set to 5, with initial breakpoints at 18°, 36°, 54° and 72°. The results show that when c_n = 8 and the breakpoints are 5.111°, 9.756°, 13.705°, 17.422°, 21.371°, 26.017°, 32.521° and 59.236°, Var_a is minimized, representing the optimal number of classifications and breakpoints, as shown in Table 5 and Figure 6.

FIGURE 6

3 Methodology

3.1 Bayesian network

A Bayesian network uses a joint probability distribution to describe probabilistic dependencies among variables (He et al., 2022). Given a set of random variables X = {X₁, …, X_n}, their joint probability P (X₁, …, X_n) can be represented as a Bayesian network B=(G, θ), where G is a DAG. In this graph, nodes correspond to random variables, and directed edges represent probabilistic dependencies. A directed edge from X_i to X_j indicates that X_i is a parent of X_j, and X_j is a child of X_i. θ represents the CPTs, which quantify the relationship between each variable and its parents. Constructing a Bayesian network involves three steps: (1) identifying relevant variables and their possible values; (2) obtaining an optimal network structure using machine learning algorithms (structure learning); and (3) computing the CPTs for each node (parameter learning) (Bakhouya et al., 2026).

3.1.1 Structure learning

Structure learning aims to determine the DAG of a Bayesian network by analyzing dependencies among variables. This task is considerably more difficult than parameter learning because the number of possible DAGs grows super-exponentially with the number of nodes. Structure learning methods fall into two categories: constraint-based and score-based. Constraint-based methods use conditional independence tests to decide whether an edge should exist between two variables. Their theoretical foundation is the Markov property of Bayesian networks, which states that each node is conditionally independent of its non-descendants given its parents. Conversely, if two variables are not directly connected by an edge, they should be conditionally independent given some other set of variables. The most classical constraint-based algorithm is the Peter-Clark algorithm (Li et al., 2021).

Score-based methods formulate structure learning as an optimization problem. They define a scoring function that measures how well a given graph structure fits the observed data and then search the space of possible graphs to maximize this score. The most widely used scoring functions are the Bayesian Information Criterion (BIC) and the Bayesian Dirichlet equivalent uniform (BDeu) score. Due to the vast space of possible DAGs, heuristic search strategies are essential; the most popular is hill climbing. To overcome local optima, hill climbing is often combined with random restarts (i.e., initiating the search from multiple different starting graphs) or with Tabu search (i.e., forbidding the algorithm from revisiting recently explored graphs) (Hung et al., 2025).

Hybrid methods typically outperform pure constraint-based methods because the subsequent scoring phase can correct errors from independence tests, and they are faster than pure score-based methods due to a drastically reduced search space (Li et al., 2022). A representative hybrid algorithm is MMHC. The MMHC algorithm consists of two stages. In the first stage, a max-min heuristic search strategy is employed. Specifically, the MaxMinHeuristic function obtains the Candidate Parents and Children (CPC) set for each variable (Men et al., 2025). This function calculates the minimum association value between the target node T and every other node, and then selects the node corresponding to the maximum of these minimum association values to be included in the CPC. The first stage terminates when, given all subsets of the CPC, the remaining nodes are all independent of T. In the second stage, the Ind(X; T|Z) function removes any nodes from the CPC that were erroneously included in the first stage. If node X and target node T are independent given a conditioning set Z, then X is removed from the CPC. The Ind(X; T|Z) function returns true if X and T are conditionally independent given Z.

3.1.2 Parameter learning

Two fundamental methods exist for parameter estimation in Bayesian networks: maximum likelihood estimation (MLE) and Bayesian estimation (Natalia et al., 2025).

3.1.2.1 Maximum likelihood estimation

Maximum likelihood estimation assumes that X₁, X₂, X₃, X₄ are samples drawn from X (Wu et al., 2020). The joint distribution of X₁, X₂, X₃, X₄ is given by Equation 6.where: θ is the parameter to be estimated, θ∈Θ, and Θ is the parameter space (the range of possible values of θ).

Let x₁, x₂, …, x_n be the observed sample values corresponding to X₁, X₂, …, X_n. The probability that the sample X₁, X₂, …, X_n takes the specific observations x₁, x₂, …, x_n is given by Equation 7.

L (θ) is the likelihood function of the sample. Maximum likelihood estimation seeks to fix the sample observations x₁, x₂, … x_n and select the parameter value ∈Θ that maximizes L(θ), as expressed in Equation 8.

3.1.2.2 Bayesian estimation

Bayesian estimation treats θ as a random variable and computes its posterior probability distribution. A prior probability distribution P(θ) is chosen to summarize prior knowledge about θ. The effect of the observed data X=(X₁, X₂, …, X_n) is then summarized by the likelihood function . The prior distribution and the likelihood function are combined using Bayes’ formula to obtain the posterior distribution (John et al., 2025). Bayes’ theorem and Bayesian estimation are shown in Equations 9, 10, respectively.

3.1.3 Improvement of the MMHC algorithm

MMPC-Tabu, Fast.iamb-Tabu and Inter.iamb-Tabu are improved algorithms based on MMHC. Their core idea is to use constraint-based algorithms to identify the CPC for each node, thereby greatly reducing the search space. Then, with the candidate set as a constraint, Tabu search is performed within the reduced space to find a DAG that optimizes the network score.

3.1.3.1 MMPC-Tabu

MMPC (max-min parents and children) is the core algorithm used in the skeleton discovery phase of MMHC. The MMPC-Tabu algorithm consists of three stages (Zhong et al., 2022):

Forward selection: The max-min heuristic is applied to progressively add candidate nodes to the CPC until no further node can pass conditional independence tests.

Backward elimination: For each node in the CPC, it is checked whether it is conditionally independent of T given the other nodes in the CPC. If yes, the node is removed from the CPC.

Tabu search: A Tabu List is maintained to record the most recently visited network structures or recently performed operations (e.g., add edge, delete edge, reverse edge).

The initial parameters used in constructing the MMPC-Tabu framework are shown in Table 6.

3.1.3.2 Fast.iamb-Tabu

Fast.iamb (fast incremental association Markov blanket) belongs to the IAMB algorithm family. Its core idea is to construct the Markov blanket incrementally. The main difference between Fast.iamb and MMPC lies in the heuristic: Fast.iamb employs a maximum association heuristic (typically using mutual information or absolute correlation coefficient) instead of the max-min heuristic. That is, at each step it simply selects the node most strongly associated with T, without considering the minimum association over all conditional subsets as in MMPC (Fallahnejad et al., 2025). The initial parameters used in constructing the Fast.iamb-Tabu framework are shown in Table 7.

3.1.3.3 Inter.iamb-tabu

Inter.iamb is also based on the IAMB framework but improves upon classic IAMB’s tendency to accumulate false positives. It adopts an interleaved procedure: after each forward-selection step that adds one variable, a backward-elimination step is immediately performed. This allows false positives to be identified and removed early, preventing them from contaminating subsequent selection steps. Consequently, Inter.iamb-Tabu is more robust than Fast.iamb-Tabu, especially when the sample size is small or when complex dependencies exist among variables (Zhou et al., 2022). The initial parameters used in constructing the Inter.iamb-Tabu framework are shown in Table 8.

3.2 Establishment of training and validation datasets

70% of the landslide samples (18,528) and 70% of the non-landslide samples (18,528) were randomly selected as the training dataset, with the remaining 30% of landslide samples (7,940) and 30% of non-landslide samples (7,940) serving as the validation dataset.

3.3 Selection of the optimal modeling method

Multiple assessment metrics were used to select the optimal modeling method, including Accuracy, Precision, Recall, F1-Score, the Receiver Operating Characteristic (ROC) curve and the area under the curve (AUC) (Sachinthaka et al., 2025). The calculation methods are shown in Equations 11–15.where: TP is the number of landslide grids correctly predicted, FP is the number of landslide grids incorrectly predicted, TN is the number of non-landslide grids correctly predicted, FN is the number of non-landslide grids incorrectly predicted.

4 Research results

The MMPC-Tabu, Fast.iamb-Tabu and Inter.iamb-Tabu algorithms were trained using the training dataset. The trained models were then applied to the validation dataset to compute output values for each sample. The modeling hardware environment consisted of a CPU i7-6700 processor, 8 GB RAM and a GTX1050 Ti-8G graphics card; the software environment was the Bayesian network learning package in R3.5.0.

In the ROC curve, the x-axis is the FPR and the y-axis is the Recall. The ROC curves for the three algorithms are shown in Figure 7. The calculated values of TP, FP, TN, FN, Accuracy, Precision, Recall, F1-Score and AUC are shown in Table 9.

FIGURE 7

Considering all assessment metrics, the Inter.iamb-Tabu algorithm achieved the best modeling performance, followed by Fast.iamb-Tabu, while MMPC-Tabu performed the worst.

The trained MMPC-Tabu, Fast.iamb-Tabu and Inter.iamb-Tabu algorithms were used to calculate the probabilities of landslide occurrence for 15726358 grids in the research area. Although the study area is prone to landslides, the landslide area accounts for only 0.168% of the total area. Based on historical hazard experience and natural conditions, the area of regions that are extremely prone to landslides generally does not exceed 10% of the total area (Gebru and Shewangzaw, 2026). On the other hand, from the perspective of landslide prevention and control, to improve the efficiency of limited resources (Bednarik et al., 2024), the delineated High and Extreme susceptible areas should not be too large. Therefore, the probabilities were classified into five levels (from low to high) so that the area percentages of Minimal, Minor, Medium, High and Extreme susceptible areas were 40%, 30%, 15%, 10% and 5%, respectively. The resulting landslide susceptibility maps are shown in Figure 8. The distribution of landslide grids across susceptible levels is shown in Table 10.

FIGURE 8

Based on the optimal Inter.iamb-Tabu algorithm, the Minimal, Minor, Medium, High and Extreme susceptible areas accounted for 0.62%, 6.14%, 11.82%, 20.17% and 61.25% of the total research area, respectively, and among the 124 landslides, the numbers falling into these susceptible levels were 1, 4, 11, 26 and 82, respectively. Topographically, the High and Extreme susceptible areas are mainly distributed in the loess hilly region, the steep front margins of thick loess layers on the northern and southern mountains, and the secondary and tertiary gullies of the Weihe River valley (e.g., Goumen Village in Gupo Town, Yijiawan in Jinshan Town and Maiduiping in Liufeng Town). Lithologically, they are primarily distributed in areas where Quaternary loess and Neogene red sandstone and claystone are exposed.

5 Discussions

5.1 Relationship between landslide and hazard-prone environment

Based on the DAG and CPTs generated by the Inter.iamb-Tabu algorithm, the relationship between landslide and element of the hazard-prone environment is analyzed, as shown in Figure 9.

FIGURE 9

The DAG comprises 10 nodes (nine hazard factor nodes and one terminal node), their corresponding CPTs, and 15 directed edges. Slope gradient, elevation, land use and distance from fault are parent nodes of landslide, exerting direct triggering effects on landslide occurrence. For example, slope gradient is classified into eight grades: Grade I (0°–5.111°), Grade II (5.111°–9.756°), Grade III (9.756°–13.705°), Grade IV (13.705°–17.422°), Grade V (17.422°–21.371°), Grade VI (21.371°–26.017°), Grade VII (26.017°–32.521°) and Grade VIII (32.521°–59.236°), with conditional probabilities of 0.0396, 0.0526, 0.0845, 0.1039, 0.1693, 0.2051, 0.3172 and 0.0278, respectively, indicating that landslides are most likely to occur when slope gradient falls within the range of 26.017°–32.521°. Profile curvature and NDVI serve as mutual feedback nodes with respect to landslides. In addition to directly triggering landslides, landslide occurrence can also drive the evolution of the hazard-prone environment, thereby altering profile curvature and NDVI. Slope aspect, distance from road and SPI exhibit mutual feedback relationships with other hazard factors and play an indirect role in landslide occurrence. For instance, distance from road can induce changes in slope gradient and slope aspect and slope aspect can subsequently influence land use, thereby triggering landslides.

5.2 Comparison of non-landslide sample selection methods

The method proposed in Section 4.2 significantly outperforms Traditional Methods 1 and 2. The underlying reason is that although the traditional methods avoid the misclassification of landslide grids, they restrict certain hazard factors to specific value ranges, thereby losing generalizability. For example, Traditional Method 1 selects non-landslide grids based on field surveys, which cannot cover deep mountainous areas characterized by high slope gradient, high elevation and high NDVI. However, high values of slope gradient, elevation and NDVI do not necessarily lead to landslides. As another example, in Traditional Method 2, areas with large surface deformation are not necessarily landslide-prone zones; rockfalls, debris flows and land subsidence can also cause significant surface deformation.

6 Conclusion

References

• 1

  AditianA.KubotaT.ShinoharaY. (2018). Comparison of GIS-based landslide susceptibility models using frequency ratio, logistic regression, and artificial neural network in a tertiary region of ambon, Indonesia. Geomorphology318, 101–111. 10.1016/j.geomorph.2018.06.006

  • CrossRef
  • Google Scholar
  • View reference in article

• 2

  AnN.TengY.YangJ. Y.LiL. (2018). Bayesian network structure learning method based on causal effect. Appl. Res. Comput.35 (12), 3609–3613.

  • Google Scholar
  • View reference in article

• 3

  AnsariA.RaoK. S.JainA. K.AnsariA. (2023). Formulation of multi-hazard damage prediction (MhDP) model for tunnelling projects in earthquake and landslide-prone regions: a novel approach with artificial neural networking (ANN). J. Earth Syst. Sci.132 (4), 164–175. 10.1007/s12040-023-02178-y

  • CrossRef
  • Google Scholar
  • View reference in article

• 4

  AzizK.MirR. A.AnsariA. (2024). Precision modeling of slope stability for optimal landslide risk mitigation in ramban road cut slopes, Jammu and Kashmir (J&K) India. Model. Earth Syst. Environ.10 (3), 3101–3117. 10.1007/s40808-023-01949-2

  • CrossRef
  • Google Scholar
  • View reference in article

• 5

  BakhouyaM.ElghaziK.RamchounH.HaddaM.MasrourT. (2026). Variational bayesian gaussian mixture last layer of convolutional neural network. SN Comput. Sci.7 (2), 131–150. 10.1007/s42979-026-04736-9

  • CrossRef
  • Google Scholar
  • View reference in article

• 6

  BednarikM.YilmazI.KralovičováL. (2024). Deterministic approach to assess landslide susceptibility and landslide activity in the central-western region of Slovakia. Bull. Eng. Geol. Environ.83, 327. 10.1007/s10064-024-03795-7

  • CrossRef
  • Google Scholar
  • View reference in article

• 7

  ChenK.FangH. L.JiangJ. (2025a). Landslide susceptibility prediction in Lin’an district, China, using ensemble learning with non-landslide sample uncertainty. Nat. Hazards121, 24347–24372. 10.1007/s11069-025-07708-z

  • CrossRef
  • Google Scholar
  • View reference in article

• 8

  ChenL.JiaY. C.HeS.WuS. J. X.WangF. (2025b). A hybrid neural network-based model for landslide susceptibility. Earth Syst. Environ., 1–20.

  • Google Scholar
  • View reference in article

• 9

  ConfortiM.IettoF. (2021). Modeling shallow landslide susceptibility and assessment of the relative importance of predisposing factors, through a GIS-based statistical analysis. Geosciences11 (8), 333–361. 10.3390/geosciences11080333

  • CrossRef
  • Google Scholar
  • View reference in article

• 10

  CostacheR.BuiD. T. (2019). Spatial prediction of flood potential using new ensembles of bivariate statistics and artificial intelligence: a case study at the putna river catchment of Romania. Sci. Total Environ.691, 1098–1118. 10.1016/j.scitotenv.2019.07.197

  • CrossRef
  • Google Scholar
  • View reference in article

• 11

  ErzagianE.FathaniT. F.WilopoW. (2025). Statistical landslide susceptibility zonation: Kulon progo Mountains, Indonesia. Environ. Earth Sci.84 (21), 601–620. 10.1007/s12665-025-12519-3

  • CrossRef
  • Google Scholar
  • View reference in article

• 12

  FallahnejadM.RezaeitabarV.KazemiM. (2025). A generalized K2 algorithm for learning bayesian network structures using ridge regression. Lobachevskii J. Math.46 (4), 1567–1579. 10.1134/s1995080225606290

  • CrossRef
  • Google Scholar
  • View reference in article

• 13

  GebruS. G.Shewangzaw (2026). Reducing subjectivity in landslide susceptibility mapping method: an integrated methodology. Geoenvironmental Disasters13, 15. 10.1186/s40677-026-00365-9

  • CrossRef
  • Google Scholar
  • View reference in article

• 14

  GuT. F.DuanP.WangM. G.LiJ.ZhangY. K. (2024). Effects of non-landslide sampling strategies on machine learning models in landslide susceptibility mapping. Sci. Rep.14 (1), 7201–7217. 10.1038/s41598-024-57964-5

  • CrossRef
  • Google Scholar
  • View reference in article

• 15

  GurungB.ChenN. S.HuG. S.KhadkaN.SapkotaN.GouliM. R.et al (2026). Integrating machine learning and numerical methods for enhanced landslide susceptibility and hazard mapping in the bhotekoshi watershed, central Nepal. J. Mt. Sci.23, 1–21. 10.1007/s11629-025-9951-2

  • CrossRef
  • Google Scholar
  • View reference in article

• 16

  HeD. L.ChengY.ZhaoR. L. (2022). Study on learning Bayesian network structure based on MMHC algorithm. J. Beijing Technol. Bus. Univ. Nat. Sci. Ed.26 (3), 43–48. (in Chinese).

  • Google Scholar
  • View reference in article

• 17

  HuangF. M.YangY.JiangB. C.ChangZ. L.ZhouC. B.JiangS. H.et al (2025). Effects of different division methods of landslide susceptibility levels on regional landslide susceptibility mapping. Bull. Eng. Geol. Environ.84 (6), 276. 10.1007/s10064-025-04281-4

  • CrossRef
  • Google Scholar
  • View reference in article

• 18

  HungE.MantziouA.ReinertG. (2025). A Bayesian mixture model for poisson network autoregression. Soc. Netw. Analysis Min.15 (1), 70–88. 10.1007/s13278-025-01485-0

  • CrossRef
  • Google Scholar
  • View reference in article

• 19

  IsmailovV. A.KhusomiddinovA. S.MukhammadkulovN. M.YadigarovE. M.KhayriddinovB. B.RakhmatovA. R.et al (2025). Assessment of landslide dynamics based on microseismic parameters. Seism. Instrum.61, 346–353. 10.3103/s0747923925700446

  • CrossRef
  • Google Scholar
  • View reference in article

• 20

  JiangH. W.ZhouH.WuJ. Y.LiuM. J.WuY. X.GuoY. F. (2026). Ensemble learning for landslide susceptibility mapping: a review of machine learning and hybrid approaches. Carbonates Evaporites41 (1), 16. 10.1007/s13146-025-01221-x

  • CrossRef
  • Google Scholar
  • View reference in article

• 21

  JohnN.LimA.SuntharS. R.ZhangL.ChengJ.BoktorJ.et al (2025). Immunotherapies in neuromyelitis optica: bayesian network meta-analysis. J. Neurology272 (9), 563–573. 10.1007/s00415-025-13279-7

  • CrossRef
  • Google Scholar
  • View reference in article

• 22

  KaushalA.GuptaA. K.SehgalV. K. (2024). A semantic segmentation framework with UNet-pyramid for landslide prediction using remote sensing data. Sci. Rep.14 (1), 30071–30094. 10.1038/s41598-024-79266-6

  • Pubmed Abstract
  • CrossRef
  • Google Scholar
  • View reference in article

• 23

  KeC. Y.SunP.ZhangS.LiR.SangK. Y. (2025). Influences of non-landslide sampling strategies on landslide susceptibility mapping: a case of tianshui city, northwest of China. Bull. Eng. Geol. Environ.84 (3), 123–142. 10.1007/s10064-025-04147-9

  • CrossRef
  • Google Scholar
  • View reference in article

• 24

  KhalajS.ToroodyF. B.AbaeiM. M.BahootoroodyA.AbbassiR. (2020). A methodology for uncertainty analysis of landslides triggered by an earthquake. Comput. Geotechnics117, 103262. 10.1016/j.compgeo.2019.103262

  • CrossRef
  • Google Scholar
  • View reference in article

• 25

  KimH.LeeJ. H.ParkH. J.HeoJ. H. (2021). Assessment of temporal probability for rainfall-induced landslides based on nonstationary extreme value analysis. Eng. Geol.294, 106372. 10.1016/j.enggeo.2021.106372

  • CrossRef
  • Google Scholar
  • View reference in article

• 26

  LiM.ZhangR.LiuK. F. (2021). A new marine disaster assessment model combining bayesian network with information diffusion. J. Mar. Sci. Eng.9 (6), 640–657. 10.3390/jmse9060640

  • CrossRef
  • Google Scholar
  • View reference in article

• 27

  LiT. T.ZhouY. Z.ZhaoY.ZhangC.ZhangX. (2022). A hierarchical object-oriented Bayesian network-based fault diagnosis method for building energy systems. Appl. Energy306 (B), 118088. 10.1016/j.apenergy.2021.118088

  • CrossRef
  • Google Scholar
  • View reference in article

• 28

  LiX. B.ZhouX. H.PengD.OuyangG. L.WuY. W.YanL. Y. (2024). Distribution and characteristics of loess landslides induced by the 1718 tongwei earthquake based on a field survey. Bull. Eng. Geol. Environ.83 (5), 163–179. 10.1007/s10064-024-03666-1

  • CrossRef
  • Google Scholar
  • View reference in article

• 29

  LiuM. X.MiJ. L.WangS. Y.XiaoS. R.LiL.XiaoY. D. (2023). Characteristics of soil organic carbon distribution in different economic forests in gangu county, Gansu province, China. Eurasian Soil Sci.56 (11), 1641–1652. 10.1134/s1064229323700242

  • CrossRef
  • Google Scholar
  • View reference in article

• 30

  LusianaN.AdliyaG. E.DeviantoL. A.HusinN. A. (2026). Rapid post-landslide vegetation regrowth detected by multi-temporal satellite imagery in the southern part of mt. Rinjani national park, lombok, Indonesia. Nat. Hazards122, 191. 10.1007/s11069-025-07964-z

  • CrossRef
  • Google Scholar
  • View reference in article

• 31

  ManiA.KumariM.BadolaR. (2024). Landslide hazard zonation (LHZ) mapping of doon valley using multi-criteria analysis method based on remote sensing and GIS techniques. Discov. Geosci.2 (1), 35–53. 10.1007/s44288-024-00044-y

  • CrossRef
  • Google Scholar
  • View reference in article

• 32

  MaoY. M.MwakapesaD. S.LiY. C.XuK. B.NanehkaranY. A.ZhangM. S. (2022). Assessment of landslide susceptibility using DBSCAN-AHD and LD-EV methods. J. Mt. Sci.19 (1), 184–197. 10.1007/s11629-020-6491-7

  • CrossRef
  • Google Scholar
  • View reference in article

• 33

  MengT.ZhangZ. Y.ZhangX.ZhangC. (2025). Bayesian network for predicting mandibular third molar extraction difficulty. BMC Oral Health25 (1), 56–64. 10.1186/s12903-025-05432-5

  • Pubmed Abstract
  • CrossRef
  • Google Scholar
  • View reference in article

• 34

  MichalowskiR. L.ParkD. (2020). Stability assessment of slopes in rock governed by the hoek-brown strength criterion. Int. J. Rock Mech. Min. Sci.127, 104217. 10.1016/j.ijrmms.2020.104217

  • CrossRef
  • Google Scholar
  • View reference in article

• 35

  NataliaY. A.CauwerH. D.NeyensT.GoniewiczK.SomvilleF.MolenberghsG. (2025). A Bayesian network analysis of aviation terrorism attack risks. J. Transp. Secur.18 (1), 19–36. 10.1007/s12198-025-00315-w

  • CrossRef
  • Google Scholar
  • View reference in article

• 36

  QianL. X.LiuN. J.HongM.DangS. Z. (2024). An improved nonlinear dynamical model for monthly runoff prediction for data scarce basins. Stoch. Environ. Res. Risk Assess.38 (10), 3771–3798. 10.1007/s00477-024-02773-5

  • CrossRef
  • Google Scholar
  • View reference in article

• 37

  SachinthakaR.KalatehjariR.BrookM. S. (2025). Evolution and critical evaluation of deterministic physically based rainfall-induced landslide susceptibility mapping: a mixed review. Nat. Hazards121 (18), 20795–20818. 10.1007/s11069-025-07634-0

  • CrossRef
  • Google Scholar
  • View reference in article

• 38

  ShangH.LiuS. H.ZhongJ. X.TsangaratosP.IliaI.ChenW.et al (2024). Application of naive bayes, kernel logistic regression and alternation decision tree for landslide susceptibility mapping in pengyang county, China. Nat. Hazards120, 12043–12079. 10.1007/s11069-024-06672-4

  • CrossRef
  • Google Scholar
  • View reference in article

• 39

  ShojaeezadehS. A.NikooM. R.MirchiA.MallakpourI.AghaKouchakA.SadeghM. (2020). Probabilistic hazard assessment of contaminated sediment in Rivers. Sci. Total Environ.703, 134875. 10.1016/j.scitotenv.2019.134875

  • Pubmed Abstract
  • CrossRef
  • Google Scholar
  • View reference in article

• 40

  SunD. L.WenH. J.WangD. Z.XuJ. H. (2020). A random forest model of landslide susceptibility mapping based on hyperparameter optimization using bayes algorithm. Geomorphology362, 107201. 10.1016/j.geomorph.2020.107201

  • CrossRef
  • Google Scholar
  • View reference in article

• 41

  SunB. D.ZhangX. Y.JiangJ. H.GongJ. H.LinD. (2025). Bayesian network structure learning by opposition-based learning. Sci. Rep.15 (1), 18447. 10.1038/s41598-025-03267-2

  • Pubmed Abstract
  • CrossRef
  • Google Scholar
  • View reference in article

• 42

  WangY.FangZ. C.HongH. Y. (2019). Comparison of convolutional neural networks for landslide susceptibility mapping in yanshan county, China. Sci. Total Environ.666, 975–993. 10.1016/j.scitotenv.2019.02.263

  • Pubmed Abstract
  • CrossRef
  • Google Scholar
  • View reference in article

• 43

  WangY.ZhouC.CaoY.MeenaS. R.FengY.WangY. (2024). Utilizing deep learning approach to develop landslide susceptibility mapping considering landslide types. Bull. Eng. Geol. Environ.83 (11), 430. 10.1007/s10064-024-03889-2

  • CrossRef
  • Google Scholar
  • View reference in article

• 44

  WangZ. F.YinC.LiJ. J. (2026). Landslide susceptibility mapping considering time-varying factors based on different models. Coatings16 (2), 207–244. 10.3390/coatings16020207

  • CrossRef
  • Google Scholar

• 45

  WuZ. N.ShenY. X.WangH. L.WuM. (2020). Urban flood disaster risk evaluation based on ontology and bayesian network. J. Hydrology583, 124596. 10.1016/j.jhydrol.2020.124596

  • CrossRef
  • Google Scholar

• 46

  WuR. Z.LiD. H.MeiH. B.LiZ. H.HuX. D.QinW. B. (2025a). Assessment of landslide susceptibility based on Bayes-optimized RUSBoost model-taking the three gorges reservoir area as an example. Bull. Eng. Geol. Environ.84 (10), 441–465. 10.1007/s10064-025-04436-3

  • CrossRef
  • Google Scholar

• 47

  WuY. Y.LiuX. L.ZhaoQ. Q.LiuH. Y.QuF.ZhangM. M. (2025b). Study on the scale dependence of the spatial distribution pattern of ecosystem service value-a case study of gangu county, China. Sci. Rep.15 (1), 18671–18689. 10.1038/s41598-025-03497-4

  • CrossRef
  • Google Scholar

• 48

  WuQ. R.XieZ.ZhaoY. F.TianM.WuY. Y.QiuQ. J.et al (2026). Integrating textual and remote sensing data for time-sensitivity of landslide susceptibility assessment. Nat. Hazards122 (2), 56–85. 10.1007/s11069-025-07831-x

  • CrossRef
  • Google Scholar

• 49

  XuS.SongY. X.LuP.MuG. Z.YangK.WangS. X. (2025). An optimized non-landslide sampling method for landslide susceptibility evaluation using machine learning models. Nat. Hazards121 (5), 5873–5900. 10.1007/s11069-024-07021-1

  • CrossRef
  • Google Scholar

• 50

  YangW. Q.GeQ. S.TaoZ. X.XuD. Y.WangY.HaoZ. X. (2026). Landslide susceptibility on the Qinghai-Tibet Plateau: key driving factors identified through machine learning. J. Geogr. Sci.36 (1), 199–218. 10.1007/s11442-026-2444-6

  • CrossRef
  • Google Scholar

• 51

  YilmazI. (2009). Landslide susceptibility mapping using frequency ratio, logistic regression, artificial neural networks and their comparison: a case study from kat landslides (Tokat-Turkey). Comput. and Geosciences35 (6), 1125–1138. 10.1016/j.cageo.2008.08.007

  • CrossRef
  • Google Scholar

• 52

  YinC.LiH. R.CheF.LiY.LiuD. (2020). Susceptibility mapping and zoning of highway landslide disasters in China. PLoS ONE15 (9), e0235780. 10.1371/journal.pone.0235780

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 53

  ZhaiS. Y.SunY.LeiJ. T.ShaoC. J. (2025). An improved information quantity method for non-landslide selection to enhance landslide susceptibility evaluation: a case study in yongfeng, south China. Nat. Hazards121 (10), 11773–11797. 10.1007/s11069-025-07261-9

  • CrossRef
  • Google Scholar

• 54

  ZhengX. S.LuW. J.JiangR. C.LiJ. H.ZhangL. M. (2025). Analysis of landslide on meizhou-dapu expressway based on satellite remote sensing. Geoenvironmental Disasters12 (1), 25–37. 10.1186/s40677-025-00331-x

  • CrossRef
  • Google Scholar

• 55

  ZhongK. H.ChenY. W.QinX. L. (2022). Sub-BN-Merge based Bayesian network structure learning algorithm. Comput. Sci.49 (S2), 64–70.

  • Google Scholar

• 56

  ZhouY. S.LiX.FaiYuenK. (2022). Holistic risk assessment of container shipping service based on bayesian network modelling. Reliab. Eng. and Syst. Saf.220, 108305. 10.1016/j.ress.2021.108305

  • CrossRef
  • Google Scholar