Nothing Special   »   [go: up one dir, main page]

Next Article in Journal
Implementing Sentinel-2 Data and Machine Learning to Detect Plant Stress in Olive Groves
Next Article in Special Issue
First Results on the Systematic Search of Land Surface Temperature Anomalies as Earthquakes Precursors
Previous Article in Journal
DOA Estimation under GNSS Spoofing Attacks Using a Coprime Array: From a Sparse Reconstruction Viewpoint
Previous Article in Special Issue
A Novel Historical Landslide Detection Approach Based on LiDAR and Lightweight Attention U-Net
You seem to have javascript disabled. Please note that many of the page functionalities won't work as expected without javascript enabled.
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Earthquake-Induced Landslide Susceptibility Assessment Using a Novel Model Based on Gradient Boosting Machine Learning and Class Balancing Methods

1
Faculty of Geosciences and Environmental Engineering, Southwest Jiaotong University, Chengdu 611756, China
2
China Railway Kunming Construction Investment Co., Ltd., Kunming 650500, China
3
Department of Railway Engineering, Sichuan College of Architectural Technology, Chengdu 610399, China
*
Author to whom correspondence should be addressed.
Remote Sens. 2022, 14(23), 5945; https://doi.org/10.3390/rs14235945
Submission received: 8 October 2022 / Revised: 15 November 2022 / Accepted: 21 November 2022 / Published: 24 November 2022
(This article belongs to the Special Issue Intelligent Perception of Geo-Hazards from Earth Observations)
Figure 1
<p>The scheme of the present study.</p> ">
Figure 2
<p>Overview map of the study area. (<b>a</b>) Geologic map at a scale of 1:200,000 overlain by the distribution of co-seismic landslides (source: geologic map in <a href="#remotesensing-14-05945-t001" class="html-table">Table 1</a>. The classification of rock groups has been modified.). (<b>b</b>) The locations and satellite imagery (source: <a href="http://goto.arcgisonline.com/maps/World_Imagery" target="_blank">http://goto.arcgisonline.com/maps/World_Imagery</a> accessed on 1 September 2022) of the study area. (<b>c</b>) Distribution of peak ground acceleration in the Hokkaido region during the 6 September 2018, earthquake (source: QuiQuake in <a href="#remotesensing-14-05945-t001" class="html-table">Table 1</a>). (<b>d</b>) Satellite imagery (source: Maxar imagery acquired on 11 September 2018) of the co-seismic landslides in the local zone in (<b>a</b>).</p> ">
Figure 3
<p>The landslide influencing factors used in this study, obtained and processed by remote sensing and GIS technology. The yellow pentagram represents the epicenter. The names of the factors in (<b>a</b>–<b>w</b>) are below the color bar.</p> ">
Figure 4
<p>Automatic extraction method for different parts of a landslide. (<b>a</b>) The division and extraction of the landslide parts. (<b>b</b>) A flowchart of our automatic extraction methodology. (<b>c</b>–<b>f</b>) The results of the landslide parts corresponding to <span class="html-italic">R<sub>d</sub></span> values of (<b>c</b>) 0.1, (<b>d</b>) 0.3, (<b>e</b>) 0.5, and (<b>f</b>) 0.9. The solid black lines represent the boundaries of the landslides. The colored cells represent the extracted parts of the landslide.</p> ">
Figure 5
<p>The influence that parameters <span class="html-italic">a</span> and <span class="html-italic">b</span> in the DCE loss function exert on landslide susceptibility prediction using the (<b>a</b>,<b>b</b>,<b>e</b>,<b>f</b>) XGB<sub>DCE</sub> and (<b>c</b>,<b>d</b>,<b>g</b>,<b>h</b>) LGB<sub>DCE</sub> models. As <span class="html-italic">R<sub>d</sub></span> increases, the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp.</p> ">
Figure 6
<p>Prediction of the susceptibilities of various landslide parts with combinations of the XGB algorithm with different class balancing methods: (<b>a</b>–<b>c</b>) XGB<sub>NB</sub>, (<b>d</b>–<b>f</b>) XGB<sub>EQS</sub>, (<b>g</b>–<b>i</b>) XGB<sub>DCE</sub>. As <span class="html-italic">R<sub>d</sub></span> increases, the part of the landslide that the model seeks to predict (the black polygon) is reduced from the entire landslide to the landslide scarp.</p> ">
Figure 7
<p>Prediction of the susceptibilities of various landslide parts with combinations of the RF algorithm with different class balancing methods: (<b>a</b>–<b>c</b>) RF<sub>NB</sub>, (<b>d</b>–<b>f</b>) RF<sub>EQS</sub>, (<b>g</b>–<b>i</b>) RF<sub>ILW</sub>. As <span class="html-italic">R<sub>d</sub></span> increases, the part of the landslide that the model seeks to predict (the black polygon) is reduced from the entire landslide to the landslide scarp.</p> ">
Figure 8
<p>The relationship between the proportion of landslides and the landslide susceptibility. From subfigures (<b>a</b>–<b>f</b>), the <span class="html-italic">R<sub>d</sub></span> increases, and the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp. Because some of the plots almost overlap (XGB-based model versus LGB-based model and EQS-based model versus ILW-based model), for greater clarity we only show representative model plots.</p> ">
Figure 9
<p>The relationship between the landslide frequency and the landslide susceptibility. The models with the plots closest to the reference lines have predicted susceptibility values that are close to the actual landslide frequency. (<b>a</b>–<b>f</b>), the <span class="html-italic">R<sub>d</sub></span> increases, and the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp. Because some of the plots almost overlap (XGB-based model versus LGB-based model and EQS-based model versus ILW-based model), for clarity we show only representative model plots.</p> ">
Figure 10
<p>The ranking of the importance of different landslide influencing factors using the (<b>a</b>) XGB<sub>NB</sub>, (<b>b</b>) XGB<sub>EQS</sub>, (<b>c</b>) XGB<sub>DCE</sub>, and (<b>d</b>) XGB<sub>ILW</sub> models. Different <span class="html-italic">R<sub>d</sub></span> (0.0–0.9) values correspond to predictions of different landslide parts in the models when landslide influencing factors are ranked by importance. Mean and STD refer to the mean and standard deviation, respectively, of the ranking of each landslide factor’s importance for all <span class="html-italic">R<sub>d</sub></span> values. Lower rankings represent more significant factor contributions. EL = elevation, diff = difference, and DIR = direction.</p> ">
Figure 11
<p>The frequencies of the various parts of the landslide with respect to the landslide influencing factors (subfigures (<b>a</b>–<b>l</b>)). In Subfigure (<b>h</b>), Un (unknown lithology), Mm (marine mudstone), Sm (siliceous mudstone), Mc (marine conglomerate), Ss (Sandstone), and Qd (Quaternary deposits) denote the lithology.</p> ">
Figure 12
<p>A combination of landslide source area prediction and landslide run-out modelling portraying the respective distribution of (<b>a</b>) the landslide impact frequency (expressed as the cumulative distribution function) and (<b>b</b>) the landslide impact probability. The landslide run-out model uses r.randomwalk [<a href="#B59-remotesensing-14-05945" class="html-bibr">59</a>], in which the motion is controlled by the landslide angle of reach. The average angle of reach (25°) and the probability distribution function used for modeling were calculated during the statistical analysis of the landslides.</p> ">
Versions Notes

Abstract

:
Predicting the susceptibility of a specific part of a landslide (SSPL) involves predicting the likelihood that the part of the landslide (e.g., the entire landslide, the source area, or the scarp) will form in a given area. When predicting SSPL, the landslide samples are far less than the non-landslide samples. This class imbalance makes it difficult to predict the SSPL. This paper proposes an advanced artificial intelligence (AI) model based on the dice-cross entropy (DCE) loss function and XGBoost (XGBDCE) or Light Gradient Boosting Machine (LGBDCE) to ameliorate the class imbalance in the SSPL prediction. We select the earthquake-induced landslides from the 2018 Hokkaido earthquake as a case study to evaluate our proposed method. First, six different datasets with 24 landslide influencing factors and 10,422 samples of a specific part of the landslides are established using remote sensing and geographic information system technologies. Then, based on each of the six datasets, four landslide susceptibility algorithms (XGB, LGB, random-forest (RF) and linear discriminant analysis (LDA)) and four class balancing methods (non-balance (NB), equal-quantity sampling (EQS), inverse landslide-frequency weighting (ILW), and DCE loss) are applied to predict the SSPL. The results show that the non-balanced method underestimates landslide susceptibility, and the ILW or EQS methods overestimate the landslide susceptibility, while the DCE loss method produces more balanced results. The prediction performance of the XGBDCE (average area under the receiver operating characteristic curve (0.970) surpasses that of RF (0.956), LGB (0.962), and LDA (0.921). Our proposed methods produce more unbiased and precise results than the existing models, and have a great potential to produce accurate general (e.g., predicting the entire landslide) and detailed (e.g., combining the prediction of the landslide source area with the landslide run-out modeling) landslide susceptibility assessments, which can be further applied to landslide hazard and risk assessments.

1. Introduction

Regional landslide susceptibility is defined as an estimation of the probability that a landslide will occur in a given area under certain environmental conditions [1,2]. A landslide susceptibility assessment provides valuable information that can be utilized in a landslide hazard and risk assessment and/or in land-use planning and traffic route decision-making in order to reduce the impacts of future landslide damage. Regional landslide susceptibility models usually focus on either predicting the susceptibility of an entire landslide [3,4,5,6] or predicting the susceptibility of one specific part of a landslide (SSPL) [7,8,9,10,11,12,13,14]. When the model focuses on predicting the SSPL, only the zones that contain that specific part (e.g., the source area or a landslide scarp) are defined as the landslide class. Compared with the prediction of the entire landslide, the prediction of SSPL can make the landslide susceptibility assessments more precise and widely applied. For example, the prediction of the landslide source area can be combined with a physical-based landslide run-out model to calculate the landslide deposition and the kinetic parameters of landslides (i.e., velocity and landslide impact frequency or probability) in a way that is more precise and consistent with physical laws. In addition, these parameters can be further used for quantitative landslide hazard and risk assessment. Therefore, it is of great significance to study the prediction of SSPL.
However, by dividing a landslide into smaller areas and only classifying one or a few of these areas as the landslide class, an extreme class imbalance arises between the number of non-landslide samples and the number of landslide samples. This class imbalance makes it more difficult to predict the SSPL than it is to predict the susceptibility of an entire landslide. In order to ameliorate this class imbalance, it is necessary to apply a class balancing method. There are three types of class balancing methods in current landslide susceptibility assessments: under-sampling techniques, over-sampling techniques, and cost-sensitive techniques. Under-sampling methods synthetically or randomly collect the same number of landslide samples and non-landslide samples using algorithms such as the equal-quantity sampling (EQS) [12,13,14], fuzzy c-means [15], easy ensemble [16], density peak clustering [17] and self-organizing-map [18] algorithms. While under-sampling significantly reduces the sample size and improves the computational efficiency, model training may be affected by the loss of data from the original samples. In over-sampling methods, algorithms such as SMOTE (synthetic minority oversampling technique) [19] and ADASYN (adaptive synthetic sampling) [20] sample or synthetically generate landslide samples to ensure that the number of landslide samples is equal to the number of non-landslide samples. While over-sampling maintains the sample balance and increases the overall sample size, it may result in higher computational costs due to the introduction of incorrect and/or redundant information into the sample data. In cost-sensitive techniques such as the direct approach, the pre-processing approach, and the post-processing approach, the class balance is maintained by explicitly accounting for the cost of model mispredictions during the training process [21]. In general, most landslide susceptibility studies employ under-sampling methods such as EQS, while over-sampling, cost-sensitive methods, inverse landslide-frequency weighting, and loss correction are hardly ever used. The EQS method sets a landslide frequency (in this article the term frequency refers to the relative frequency) of 50% in the model training and validation sets. Because the actual landslide frequency is much lower than 50%, the EQS method is not a true representation of real-world conditions, and thus yields a biased prediction of landslide susceptibility.
The susceptibility algorithm is important in predicting the SSPL, or the susceptibility of an entire landslide. Early research on this topic focused mainly on susceptibility algorithms reliant on data-driven statistical methods such as linear discriminant analysis (LDA) [2,22], logistic regression [8,23], and the weight of evidence method [24,25]. However, the landslide influencing factors used in these algorithms must be conditionally independent and linearly independent, and it is difficult to solve complicated nonlinear problems using such algorithms. In the past ten years, machine learning techniques such as decision tree (DT) [26,27], random forest (RF) [5,9,28], support vector machine (SVM) [29,30,31,32], artificial neural network (ANN) [33,34,35], and deep neural network [36] have been widely used in landslide susceptibility assessments. While these machine learning techniques can deal with nonlinear problems, they struggle with large datasets, which limits the sample size and resolution of the resulting landslide susceptibility assessment. Recently, several new breakthrough ensemble learning methods, such as the XGBOOST (XGB) [37] and LightGBM (LGB) [38] methods, have been invented. XGB and LGB can handle large datasets and have more flexible and accurate objective functions. These ensemble machine learning methods can provide new solutions for landslide susceptibility prediction.
Although a great deal of research has been carried out on prediction of the SSPL, several problems remain unsolved.
(1) The EQS class balancing method is commonly used in current research. However, this method sets a landslide frequency of 0.5 in the model training and validation sets, which is not a true representation of real-world conditions. At present, a more reasonable class balancing method has not been proposed. (2) No studies in the literature have systematically compared and evaluated the application of different class balancing methods in SSPL prediction. (3) The commonly used linear and nonlinear landslide susceptibility algorithms have many limitations. Despite the numerous advantages of using the XGB and LGB algorithms, no studies have employed these frameworks in prediction of the SSPL.
The goals of this study were to develop a landslide susceptibility model that combines the XGB, LGB, and dice cross-entropy (DCE) loss function to solve the class imbalance problem and achieve better model performance in SSPL prediction. We explored how the choice of class balancing methods (EQS, inverse landslide-frequency weighting (ILW), DCE loss, or non-balance (NB) methods) and susceptibility algorithms (XGB, LGB, RF, or LDA) affects the prediction of the SSPL. In this study, the landslides triggered by the 2018 Hokkaido earthquake (MW 6.6) were taken as a case study, and the prediction of the susceptibility of six different landslide parts was considered.

2. Materials and Methods

The main framework of this article is presented in Figure 1. In Part I, we describe our study area and construct our dataset. Our high-quality dataset consists of 10,422 earthquake-induced landslides in the 2018 Hokkaido earthquake and 24 landslide influencing factors. We discuss our methodology in Parts II and III. Our technique automatically extracts six different parts of the landslide (e.g., the entire landslide, the landslide source area, and the landslide scarp) and compiles separate inventories for them. Each type of landslide inventory is combined with 24 landslide factors and then separated into six datasets. We created 13 landslide susceptibility models using combinations of four susceptibility algorithms (XGB, LGB, linear discriminant analysis (LDA), and random forest (RF)) and four class balancing methods (non-balance (NB), equal-quantity sampling (EQS), inverse landslide-frequency weighting (ILW), and DCE loss). We input all six datasets into each of the 13 models and generated 78 sets of results. In each case, the modified group cross-validation method (MGCV) was used for the model training, prediction, and validation. Part IV includes the Results and Discussion sections. We used machine learning metrics and landslide statistics to evaluate and analyze the results of the landslide susceptibility.

2.1. Study Area

Our study area (Figure 2) is located in the eastern Iburi area of central Hokkaido (bounded by latitudes of 42.49°N and 43.02°N and longitudes of 141.69°E and 142.36°E; Figure 2a). The terrain in this area is low and hilly, with flatter profiles in the west and steeper topography in the east. The altitude ranges from 1 m to 873 m and the slope degree ranges from 0.1° to 57.7°. The minimum and maximum monthly average temperatures are −7 °C and 20.4 °C, respectively. Most (42%) of the yearly precipitation falls between July and September, with a monthly average precipitation of 142 mm.
The lithology of our study area (Figure 2a) is predominantly (63.4%) marine sedimentary rocks such as marine siliceous mudstones, marine conglomerates, marine mudstones, and marine sandstones. The igneous rocks (5.2%) in our study area mainly consist of rhyolite, andesite, basalt, and diorite, while metasandstone chlorite, pelitic schist, and mafic schist are the main metamorphic rocks (1.1%) found in this area. The remaining strata are mainly Quaternary deposits. In terms of its structural features, our study area is located in an active compression collision zone created by the tectonic action of the Amur Plate, the Philippine Sea Plate, the Okhotsk Plate, and the Pacific Plate. Per the USGS seismic database, a MW 6.6 earthquake (epicenter location of 42°41′10″N 141°55′44″E and focal depth of 35 km) struck central Hokkaido on 6 September 2018. According to our landslide inventory, this earthquake triggered at least 10,422 landslides.

2.2. Basic Data

Before carrying out a landslide susceptibility assessment, it is necessary to prepare the basic data such as the landslide inventory and the landslide influencing factors. These data were gained and processed through the remote sensing technology, geographic information system (GIS) technology, and field surveys. The dataset for preparing the landslide inventory and the landslide influencing factors are presented in Table 1. By comprehensively considering the data quality, sample size, and calculation cost, we selected the grid cells (referred to as pixels) as the mapping units and resampled all of the basic data to a resolution of 30 m.
We prepared a landslide inventory containing 10,422 co-seismic landslides from the Hokkaido earthquake on 6 September 2018. Our landslide inventory was created based on the landslide data (3307 landslides) obtained from the Geospatial Information Authority of Japan and the landslide inventory of [39] (5625 landslides, covering an area of 46.3 km2). We fully expanded and revised the original landslide inventories through visual interpretation of high-resolution Google Earth images and UAV photos taken before and after the earthquake (Table 1) and through digitization. We interpreted the landslides by starting at the area with the highest seismic intensity and working our way outwards until we found no further landslides; the total area of the interpretation was greater than 2000 km2. The interpretation and digitization of each landslide covered the range from the source area to the deposition area (Figure 2d).
The landslides in this inventory have a total area of 49.97 km2, a total volume of 62 × 106 m3 (based on the volume formula for soil landslides proposed by Larsen, et al. [40]), and a density of 5.2 landslides/km2. The quartiles and the average of the area of the individual landslides are 1378 m2, 2763 m2, 5575 m2, 170,255 m2, and 4795 m2. The landslides mainly occurred in the Quaternary deposits and the mudstone layer, and can be classified as flow-slides, debris avalanches, debris slides, and rock slides. The study area had experienced prolonged rainfall and strong earthquakes, and numerous landslides liquefied as a result.
In our landslide susceptibility assessment, we focused on 24 landslide influencing factors (Table 2 and Figure 3), including the seismology, geomorphology, hydrology, geology, vegetation, and rainfall of our study area. PGA, PGVA, and distancefocus are directly related to both the amount of ground shaking and earthquake intensity. The distanceridge factor provides insight into the seismic site effects. Sites that are located close to mountain ridges generally experience seismic amplification, as after propagating through a mountain, seismic waves are reflected, diffracted, and channeled toward the mountain ridge, where the waves interfere constructively with one another [41]. By summarizing the relationship between the direction of seismic wave propagation and the slope attitude [42], factors such as the epicentral direction, angleES, angleFH, and angleFS influence the occurrence of landslides.
In terms of topography and geomorphology, slopeMAX and the elevation difference are used to describe the steepness of the terrain. Landslides rarely develop in gentle terrain because it is difficult to generate effective free faces and the potential energy required to trigger a landslide. The elevation ratio describes the position of a point with respect to the slope (e.g., at the bottom of the slope, on the shoulder of the slope, etc.). Differences in the slope stress state and the thickness of the weathering layer or unloading layer affect the slope stability as well. In terms of the hydrological conditions, we considered river factors such as distanceSR, distanceMR, and SPI [44]. We accounted for these variables because a river can alter the stress state (e.g., by undercutting unloading or via the seepage force) and/or the physical properties (e.g., by softening the rock or soil) of a slope. According to the influence a river of a given size exerts on a slope, that river is classified as either a major river or a minor river according to the Strahler stream order. The meteorological factors include the cumulative precipitation in the four days preceding the earthquake. Because the earthquake occurred during the rainy season, we inferred that abundant rainfall facilitated the onset of multiple landslides during the earthquake by saturating and softening the slope rocks and soils.
Other geological factors include the lithological and structural features of a given area. The lithology affects the material composition, strength, structure, and weathering resistance of the slope. For example, mudstone is characterized by low strength, low structural integrity, and weak intercalation. Structural factors include distancefault, distancefold, and fault density. Slopes with faults and folds often develop additional joints, fissures, weathered zones, and groundwater paths; naturally, these features increase the slope instability.

2.3. Methodology

2.3.1. Automatic Extraction of Different Parts of Landslides

In order to establish the landslide susceptibility model, we developed an algorithm that automatically and efficiently extracts the different parts of the landslide. This algorithm assumes that each part of a landslide has a different relative elevation [45]. As is shown in Figure 4a, a landslide can be divided into two parts at a certain elevation Em, namely, the upper part (i.e., landslide scarp or landslide source area) as the landslide class and the lower part (i.e., landslide accumulation area) as the non-landslide class; this allows us to specifically extract the areas where the landslide initiates, such as the source area. To facilitate this process, we can define the relative elevation ratio Rd (0–100%) as follows:
R d = E m E min E max E min = H d H
where Emax and Emin (m) are the highest and lowest elevation of the landslide, respectively, and Hd and H (m) represent the height difference of the accumulation zone and the height of the entire landslide, respectively. In order to accurately represent the landslide feature, we chose to define the landslide parts as different sections of a surface rather than as a single point (e.g., the centroid). The grid cells are used as the mapping units, with the landslide class and non-landslide class cells assigned values of one and zero, respectively.
Using our algorithm (Figure 4b), we were able to automatically extract and identify 62,532 parts representing six different landslide features (Figure 4c–f). Instead of extracting the landslide accumulation area separately, we extracted the area containing the landslide source area for landslide susceptibility prediction, as the landslide accumulation area is not the location where the landslide initiates. We compared and analyzed the extraction results with satellite images of the landslides. An Rd value of 0% corresponds to an entire landslide, from the source area to the accumulation area. An Rd of 10–30% encompasses the parts of the landslide that are not the accumulation area, such as the source and the transport area of a rock avalanche. An Rd value of 50–70% corresponds to the source area. An Rd value of 90% corresponds to the landslide scarp. In the automatic landslide extraction algorithm, the accuracy of the landslide polygon boundaries is related to the output resolution, while the accuracy of landslide cell extraction and identification is affected by the DEM quality. A statistical analysis of each part of the landslide is shown in Table 3. Here, we use realistic landslide frequency values, rather than the ideal value of 0.5 that has often been used in previous studies.

2.3.2. Validation Method

In this study, we devised and implemented the modified group cross-validation method (MGCV) for validation of the landslide susceptibility models. In this method, we treat each landslide or landslide part as a whole; in this way, we ensure that the cells in a given landslide can only be assigned to the training set or to the validation set. The process of the MGCV is as follows: (1) MGCV randomly divides all of the landslides (or landslide parts) in the inventory into KG groups (KG = 4 in this study), where each independent group is composed of the same number of samples. It then randomly selects non-landslide cells and distributes them evenly amongst the four groups. This ensures that the training set, validation set, and total sample set have the same realistic landslide frequency. (2) Then, taking the i-th group of samples as the model validation set, the remaining KG − 1 groups are combined into a model training set, the model landslide susceptibility prediction is compared to the observed landslides of the validation set, and the evaluation metric is calculated. (3) Finally, step 2 is repeated KG times to complete the cross-validation process, with each group of samples used as the validation set exactly once. The model performance is the average of the metric scores calculated in KG iterations.

2.3.3. Evaluation Metrics

The receiver operating characteristic curve (ROC) [46] is an evaluation metric frequently used in landslide susceptibility assessments. The area under the ROC (AUC) is used as a quantitative indicator of model performance; a higher AUC corresponds to a better model performance. This study employs the ROC as a baseline reference for comparison with other studies, not as a core indicator of model performance. In fact, the ROC only measures the relative ranks of landslide susceptibility, rather than the size of the landslide susceptibility value.
The Jaccard index (Jaccard), which ranges from 0 to 1, measures the similarity between the landslide susceptibility and the presence or absence of a landslide. Larger Jaccard values correspond to better model performance. In this paper, Jaccard is the geometric mean of JaccardP and JaccardN, which measure the similarity in negative samples (i.e., the non-landslide cells) and positive samples (i.e., the landslide cells), respectively:
JaccardP = TP / ( TP + FP + FN )
JaccardN = TN / ( TN + FN + FP )
By equally considering the contributions of both the positive and negative samples to the final prediction result, the Jaccard value avoids placing too much influence on classes that have proportionally more samples in a given dataset. More traditional metrics (i.e., accuracy score) can be biased by such classes. The F1 score is ideal for tasks in which a class imbalance may arise, as it considers both the recall and precision:
F 1 = 2 Precision × Recall Precision + Recall
where Precision = TP/(TP + FP) and Recall = TP/(TP + FN). The true positive (TP), false positive (FP), true negative (TN), and false negative (FN) are calculated using the confusion matrix approach [46], where the terms positive and negative refer to the presence or absence of landslides, respectively.

2.3.4. Landslide Susceptibility Algorithms

XGBoost (XGB) [37] is a gradient boosting framework that supports GPU and parallel computing. Because XGB is flexible, efficient, and accurate, it can handle large datasets such as those found in landslide susceptibility assessments. XGB is a tree-based ensemble machine learning approach. The generation of an XGB model begins with a single classification and regression tree (CART); in each model iteration, a new CART is created using the current error, then that CART is then added to the output of the existing model. The XGB model grows iteratively in this fashion until the convergence criteria are satisfied. The output of an XGB model is the combination of the outputs of K CARTs:
y ^ i = ϕ ( x i ) = t = 1 K f t ( x i ) ,
where i represents the i-th sample, xi is the landslide influencing factors, y ^ i is the landslide susceptibility prediction, K is the number of iterations, and ft(xi) is the landslide susceptibility prediction generated by CART Tt added in the t-th iteration. Here, Tt is determined mainly by minimizing L(t), which is the objective function
L ( t ) = i = 1 n l ( y i , y ^ i ( t ) ) + R ( f t ) ,
where y ^ i ( t ) = k = 1 t f k ( x i ) = y ^ i ( t 1 ) + f t ( x i ) is the output of the model in the t-th iteration; yi represents whether landslides occur (1) or do not occur (0), l is the loss function, and R(ft) is the regularization term of the model. By further expanding i = 1 n l ( y i , y ^ i ( t ) ) using the second-order Taylor expansion, we account for the structure of the CART in the objective function and rewrite the objective function as [37]:
L ( t ) = 1 2 j = 1 V ( i I j g i ) 2 i I j h i + α + β V
Where g i = l ( y i , y ^ i ( t 1 ) ) / y ^ i ( t 1 ) and h i = 2 l ( y i , y ^ i ( t 1 ) ) / ( y ^ i ( t 1 ) ) 2 , V is the total number of leaves in Tt, Ij is a sample group belonging to leaf node j, and α and β are regularization coefficients. The optimal structure of Tt and ft(xi) can be determined based on Equation (7) and the specific growth mode of the CARTs. This y ^ i ( t ) calculation completes the tth iteration.
LightGBM (LGB) [38] is a gradient boosting framework that is similar to XGB. LGB was specifically designed to expedite the slower parts of XGB. Gradient-based one-sided sampling reduces the training cost while preserving the accuracy of the model. The histogram algorithm divides the continuous landslide influencing factors into discrete bins and uses those bins to create new CARTs, which greatly reduces the costs of computation and storage in model training. Exclusive feature bundling binds mutually exclusive landslide factors to reduce the factors involved in model training, which significantly improves training speed without affecting model accuracy.
Linear discriminant analysis (LDA) and random forest (RF) have been widely used in studies of debris flows and landslides. Detailed descriptions of LDA and RF can be found in [47] and the [48], respectively. None of the aforementioned algorithms considers hyperparameter optimization, because LDA, RF, XGB, and LGB all have good initial hyperparameters. The main hyperparameters are described as follows. In LDA, the solver is singular value decomposition and the threshold for a singular value is 1 × 10−4. In RF, the number of trees is 100, the split criterion is Gini impurity, and the depth of the tree is unlimited. In XGB, the number of trees is 100, the learning rate is 0.3, the subsample ratio is 1, and the depth of the tree is 6. In LGB, the number of trees is 100, the learning rate is 0.1, the subsample ratio is 1, and the depth of the tree is unlimited. The above algorithms and the DCE loss function (described below) were coded and calculated using Python’s scientific computing libraries (i.e., NumPy, XGBoost and LightGBM).

2.3.5. DCE Loss Function

This study proposed two new XGB and LGB frameworks that employ the DCE loss function. We synthesize the DCE loss using the dice loss [49] and weighted cross-entropy loss algorithms, where the dice loss is actually equivalent to the F1 score. Let l ( y , y ^ ) be defined as an n-dimensional vector l ( y , y ^ ) = [ l 1 , l 2 , , l n ] where li is the DCE loss of the i-th sample. The relevant formulas are as follows:
l i ( y i , y ^ i ) = b l C i + ( 1 b ) l D i
l C i = [ ( 1 a ) y i log p i + a ( 1 y i ) log ( 1 p i ) ]
l D i = n l ( 1 2 < p , y > + ε p i + y i + ε )
where l C i and l D i are the weighted cross-entropy loss and dice loss of the i-th sample, respectively, and a and (1 − a) are the weights of the non-landslide and landslide samples in l C i respectively. It is possible to address the class imbalance by using a smaller a weight in the training process; this choice reduces the contribution of non-landslide samples and increases the contribution of landslide samples. The parameters b and (1 − b) are the weights of l C i and l D i in l i . The parameters a and b must be tuned in order to identify the optimal value. The n l parameter is the coefficient that adjusts the relative sizes of l C i and l D i , while pi represents the landslide susceptibility prediction, where p i = 1 / ( 1 + e y ^ i ) . The ε parameter is the smoothing coefficient, which we set to 1. Parameters p, y ^ , and y are n-dimensional vectors, p = [ p 1 , p 2 , , p n ] , y = [ y 1 , y 2 , , y n ] , and y ^ = [ y ^ 1 , y ^ 2 , , y ^ n ] . The first- and second-order gradients of l ( y , y ^ ) with respect to y ^ are the Jacobi and Hessian matrices, respectively. From the diagonal elements of these matrices, we can identify the first-order gradient and its ith term:
y ^ l ( y , y ^ ) = [ l 1 y ^ 1 , l 2 y ^ 2 , , l n y ^ n ]
where
l i y ^ i = y i ( p i + a 1 ) + a p i ( 1 2 y i ) + n p i ( 1 p i ) s t 2 s = 2 < p , y > + 1 2 y i t t = p i + y i + 1
In addition, we can identify the second-order gradient and its i-th term:
y ^ 2 l ( y , y ^ ) = [ 2 l 1 y ^ 1 y ^ 1 , 2 l 2 y ^ 2 y ^ 2 , , 2 l 1 y ^ n y ^ n ]
2 l i y ^ i y ^ i = p i ( 1 p i ) { a ( 1 2 y i ) + y i + n s [ ( 1 2 p i ) t 2 p i ( 1 p i ) ] t 3 }
The results of Equations (12) and (14) are substituted into the gi and hi terms of Equation (7), respectively, for the t-th iteration of model training.

2.3.6. Class Balancing Method

We explored four class balancing methods: the equal quantity sampling (EQS), inverse landslide frequency weighting (ILW), loss correction, and non-balance (NB) methods. EQS creates a class-balanced dataset by combining the landslide samples and an equal number of non-landslide samples obtained through random under-sampling. The landslide frequency of the dataset produced using this method is 0.5. In order to ensure that the landslide frequency in the validation set is realistic (rather than 0.5), the EQS in this study was used to generate only the training sets, instead of both training and validation sets as in previous studies. ILW addresses the class imbalance by using the landslide prior probabilities of the samples to amplify the weight of the landslide samples in the model training process. The weights of the landslide and non-landslide samples are set to (Nl + Nf)/Nl and (Nl + Nf)/Nf, respectively, where Nl and Nf are the number of landslide samples and non-landslide samples in the training set, respectively.
When using the loss correction method, loss functions are created or corrected in order to improve the model training when datasets are class imbalanced. In this study, we employed the DCE loss correction (described in Section 2.3.6), which is specifically applicable to the XGB and LGB. For comparison, we considered a model that does not include a class balancing method, which we name the non-balance (NB) model.
These four methods for balancing classes were used in the training of the MGCV method. The four landslide susceptibility algorithms were combined with four class balancing methods to create 13 landslide susceptibility models. We refer to a given landslide susceptibility model using nomenclature based on the algorithm and class balancing method used for that model; for example, XGBNB represents the model combining XGBoost (the susceptibility algorithm) and NB (the class balancing method).

3. Results

3.1. Model Performance Evaluation

As discussed previously, we used 13 different models to predict the landslide susceptibility of various landslide parts. As is shown in Table 4, for every model, increasing Rd values correspond to decreasing PRE, REC, F1, JACP, JAC, and Jaccard values. This means that when landslide parts are reduced in size, the prediction performance of all models deteriorates. The AUC only decreased by up to 5% (i.e., RFILW) with increasing Rd; this means that this metric is insensitive to changes in the landslide part.
The evaluation metric scores for each model varied considerably with the class balancing method, especially at higher values of Rd. We grouped the models according to their balancing methods and calculated the average values of the evaluation metric scores (considering different Rd). The EQS (or ILW) method yielded the highest recall (mean recall: 0.92) and the lowest precision (mean precision: 0.21). The NB method produced the highest precision (0.52) and the lowest recall (0.18). The precision (0.38) and recall (0.60) values for the DCE loss (XGBDCE and LGBDCE) models were more balanced.
The models, again classified by their class balancing methods, were comprehensively evaluated based on the F1 and Jaccard metrics. The DCE loss method produced the best average F1 (0.465) and Jaccard (0.540) scores, while the NB method produced the lowest average F1 (0.258) and Jaccard (0.365) scores. The EQS and ILW methods produced F1 (0.311 and 0.293) and Jaccard (0.406 and 0.381) scores that fell in between these two extremes. We conducted statistical analyses on the model results based on their primary algorithms. The average F1 scores for the algorithms, from best to worst, are XGB (0.369), LGB (0.337), RF (0.304), and LDA (0.163) (Table 4).
As mentioned previously, a and b in the XGBDCE and LGBDCE models represent specific model parameters. Here, we investigate how these values affect the models themselves and their scoring metrics (Figure 5). As is shown in Figure 5a,e, the highest AUC and F1 scores occurred when a = 0.2. A value of index a that is too high or too low proves detrimental to the performance of the model. When b = 0, the DCE loss function reduces into the dice loss function, resulting in very low F1 and AUC values. When b > 0, the AUC metric remains largely unchanged (Figure 5f). When b = 1, the DCE loss function reduces into the weighted cross-entropy loss function, causing a sharp reduction in the F1 score. This suggests that the model does not perform as well as the DCE loss function does when either the Dice loss or cross-entropy loss is used separately as an individual loss function. The best AUC and F1 scores correspond to a b value of 0.05. The LGBDCE model produced F1 and AUC trends (Figure 5c,d,g,h) that are similar to those in XGBDCE.

3.2. Landslide Susceptibility Mapping

Using the MGCV method, we used four iterations of each model to generate cross predictions of landslide susceptibility for our entire study area (Figure 6 and Figure 7). Because the results produced by the LGB and XGB algorithms are relatively similar, we chose to highlight the XGB and RF algorithm results in our figures. The statistical summary of the cells in each model (Figure 6 and Figure 7) are summarized in Table 5.
The models using either the EQS or ILW class balancing methods included the largest number of cells in the high susceptibility class. However, the NB method models had the largest number of cells with a susceptibility of 0.0–0.1 (>92% of the total), indicating that the landslide susceptibility determined using the NB method is low overall. As the Rd value increase from 0.0 to 0.9, the high-susceptibility area (0.6–1.0) predicted by the XGBDCE model changed from the entire landslide area (Figure 6g) to the landslide source area (Figure 6h) and the landslide scarp (Figure 6i), demonstrating that the XGBDCE model is able to precisely adjust the spatial distribution of the landslide susceptibility according to the target landslide part. With the EQS balancing method, when Rd increases from 0.0 to 0.9, the number of cells in the high-susceptibility area (0.6–1.0) only decreases by 10% (Table 5). Figure 6d–f and Figure 7d–f show that the high-susceptibility areas on the hillside are much larger than the actual landslide areas; therefore, the EQS method cannot precisely predict a specific part of the landslides. For the NB method, Figure 6a–c and Figure 7a–c show that the XGBNB and RFNB models are able to accurately adjust the spatial distribution of the landslide susceptibility according to the target landslide part. However, in these models, the susceptibility values within the dense landslide area were low.
We quantified the relationship between the proportion of landslides and the susceptibility on the landslide susceptibility maps (Figure 8). In the NB-based model, the proportion of landslides is negatively correlated with landslide susceptibility. However, in the EQS-based and DCE-based models the proportion of landslides increases with the landslide susceptibility, indicating that the landslide proportion in these models is more reasonable than that in the NB model. The proportion of landslides varies with increasing Rd. In the NB-based model, the proportion of landslides in the 0.0–0.2 susceptibility class increased very rapidly (average value of 55%). In the XGBDCE models, only when Rd reaches 0.9 does the proportion of landslides decrease as the susceptibility value increases. In the EQS-based models, the proportion of landslides changed relatively little. The above results indicate that the proportion of landslides in EQS-based and DCE-based models is more stable than that in NB-based models when the target landslide part changes.
Our findings on the relationship between the landslide frequency and the susceptibility on the landslide susceptibility maps are summarized in Figure 9. The landslide frequency differs significantly amongst the various models. The landslide frequency in the NB-based models is closest to the reference line (the actual landslide frequency), while that in the EQS-based models is very low and furthest away from the reference line. The landslide frequency of the DCE-based models falls in between the above two extremes.
As is shown in Figure 9, the landslide frequency decreased as Rd increased in all of the models. In the EQS model, the largest landslide frequency reduction occurred in the 0.6–1.0 susceptibility class (on average 6.5 times lower). In this susceptibility class, the model predicted that many landslide cells were actually non-landslide cells (i.e., false positives). For the NB-based model, when Rd ≥ 0.7, the landslide frequency in the high-susceptibility class is even lower than that in the low susceptibility class. The frequency of landslides in the XGBDCE model decreased by a factor of only 1.5 on average. Therefore, the model was maximally stable when there were changes to the target landslide parts.

3.3. The Analysis of Landslide Influencing Factors

Next, we explored how various landslide influencing factors affect the prediction of different landslide parts in these models. We quantified the influence exerted by each factor using the reduction of the loss function in the XGB-based model (Figure 10). The rank, Mean, and STD of the landslide influencing factors vary greatly depending on the model (Figure 10). Because the XGBDCE model generated the most accurate predictions, we focus on the results for this model here. According to the ranked mean of the landslide influencing factors (Figure 10c), the PGA, curvatureSTD, angleFS, angleFH, precipitation, distancefocus, and lithology are the primary factors controlling landslide occurrence. It is worth noting that earthquakes and rainfall events are factors that can trigger landslides. As shown in Figure 11a, there is a significant positive correlation between the PGA and landslide frequency. However, contradicting known patterns, more precipitation corresponds to a lower landslide frequency, that is, the maximum landslide frequency (0.07) was observed in the low rainfall (5–17 mm) area (Figure 11f). Clearly, the landslides in this study area were more heavily influenced by ground motion than by the amount of rainfall.
Figure 10c shows that certain factors, such as PGA, curvatureSTD, angleFS, angleFH, precipitation, and distancefocus, have relatively low-ranking standard deviations. Additionally, because there is no comparative change in the landslide frequency for a given factor (Figure 11) for increasing values of Rd, we infer that these factors are not sensitive to changes in the parts of the landslide. Conversely, the ranked standard deviations of factors such as slopeMAX, distanceSR, distanceridge, elevation ratio, and stream power index (SPI) are significantly larger (Figure 10c). The relationships between these factors and landslide frequencies vary for the different parts of the landslide. For higher values of Rd, the landslide frequency increases with the elevation ratio and decreases with the SPI and distanceridge metrics (Figure 11i,k,l). Therefore, the source areas of the landslides are most likely to be found on the upper slopes of the hills, closer to the ridges. The opposite trend is observed for smaller values of Rd (Figure 11).
The mechanism of the earthquake-induced landslides in this study area was analyzed based on the main landslide influencing factors, as follows. The type of the soil and rock determines the strength of the slope. The landslide frequency of siliceous mudstone and marine conglomerate is significantly higher than that of other lithologies (Figure 11h), mainly because these two lithologies have low strength and weak regoliths. The slip surfaces of most landslides in the study area are shallow and nearly planar; therefore, the slip surface should mainly follow the interface of the bedrock and the overburden or interface between different layers in the overburden (mostly composed of the residual soil, unweathered volcanic ash, weathered pyroclastic soil, and pumice [50]). These weak interfaces are low in yield acceleration, and are more likely to deform under the seismic load. Precipitation plays a key role in landslide occurrence (Figure 10c). Prolonged rainfall before the earthquake (the cumulative precipitation in August and September was 200–300 mm) increased the water content of the soil and rock, softening them.
The seismic factors are the most important landslide influencing factors (Figure 10c). Under strong ground motion, the pore pressure of the soil increased rapidly (most of the soils in the study area are fine grained with impeded drainage), leading to a significant decrease in the undrained shear strength (and the yield acceleration) of the soil [51], with soil layers even being liquefied. Each time the ground acceleration exceeded the yield acceleration, the slope mass slid along the aforementioned weak interfaces. In the process, the joints and cracks of the slope propagated, and the deformation modulus and the strength of the slope decreased [52,53]. When the permanent displacement of the slope reaches a certain threshold (with higher PGA values associated with larger permanent displacement [54]), landslides are very likely to occur. Most landslides in the study area initiated and slid along the valley channel for a relatively long distance. The slip surfaces on the slope are smooth, without residual sliding mass [39]. These phenomena indicate that the water content of the sliding mass was relatively high.

4. Discussion

4.1. The Applicability of the Automatic Extraction Method

We devised an automatic extraction method to efficiently identify and extract the various landslide parts (Figure 4). This extraction method only requires DEM data and the polygonal landslide inventory. Using our technique, an ordinary computer can extract 360,000 landslide parts in 6 s. This automatic method is applicable to both the new landslide inventories and older landslide inventories that are missing data or are use low quality data. In contrast, manual extraction methods require experienced specialists, high-quality imagery, and DEMs. For large study areas, it may take months to prepare the data using manual extraction methods. While automatic landslide extraction methods that are based on relative elevation ratios (Rd) are empirical, their precision is sufficient for regional landslide studies, in which the resolution is usually >10 m. Within this resolution, landslides are represented by a small number of cells. Even when more precisely determining the boundaries of each part of the landslides using manual methods and then extracting these parts, the extraction results for most of the landslides change only slightly. Due to the inherent subjectivity and additional constraints present in the manual extraction method, automatic extraction methods should always match the accuracy of manual extraction methods at the regional scale.

4.2. Comparison and Prospect of Landslide Learning Algorithms

We mainly measured landslide susceptibility model performance using F1 score instead of the commonly used AUC of the ROC, as the number of landslide and non-landslide samples in our dataset was extremely imbalanced. In such cases, the F1 score is far more stringent than the AUC of the ROC. Therefore, in this study the F1 score is far less than the AUC of the ROC (Table 4), which is greater than 0.9 (indicates an excellent model performance) [55] for all of the models.
In order to enhance the standards of landslide susceptibility evaluation, the landslide susceptibility algorithm should meet multiple criteria. Compared with more traditional methods, the XGB and LGB algorithms all enable the user to adjust the loss functions based on the specific application. In this study, we devise the XGBDCE and LGBDCE models, which are based on DCE loss algorithm, to address the extreme class imbalance that is inherent in landslide susceptibility modeling. On average, the F1 values of these models are 59% (XGBDCE) and 49% (LGBDCE) higher than that of the RFEQS (the optimum RF model). In terms of the computational efficiency, we experimented with RF and ANN as well. The computational costs of the RF and ANN algorithms are 66 and 114 times higher, respectively, than that of the XGB and LGB models. An increasing demand for large datasets is anticipated in the future. For example, it is necessary to analyze landslide susceptibility using rapidly updated data (i.e., meteorological data) or high-resolution data at larger regional scales. Because the XGB and LGB algorithms are more accurate in prediction and more efficient in computation, they have wide applicability in future landslide susceptibility assessment.

4.3. Influence and Applicability of Class Balancing Methods

We explored how landslide susceptibility prediction is affected by different class balancing methods. While the NB-based models, which did not utilize any class balancing method, yielded the highest precision, they had the lowest recall ratio (Table 4). Most of the landslides predicted by these NB models (approximately 76%) are located in the lowest susceptibility area (0.0–0.2) (Figure 8). As such, we conclude that these NB models underestimate landslide susceptibility. This underestimation is due to most of the samples being in the non-landslide class; however, the weights of the landslide and non-landslide class are the same in the objective function of NB method, and therefore predicting the samples as the non-landslide class (i.e., with a low landslide susceptibility) is beneficial to the optimization of the objective function of the NB method. Therefore, the NB model is suitable for small-scale earthquake events with low landslide frequencies, in which the underestimation of the landslide susceptibility does not result in serious consequences. Moreover, the more precise values found in these NB models can accurately predict small potential landslide areas. However, when landslides are triggered by large earthquakes (e.g., the Wenchuan earthquake [56] and Nepal earthquake [57]), the NB-based models may fail to predict most landslides.
The EQS and ILW class balancing methods, on the other hand, have the highest recall values and the lowest precision values (Table 4). The landslide frequency of the high landslide susceptibility area (0.8–1.0) in these models is very low. Therefore, we conclude that these methods result in significant overestimation of the landslide susceptibility. This overestimation occurs because these techniques involve an inherent amplification of the landslide frequency or the landslide-class weight in the training dataset. For example, when using the EQS-based model, the landslide frequency in the training dataset was set to 0.5, while the landslide frequency of historically strong earthquakes such as the 2018 Hokkaido earthquake (Mw 6.6), the 2008 Wenchuan earthquake (Mw 7.9) [56], and the 2015 Nepal earthquake (Mw 7.8) [57] were only 0.047, 0.12, and 0.01, respectively. Therefore, we conclude that it is not appropriate to use the EQS and ILW class balancing methods in models used to predict landslide susceptibility in areas affected by small or medium scale earthquake events. Due to overestimation, regardless of which part of the landslide is used to create the model, these models can only predict the overall landslide susceptibility of entire hillsides, and are unable to predict the susceptibility of smaller landslide features (Figure 6 and Figure 7).
The prediction performance of the DCE loss-based models is far better than that of the other class balancing methods (Table 4). In addition, by using a specific weighted F1 score as the optimization objective, the DCE loss method achieves more balanced precision and recall, thereby solving the problem of overestimating or underestimating the landslide susceptibility that occurs in other class balancing methods. In order to implement the DCE loss function, the user must determine the values of parameters a and b. Although the maximum number of non-landslide cells is 136 times that of landslide cells (Figure 5), the optimal weighted ratio of positive to negative is only four. Therefore, the a value does not need to refer to the prior landslide probability (i.e., 1/136 when Rd = 0.9) of an event. Under different Rd values, the pattern of the F1 score and AUC curve is consistent (Figure 5). These trends demonstrate that reasonable values of a and b can be applied to a dataset of different parts of a landslide as well as a dataset of different frequencies (e.g., earthquake events of different intensities).
Model performance decreased significantly when shifting the focus from predicting the susceptibility of an entire landslide to predicting the susceptibility of a specific landslide feature (Table 4). The main reason for these results is that predicting low-frequency events (a specific landslide feature) is more difficult than predicting high-frequency events (an entire landslide). Despite this, studies have found [13,23] landslide susceptibility models that perform better when focused solely on landslide scarps. These contradictory results can be explained by the fact that the frequencies of the entire landslide and the landslide scarp in these studies [13,58] were both set to 0.5, a value that is dozens of times larger than the actual landslide frequency. Therefore, for that specific study, it became easier to predict the susceptibility of landslide scarps and the prediction performance was significantly overestimated.
While standard models can only be applied to predict the susceptibility of entire landslides, the XGBDCE and LGBDCE algorithms can be applied to predict the susceptibility of smaller landslide features. The susceptibilities predicted for an entire landslide by the XGBDCE and LGBDCE models (Figure 6g) can be applied across large geographic regions. With regard to predicting a specific part of a landslide, such as the landslide source area (Figure 6h,i), these prediction results can be combined with a physical-based landslide run-out model to calculate the landslide deposition area and the kinetic parameters of landslides in a way that is more precise and consistent with physical laws than the traditional methods (Figure 12). Because these algorithms are both flexible and accurate, they can be applied to large-scale regions to increase the precision of landslide hazard or risk assessments. As shown in Figure 12a,b, we calculated the landslide impact frequency and impact probability parameters for landslide hazard and risk assessment, respectively. The spatial distribution of these parameters is in agreement with the observed landslides.

4.4. Limitations

First, the predicted susceptibility of landslide scarps (Rd = 0.9) was not ideal. An attempt to improve the model performance which involved the addition of 30 landslide influencing factors related to the microtopographical features and earthquakes was unsuccessful. There is a relationship between subpar model performance and the lack of high-quality data. At a resolution of 30 m, it is difficult to gather high-quality data that exactly describe the characteristics of a landslide scarp. However, data with resolution higher than 30m were unavailable in this study area. Earthquakes and rainfall have been proven to be the two most significant factors driving the triggering of landslides (Figure 10c). However, the available seismic data and rainfall data do not account for the localized characteristics of the rainfall or earthquakes (i.e., the site effects).
Second, certain factors that can improve the performance of landslide susceptibility models were not considered in the analysis of the earthquake-induced landslide mechanism in Section 3.3. (1) The undrained strength of the overburden and of the interface between the bedrock and overburden. Most of the landslides in the study area are shallow translational landslides, and therefore these strength factors are helpful in predicting this type of landslide. (2) The parameters that reflects the physical properties of the soil and rock under the rainfall conditions (e.g., degree of saturation, relative density, coefficient of permeability). The occurrence of earthquake-induced landslides in the study area is related to the prolonged rainfall before the earthquake. These factors could be considered in the model to better measure the effect of rainfall. (3) The profile characteristics of the soil and rock, such as the type, depth, and sequence. These factors influence the characteristics of the slip surface of landslides and the site effect of ground motions. However, at present, data for the above factors in the study area are unavailable. Due to the large spatial variability of these geotechnical factors [32], it is necessary to carry out a large number of field investigations and experiments in the future.
Third, our method is not able to predict landslide occurrence prior to extreme triggering events such as strong earthquakes or storms. The spatial distribution of landslides largely depends on the intensity of the triggering event (i.e., the intensity or PGA of the earthquake). However, our method only establishes the relationship between the intensity of a single event and the occurrence of landslides. As such, this method is not yet suitable for landslide prediction related to future triggering events of different sizes. In future work, we intend to use more landslide inventories to thoroughly investigate the relationship between seismic event intensity and landslides in order to improve the generalization ability of our models.

5. Conclusions

This paper proposes a landslide susceptibility model which combines the XGB (or LGB) and the DCE loss function for predicting the susceptibility of a specific part of a landslide (SSPL). We systematically explored how different susceptibility algorithms (XGB, LGB, and RF) and class balancing methods (EQS, ILW, NB, and DCE) affect the prediction of the SSPL. For our study, we prepared six different datasets with 24 landslide influencing factors and 10,422 samples of a specific part of the landslides from the Hokkaido earthquake on 6 September 2018, and established samples of realistic landslide frequencies. Our results demonstrate that XGB and LGB outperform RF and LDA.
The class balancing method heavily influences the predictive performance of different landslide susceptibility models. The low recall of NB methods results in a tendency to underestimate the landslide susceptibility of a given area. NB models are suitable for low-intensity events (i.e., heavy rainfall and small earthquakes) that can trigger landslides. The EQS and ILW methods result in low precision values and significant overestimation of landslide susceptibility. Those two methods should only be applied in large-scale landslide assessments with low resolutions. Predicting a smaller part of a landslide is more difficult than predicting an entire landslide, and a model with a reasonable class balancing method is required. The proposed XGBDCE and LGBDCE models produce more balanced precision and recall values and predict the SSPL more accurately. When combined with the use of an automatic extraction algorithm on different parts of the landslide, these models are applicable for assessing landslide susceptibility, hazard, and risk at various scales and levels of landslide frequency.
In future research, higher quality data should be collected and applied to further improve the predictive performance of the models by focusing on smaller landslide features. A second possible improvement is to include more earthquake-induced landslide inventories and landslide influencing factors in future iterations. Ideally, our model should be able to evaluate landslide susceptibility either prior to a seismic event or in near real-time after the earthquake has begun.

Author Contributions

Conceptualization, S.Z., Y.W. and G.W.; methodology and software, S.Z.; Validation and Visualization, S.Z. and Y.W.; Writing—Original Draft Preparation, S.Z.; Writing—Review and Editing, S.Z., Y.W. and G.W. All of the authors contributed to this article. All authors have read and agreed to the published version of the manuscript.

Funding

This research was supported by the National Natural Science Foundation of China (Grant No. 41602293) and the Science and Technology Development Program of the China Railway Corporation (Grant No. 2010G004-I).

Data Availability Statement

The available data in this study and their access addresses have been summarized in the basic data section.

Acknowledgments

The authors sincerely thank Liming Zhang for his suggestions on revising the paper, as well as the editors and reviewers of the Journal of Remote Sensing.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Guzzetti, F.; Reichenbach, P.; Cardinali, M.; Galli, M.; Ardizzone, F. Probabilistic landslide hazard assessment at the basin scale. Geomorphology 2005, 72, 272–299. [Google Scholar] [CrossRef]
  2. Rossi, M.; Guzzetti, F.; Reichenbach, P.; Mondini, A.C.; Peruccacci, S. Optimal landslide susceptibility zonation based on multiple forecasts. Geomorphology 2010, 114, 129–142. [Google Scholar] [CrossRef]
  3. Pradhan, B. A comparative study on the predictive ability of the decision tree, support vector machine and neuro-fuzzy models in landslide susceptibility mapping using GIS. Comput. Geosci. 2013, 51, 350–365. [Google Scholar] [CrossRef]
  4. Tsangaratos, P.; Ilia, I. Comparison of a logistic regression and Naïve Bayes classifier in landslide susceptibility assessments: The influence of models complexity and training dataset size. CATENA 2016, 145, 164–179. [Google Scholar] [CrossRef]
  5. Liu, Y.; Zhang, W.; Zhang, Z.; Xu, Q.; Li, W. Risk Factor Detection and Landslide Susceptibility Mapping Using Geo-Detector and Random Forest Models: The 2018 Hokkaido Eastern Iburi Earthquake. Remote Sens. 2021, 13, 1157. [Google Scholar] [CrossRef]
  6. Youssef, A.M.; Pourghasemi, H.R. Landslide susceptibility mapping using machine learning algorithms and comparison of their performance at Abha Basin, Asir Region, Saudi Arabia. Geosci. Front. 2021, 12, 639–655. [Google Scholar] [CrossRef]
  7. Sorbino, G.; Sica, C.; Cascini, L. Susceptibility analysis of shallow landslides source areas using physically based models. Nat. Hazards 2010, 53, 313–332. [Google Scholar] [CrossRef]
  8. Conoscenti, C.; Ciaccio, M.; Caraballo-Arias, N.A.; Gómez-Gutiérrez, Á.; Rotigliano, E.; Agnesi, V. Assessment of susceptibility to earth-flow landslide using logistic regression and multivariate adaptive regression splines: A case of the Belice River basin (western Sicily, Italy). Geomorphology 2015, 242, 49–64. [Google Scholar] [CrossRef]
  9. Trigila, A.; Iadanza, C.; Esposito, C.; Scarascia-Mugnozza, G. Comparison of Logistic Regression and Random Forests techniques for shallow landslide susceptibility assessment in Giampilieri (NE Sicily, Italy). Geomorphology 2015, 249, 119–136. [Google Scholar] [CrossRef]
  10. Persichillo, M.G.; Bordoni, M.; Meisina, C.; Bartelletti, C.; Barsanti, M.; Giannecchini, R.; D’Amato Avanzi, G.; Galanti, Y.; Cevasco, A.; Brandolini, P.; et al. Shallow landslides susceptibility assessment in different environments. Geomat. Nat. Hazards Risk 2017, 8, 748–771. [Google Scholar] [CrossRef]
  11. Spinetti, C.; Bisson, M.; Tolomei, C.; Colini, L.; Galvani, A.; Moro, M.; Saroli, M.; Sepe, V. Landslide susceptibility mapping by remote sensing and geomorphological data: Case studies on the Sorrentina Peninsula (Southern Italy). GISci. Remote Sens. 2019, 56, 940–965. [Google Scholar] [CrossRef]
  12. Bordoni, M.; Galanti, Y.; Bartelletti, C.; Persichillo, M.G.; Barsanti, M.; Giannecchini, R.; Avanzi, G.D.A.; Cevasco, A.; Brandolini, P.; Galve, J.P.; et al. The influence of the inventory on the determination of the rainfall-induced shallow landslides susceptibility using generalized additive models. CATENA 2020, 193, 104630. [Google Scholar] [CrossRef]
  13. Dou, J.; Yunus, A.P.; Merghadi, A.; Shirzadi, A.; Nguyen, H.; Hussain, Y.; Avtar, R.; Chen, Y.; Pham, B.T.; Yamagishi, H. Different sampling strategies for predicting landslide susceptibilities are deemed less consequential with deep learning. Sci. Total Environ. 2020, 720, 137320. [Google Scholar] [CrossRef] [PubMed]
  14. Ozturk, U.; Pittore, M.; Behling, R.; Roessner, S.; Andreani, L.; Korup, O. How robust are landslide susceptibility estimates? Landslides 2021, 18, 681–695. [Google Scholar] [CrossRef]
  15. Wang, L.-J.; Sawada, K.; Moriguchi, S. Landslide susceptibility analysis with logistic regression model based on FCM sampling strategy. Comput. Geosci. 2013, 57, 81–92. [Google Scholar] [CrossRef]
  16. Gupta, S.K.; Shukla, D.P. Data Imbalance in Landslide Susceptibility Zonation: A Case Study of Mandakini River Basin, Uttarakhand, India. In Proceedings of the IGARSS 2020 IEEE International Geoscience and Remote Sensing Symposium, Waikoloa, HI, USA, 26 September–2 October 2020; pp. 5234–5237. [Google Scholar] [CrossRef]
  17. Zhang, Y.; Yan, Q. Landslide Susceptibility Prediction Based on High-Trust Non-Landslide Point Selection. ISPRS Int. J. Geo-Inf. 2022, 11, 398. [Google Scholar] [CrossRef]
  18. Huang, F.; Yin, K.; Huang, J.; Gui, L.; Wang, P. Landslide susceptibility mapping based on self-organizing-map network and extreme learning machine. Eng. Geol. 2017, 223, 11–22. [Google Scholar] [CrossRef]
  19. Wang, Y.; Wu, X.; Chen, Z.; Ren, F.; Feng, L.; Du, Q. Optimizing the Predictive Ability of Machine Learning Methods for Landslide Susceptibility Mapping Using SMOTE for Lishui City in Zhejiang Province, China. Int. J. Environ. Res. Public Health 2019, 16, 368. [Google Scholar] [CrossRef] [Green Version]
  20. Zhang, S.; Yu, P. Seismic landslide susceptibility assessment based on ADASYN-LDA model. IOP Conf. Ser. Earth Environ. Sci. 2020, 525, 012087. [Google Scholar] [CrossRef]
  21. Ling, C.X.; Sheng, V.S. Cost-Sensitive Learning. In Encyclopedia of Machine Learning; Sammut, C., Webb, G.I., Eds.; Springer: Boston, MA, USA, 2010; pp. 231–235. [Google Scholar]
  22. Carrara, A.; Crosta, G.; Frattini, P. Comparing models of debris-flow susceptibility in the alpine environment. Geomorphology 2008, 94, 353–378. [Google Scholar] [CrossRef]
  23. Regmi, N.R.; Giardino, J.R.; McDonald, E.V.; Vitek, J.D. A comparison of logistic regression-based models of susceptibility to landslides in western Colorado, USA. Landslides 2014, 11, 247–262. [Google Scholar] [CrossRef]
  24. van Westen, C.J.; Rengers, N.; Soeters, R. Use of Geomorphological Information in Indirect Landslide Susceptibility Assessment. Nat. Hazards 2003, 30, 399–419. [Google Scholar] [CrossRef]
  25. Thiery, Y.; Malet, J.P.; Sterlacchini, S.; Puissant, A.; Maquaire, O. Landslide susceptibility assessment by bivariate methods at large scales: Application to a complex mountainous environment. Geomorphology 2007, 92, 38–59. [Google Scholar] [CrossRef] [Green Version]
  26. Xia, D.; Tang, H.; Sun, S.; Tang, C.; Zhang, B. Landslide Susceptibility Mapping Based on the Germinal Center Optimization Algorithm and Support Vector Classification. Remote Sens. 2022, 14, 2707. [Google Scholar] [CrossRef]
  27. Tanyu, B.F.; Abbaspour, A.; Alimohammadlou, Y.; Tecuci, G. Landslide susceptibility analyses using Random Forest, C4.5, and C5.0 with balanced and unbalanced datasets. CATENA 2021, 203, 105355. [Google Scholar] [CrossRef]
  28. Ageenko, A.; Hansen, L.C.; Lyng, K.L.; Bodum, L.; Arsanjani, J.J. Landslide Susceptibility Mapping Using Machine Learning: A Danish Case Study. ISPRS Int. J. Geo-Inf. 2022, 11, 324. [Google Scholar] [CrossRef]
  29. Chang, Z.; Du, Z.; Zhang, F.; Huang, F.; Chen, J.; Li, W.; Guo, Z. Landslide Susceptibility Prediction Based on Remote Sensing Images and GIS: Comparisons of Supervised and Unsupervised Machine Learning Models. Remote Sens. 2020, 12, 502. [Google Scholar] [CrossRef] [Green Version]
  30. Huang, F.; Cao, Z.; Guo, J.; Jiang, S.-H.; Li, S.; Guo, Z. Comparisons of heuristic, general statistical and machine learning models for landslide susceptibility prediction and mapping. CATENA 2020, 191, 104580. [Google Scholar] [CrossRef]
  31. Yao, K.; Yang, S.; Wu, S.; Tong, B. Landslide Susceptibility Assessment Considering Spatial Agglomeration and Dispersion Characteristics: A Case Study of Bijie City in Guizhou Province, China. ISPRS Int. J. Geo-Inf. 2022, 11, 269. [Google Scholar] [CrossRef]
  32. Jiang, S.-H.; Huang, J.; Huang, F.; Yang, J.; Yao, C.; Zhou, C.-B. Modelling of spatial variability of soil undrained shear strength by conditional random fields for slope reliability analysis. Appl. Math. Model. 2018, 63, 374–389. [Google Scholar] [CrossRef]
  33. Bui, D.T.; Tsangaratos, P.; Nguyen, V.-T.; Liem, N.V.; Trinh, P.T. Comparing the prediction performance of a Deep Learning Neural Network model with conventional machine learning models in landslide susceptibility assessment. CATENA 2020, 188, 104426. [Google Scholar] [CrossRef]
  34. Huang, F.; Cao, Z.; Jiang, S.-H.; Zhou, C.; Huang, J.; Guo, Z. Landslide susceptibility prediction based on a semi-supervised multiple-layer perceptron model. Landslides 2020, 17, 2919–2930. [Google Scholar] [CrossRef]
  35. Bravo-López, E.; Fernández Del Castillo, T.; Sellers, C.; Delgado-García, J. Landslide Susceptibility Mapping of Landslides with Artificial Neural Networks: Multi-Approach Analysis of Backpropagation Algorithm Applying the Neuralnet Package in Cuenca, Ecuador. Remote Sens. 2022, 14, 3495. [Google Scholar] [CrossRef]
  36. Huang, F.; Zhang, J.; Zhou, C.; Wang, Y.; Huang, J.; Zhu, L. A deep learning algorithm using a fully connected sparse autoencoder neural network for landslide susceptibility prediction. Landslides 2020, 17, 217–229. [Google Scholar] [CrossRef]
  37. Chen, T.; Guestrin, C. XGBoost: A Scalable Tree Boosting System. In Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, San Francisco, CA, USA, 13–17 August 2016; pp. 785–794. [Google Scholar] [CrossRef] [Green Version]
  38. Ke, G.; Meng, Q.; Finley, T.; Wang, T.; Chen, W.; Ma, W.; Ye, Q.; Liu, T.-Y. LightGBM: A highly efficient gradient boosting decision tree. Adv. Neural Inf. Process. Syst. 2017, 30, 3149–3157. [Google Scholar] [CrossRef]
  39. Zhang, S.; Li, R.; Wang, F.; Iio, A. Characteristics of landslides triggered by the 2018 Hokkaido Eastern Iburi earthquake, Northern Japan. Landslides 2019, 16, 1691–1708. [Google Scholar] [CrossRef]
  40. Larsen, I.J.; Montgomery, D.R.; Korup, O. Landslide erosion controlled by hillslope material. Nat. Geosci. 2010, 3, 247–251. [Google Scholar] [CrossRef]
  41. Geli, L.; Bard, P.-Y.; Jullien, B. The effect of topography on earthquake ground motion: A review and new results. Bull. Seismol. Soc. Am. 1988, 78, 42–63. [Google Scholar] [CrossRef]
  42. David, P. Earthquake Induced Landslides Lessons from Taiwan and Pakistan; Chengdu University of Technology: Chengdu, China, 2008. [Google Scholar]
  43. Iwahashi, J.; Pike, R.J. Automated classifications of topography from DEMs by an unsupervised nested-means algorithm and a three-part geometric signature. Geomorphology 2007, 86, 409–440. [Google Scholar] [CrossRef]
  44. Zakerinejad, R.; Maerker, M. An integrated assessment of soil erosion dynamics with special emphasis on gully erosion in the Mazayjan basin, southwestern Iran. Nat. Hazards 2015, 79, 25–50. [Google Scholar] [CrossRef]
  45. Mergili, M.; Chu, H.J. Integrated statistical modelling of spatial landslide probability. Nat. Hazards Earth Syst. Sci. Discuss. 2015, 2015, 5677–5715. [Google Scholar] [CrossRef] [Green Version]
  46. Fawcett, T. An introduction to ROC analysis. Pattern Recognit. Lett. 2006, 27, 861–874. [Google Scholar] [CrossRef]
  47. Hastie, T.; Friedman, J.; Tibshirani, R. Linear Methods for Classification. In The Elements of Statistical Learning: Data Mining, Inference, and Prediction; Springer: New York, NY, USA, 2009; pp. 79–113. [Google Scholar]
  48. Breiman, L. Random forests. Mach. Learn. 2001, 45, 5–32. [Google Scholar] [CrossRef] [Green Version]
  49. Milletari, F.; Navab, N.; Ahmadi, S. V-Net: Fully Convolutional Neural Networks for Volumetric Medical Image Segmentation. In Proceedings of the 2016 Fourth International Conference on 3D Vision (3DV), Stanford, CA, USA, 25–28 October 2016; pp. 565–571. [Google Scholar]
  50. Kanda, T.; Takata, Y.; Kohyama, K.; Ohkura, T.; Maejima, Y.; Wakabayashi, S.; Obara, H. New Soil Maps of Japan based on the Comprehensive Soil Classification System of Japan First Approximation and its Application to the World Reference Base for Soil Resources 2006. Jpn. Agric. Res. Q. JARQ 2018, 52, 285–292. [Google Scholar] [CrossRef] [Green Version]
  51. Katsenis, L.C.; Stamatopoulos, C.A.; Panoskaltsis, V.P.; Di, B. Prediction of large seismic sliding movement of slopes using a 2-body non-linear dynamic model with a rotating stick-slip element. Soil Dyn. Earthq. Eng. 2020, 129, 105953. [Google Scholar] [CrossRef]
  52. Ambraseys, N.; Srbulov, M. Earthquake induced displacements of slopes. Soil Dyn. Earthq. Eng. 1995, 14, 59–71. [Google Scholar] [CrossRef]
  53. Rathje, E.M.; Bray, J.D. Nonlinear Coupled Seismic Sliding Analysis of Earth Structures. J. Geotech. Geoenviron. Eng. 2000, 126, 1002–1014. [Google Scholar] [CrossRef]
  54. Ji, J.; Wang, C.-W.; Gao, Y.; Zhang, L. Probabilistic investigation of the seismic displacement of earth slopes under stochastic ground motion: A rotational sliding block analysis. Can. Geotech. J. 2021, 58, 952–968. [Google Scholar] [CrossRef]
  55. Youssef, A.M.; Pourghasemi, H.R.; Pourtaghi, Z.S.; Al-Katheeri, M.M. Landslide susceptibility mapping using random forest, boosted regression tree, classification and regression tree, and general linear models and comparison of their performance at Wadi Tayyah Basin, Asir Region, Saudi Arabia. Landslides 2016, 13, 839–856. [Google Scholar] [CrossRef]
  56. Xu, C.; Xu, X.; Yao, X.; Dai, F. Three (nearly) complete inventories of landslides triggered by the May 12, 2008 Wenchuan Mw 7.9 earthquake of China and their spatial distribution statistical analysis. Landslides 2014, 11, 441–461. [Google Scholar] [CrossRef] [Green Version]
  57. Roback, K.; Clark, M.K.; West, A.J.; Zekkos, D.; Li, G.; Gallen, S.F.; Chamlagain, D.; Godt, J.W. The size, distribution, and mobility of landslides caused by the 2015 Mw7.8 Gorkha earthquake, Nepal. Geomorphology 2018, 301, 121–138. [Google Scholar] [CrossRef]
  58. Regmi, A.D.; Dhital, M.R.; Zhang, J.-Q.; Su, L.-J.; Chen, X.-Q. Landslide susceptibility assessment of the region affected by the 25 April 2015 Gorkha earthquake of Nepal. J. Mt. Sci. 2016, 13, 1941–1957. [Google Scholar] [CrossRef]
  59. Mergili, M.; Krenn, J.; Chu, H.J. r. randomwalk v1, a multi-functional conceptual tool for mass movement routing. Geosci. Model Dev. 2015, 8, 4027–4043. [Google Scholar] [CrossRef]
Figure 1. The scheme of the present study.
Figure 1. The scheme of the present study.
Remotesensing 14 05945 g001
Figure 2. Overview map of the study area. (a) Geologic map at a scale of 1:200,000 overlain by the distribution of co-seismic landslides (source: geologic map in Table 1. The classification of rock groups has been modified.). (b) The locations and satellite imagery (source: http://goto.arcgisonline.com/maps/World_Imagery accessed on 1 September 2022) of the study area. (c) Distribution of peak ground acceleration in the Hokkaido region during the 6 September 2018, earthquake (source: QuiQuake in Table 1). (d) Satellite imagery (source: Maxar imagery acquired on 11 September 2018) of the co-seismic landslides in the local zone in (a).
Figure 2. Overview map of the study area. (a) Geologic map at a scale of 1:200,000 overlain by the distribution of co-seismic landslides (source: geologic map in Table 1. The classification of rock groups has been modified.). (b) The locations and satellite imagery (source: http://goto.arcgisonline.com/maps/World_Imagery accessed on 1 September 2022) of the study area. (c) Distribution of peak ground acceleration in the Hokkaido region during the 6 September 2018, earthquake (source: QuiQuake in Table 1). (d) Satellite imagery (source: Maxar imagery acquired on 11 September 2018) of the co-seismic landslides in the local zone in (a).
Remotesensing 14 05945 g002
Figure 3. The landslide influencing factors used in this study, obtained and processed by remote sensing and GIS technology. The yellow pentagram represents the epicenter. The names of the factors in (aw) are below the color bar.
Figure 3. The landslide influencing factors used in this study, obtained and processed by remote sensing and GIS technology. The yellow pentagram represents the epicenter. The names of the factors in (aw) are below the color bar.
Remotesensing 14 05945 g003
Figure 4. Automatic extraction method for different parts of a landslide. (a) The division and extraction of the landslide parts. (b) A flowchart of our automatic extraction methodology. (cf) The results of the landslide parts corresponding to Rd values of (c) 0.1, (d) 0.3, (e) 0.5, and (f) 0.9. The solid black lines represent the boundaries of the landslides. The colored cells represent the extracted parts of the landslide.
Figure 4. Automatic extraction method for different parts of a landslide. (a) The division and extraction of the landslide parts. (b) A flowchart of our automatic extraction methodology. (cf) The results of the landslide parts corresponding to Rd values of (c) 0.1, (d) 0.3, (e) 0.5, and (f) 0.9. The solid black lines represent the boundaries of the landslides. The colored cells represent the extracted parts of the landslide.
Remotesensing 14 05945 g004
Figure 5. The influence that parameters a and b in the DCE loss function exert on landslide susceptibility prediction using the (a,b,e,f) XGBDCE and (c,d,g,h) LGBDCE models. As Rd increases, the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp.
Figure 5. The influence that parameters a and b in the DCE loss function exert on landslide susceptibility prediction using the (a,b,e,f) XGBDCE and (c,d,g,h) LGBDCE models. As Rd increases, the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp.
Remotesensing 14 05945 g005
Figure 6. Prediction of the susceptibilities of various landslide parts with combinations of the XGB algorithm with different class balancing methods: (ac) XGBNB, (df) XGBEQS, (gi) XGBDCE. As Rd increases, the part of the landslide that the model seeks to predict (the black polygon) is reduced from the entire landslide to the landslide scarp.
Figure 6. Prediction of the susceptibilities of various landslide parts with combinations of the XGB algorithm with different class balancing methods: (ac) XGBNB, (df) XGBEQS, (gi) XGBDCE. As Rd increases, the part of the landslide that the model seeks to predict (the black polygon) is reduced from the entire landslide to the landslide scarp.
Remotesensing 14 05945 g006
Figure 7. Prediction of the susceptibilities of various landslide parts with combinations of the RF algorithm with different class balancing methods: (ac) RFNB, (df) RFEQS, (gi) RFILW. As Rd increases, the part of the landslide that the model seeks to predict (the black polygon) is reduced from the entire landslide to the landslide scarp.
Figure 7. Prediction of the susceptibilities of various landslide parts with combinations of the RF algorithm with different class balancing methods: (ac) RFNB, (df) RFEQS, (gi) RFILW. As Rd increases, the part of the landslide that the model seeks to predict (the black polygon) is reduced from the entire landslide to the landslide scarp.
Remotesensing 14 05945 g007
Figure 8. The relationship between the proportion of landslides and the landslide susceptibility. From subfigures (af), the Rd increases, and the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp. Because some of the plots almost overlap (XGB-based model versus LGB-based model and EQS-based model versus ILW-based model), for greater clarity we only show representative model plots.
Figure 8. The relationship between the proportion of landslides and the landslide susceptibility. From subfigures (af), the Rd increases, and the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp. Because some of the plots almost overlap (XGB-based model versus LGB-based model and EQS-based model versus ILW-based model), for greater clarity we only show representative model plots.
Remotesensing 14 05945 g008
Figure 9. The relationship between the landslide frequency and the landslide susceptibility. The models with the plots closest to the reference lines have predicted susceptibility values that are close to the actual landslide frequency. (af), the Rd increases, and the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp. Because some of the plots almost overlap (XGB-based model versus LGB-based model and EQS-based model versus ILW-based model), for clarity we show only representative model plots.
Figure 9. The relationship between the landslide frequency and the landslide susceptibility. The models with the plots closest to the reference lines have predicted susceptibility values that are close to the actual landslide frequency. (af), the Rd increases, and the part of the landslide that the model seeks to predict is reduced from the entire landslide to the landslide scarp. Because some of the plots almost overlap (XGB-based model versus LGB-based model and EQS-based model versus ILW-based model), for clarity we show only representative model plots.
Remotesensing 14 05945 g009
Figure 10. The ranking of the importance of different landslide influencing factors using the (a) XGBNB, (b) XGBEQS, (c) XGBDCE, and (d) XGBILW models. Different Rd (0.0–0.9) values correspond to predictions of different landslide parts in the models when landslide influencing factors are ranked by importance. Mean and STD refer to the mean and standard deviation, respectively, of the ranking of each landslide factor’s importance for all Rd values. Lower rankings represent more significant factor contributions. EL = elevation, diff = difference, and DIR = direction.
Figure 10. The ranking of the importance of different landslide influencing factors using the (a) XGBNB, (b) XGBEQS, (c) XGBDCE, and (d) XGBILW models. Different Rd (0.0–0.9) values correspond to predictions of different landslide parts in the models when landslide influencing factors are ranked by importance. Mean and STD refer to the mean and standard deviation, respectively, of the ranking of each landslide factor’s importance for all Rd values. Lower rankings represent more significant factor contributions. EL = elevation, diff = difference, and DIR = direction.
Remotesensing 14 05945 g010
Figure 11. The frequencies of the various parts of the landslide with respect to the landslide influencing factors (subfigures (al)). In Subfigure (h), Un (unknown lithology), Mm (marine mudstone), Sm (siliceous mudstone), Mc (marine conglomerate), Ss (Sandstone), and Qd (Quaternary deposits) denote the lithology.
Figure 11. The frequencies of the various parts of the landslide with respect to the landslide influencing factors (subfigures (al)). In Subfigure (h), Un (unknown lithology), Mm (marine mudstone), Sm (siliceous mudstone), Mc (marine conglomerate), Ss (Sandstone), and Qd (Quaternary deposits) denote the lithology.
Remotesensing 14 05945 g011
Figure 12. A combination of landslide source area prediction and landslide run-out modelling portraying the respective distribution of (a) the landslide impact frequency (expressed as the cumulative distribution function) and (b) the landslide impact probability. The landslide run-out model uses r.randomwalk [59], in which the motion is controlled by the landslide angle of reach. The average angle of reach (25°) and the probability distribution function used for modeling were calculated during the statistical analysis of the landslides.
Figure 12. A combination of landslide source area prediction and landslide run-out modelling portraying the respective distribution of (a) the landslide impact frequency (expressed as the cumulative distribution function) and (b) the landslide impact probability. The landslide run-out model uses r.randomwalk [59], in which the motion is controlled by the landslide angle of reach. The average angle of reach (25°) and the probability distribution function used for modeling were calculated during the statistical analysis of the landslides.
Remotesensing 14 05945 g012
Table 1. Data set sources and meta-data.
Table 1. Data set sources and meta-data.
Basic DataSourceDateResolution
Digital elevation model (DEM)Geospatial Authority Institute of Japan
(https://fgd.gsi.go.jp/download/menu.php accessed on 1 September 2022)
1 October 201610 m
QuiQuakeGeological Survey of Japan, AIST
(https://gbank.gsj.jp/QuiQuake/QuakeMap/index.en.html accessed on 1 September 2022)
6 September 2018250 m
ShakemapUnited States Geological Survey (USGS)
(https://earthquake.usgs.gov/data/shakemap/ accessed on 1 September 2022)
6 September 2018
CHIRPSClimate Hazards Center
(https://www.chc.ucsb.edu/data/chirps accessed on 1 September 2022)
3 September 2018–6 September 20185000 m
MODIS Vegetation IndexNASA (https://earthdata.nasa.gov/ accessed on 1 September 2022)13 August 2018–28 August 2018250 m
Geological mapGeological Survey of Japan, AIST
(https://gbank.gsj.jp/seamless/v2.html accessed on 1 September 2022)
22 January 20211:200,000
Google earth imageryGoogle earth pro30 September 2016–10 July 20200.2 m
Aerial photosGeospatial Authority Institute of Japan
(https://kmlnetworklink.gsi.go.jp/kmlnetworklink/index.html accessed on 1 September 2022)
6 September 2018–13 September 20180.2 m
QuiQuake = Quick estimation system for earthquake maps triggered by observation records; CHIRPS = Climate Hazards Group InfraRed Precipitation with Station data.
Table 2. Description of factors influencing landslides.
Table 2. Description of factors influencing landslides.
TypeFactorBasic DataRangeUnitResolution
SeismicPeak ground acceleration (PGA)QuiQuake11.0–131.6g%250 m
PGVA, the product of PGV (peak ground velocity) and PGA4.4–9.4250 m
The Euclidean distance to the focus (distancefocus)Shakemap35.0–52.2km30 m
Epicentral direction0.0–180.0°30 m
The Euclidean distance to the nearest ridge (distanceridge)DEM0.0–488.4m30 m
The angle between epicentral direction and the slope aspect (angleES)DEM and
QuiQuake
0.0–360.0°30 m
The angle between the horizontal and the line from calculated cell to focus (angleFH). AngleFH represents the direction of seismic wave propagation at a location, which would influence the occurrence of landslides.42.9–90.0°30 m
The sum of angleFH and the slope degree (angleFS)43.7–126.8°30 m
GeomorphologicThe maximum slope in the neighbourhood of the calculated cell (slopeMAX)DEM0.0–57.7°30 m
The variation of the slope aspects in the neighbourhood of the calculated cell (aspectVAR)0.0–1.030 m
The ratio of the elevation to the maximum elevation in the neighbourhood of the calculated cell (elevation ratio)0.0–1.030 m
Elevation difference0.0–233.0m30 m
The percentage of convex cells in the neighbourhood of the calculated cell (surface convexity) [43]. Surface convexity describes the shape (convex, concave, flat) of the slope, which affects the stability of the slope under earthquake shaking.1.2–74.430 m
The standard deviation of the curvature in the neighbourhood of the calculated cell (curvatureSTD)0.0–4.130 m
HydrologicalThe shortest Euclidean distance to minor rivers (distanceSR)DEM0.0–1855.6m30 m
The shortest Euclidean distance to major rivers (distanceMR)0.0–10.2km30 m
Stream power index (SPI)−13.8–15.630 m
ClimaticCumulative precipitation in the 4 days before an earthquake (precipitation)CHIRPS0.0–39.1mm5000 m
Vegetation coverEnhanced vegetation index (EVI)MODIS Vegetation Index Products−3879.0–9748.0250 m
GeologicalLithologyGeological map
The Euclidean distance to the nearest fault (distancefault)0.0–13.2km30 m
The Euclidean distance to the nearest fold (distancefold)0.0–12.3km30 m
Fault density (LF × WF/AF). LF, WF, and AF are the total fault length, fault width, and area of the statistical zone, respectively.0.0–2.230 m
Fold density (LO × WO/AO). LO, WO, and AO are the total fold length, fold width, and area of the statistical zone, respectively.0.0–1.630 m
Table 3. Inventories of different landslide parts.
Table 3. Inventories of different landslide parts.
RdNumber of Landslide PartsNumber of Landslide CellsNumber of Non-landslide CellsRatio of Landslide to Non-LandslideArea of Landslide CellsTotal Study Area (km2)Landslide Frequency
Sum (km2)Mean (m2)
0.010,42210,48262,101,4791:2094.349052.331985.670.048
0.110,42283,6632,122,6421:2575.307224.781985.670.038
0.310,42266,1202,140,1851:3259.515709.841985.670.030
0.510,42248,8642,157,4411:4443.984219.691985.670.022
0.710,42231,1522,175,1531:7028.042690.161985.670.014
0.910,42216,0612,190,2441:13614.451386.961985.670.007
Landslide frequency: the ratio of the sum of the areas of the landslide parts to the total study area.
Table 4. Model prediction performance of different landslide features for every model.
Table 4. Model prediction performance of different landslide features for every model.
ModelRd = 0.0Rd = 0.1
PRERECF1JACPJACNJACAUCFreqPRERECF1JACPJACNJACAUCFreq
XGBNB0.6990.4680.5610.3900.9640.6130.9740.0480.6460.3540.4570.2970.9680.5360.9710.038
XGBEQS0.3500.9300.5090.3410.9110.5570.9720.0480.2900.9320.4420.2840.9080.5080.9700.038
XGBILW0.3660.9160.5230.3540.9170.5700.9720.0480.3070.9160.4600.2990.9160.5230.9700.038
XGBDCE0.5200.7710.6210.4510.9540.6550.9730.0480.4640.7270.5660.3950.9570.6150.9710.038
LGBNB0.6390.4290.5120.3440.9600.5750.9680.0480.5940.2900.3890.2420.9650.4830.9660.038
LGBEQS0.3200.9420.4770.3130.8970.5300.9710.0480.2660.9450.4150.2620.8950.4840.9690.038
LGBILW0.3220.9410.4800.3160.8980.5320.9710.0480.2680.9440.4180.2640.8970.4860.9690.038
LGBDCE0.4710.7970.5920.4210.9460.6310.9700.0480.4190.7600.5400.3700.9490.5920.9680.038
RFNB0.7510.4140.5340.3640.9650.5930.9690.0480.6560.2740.3870.2400.9670.4810.9630.038
RFEQS0.3870.9120.5430.3730.9240.5870.9720.0480.3200.9150.4750.3110.9200.5350.9690.038
RFILW0.7510.3790.5040.3370.9640.5700.9630.0480.6520.2470.3580.2180.9660.4590.9570.038
LDANB0.4260.1780.2510.1440.9490.3690.9250.0480.3670.1290.1910.1060.9580.3180.9220.038
LDAEQS0.1940.9490.3220.1920.8020.3930.9360.0480.1570.9500.2700.1560.7980.3530.9320.038
ModelRd= 0.3Rd= 0.5
PRERECF1JACPJACNJACAUCFreqPRERECF1JACPJACNJACAUCFreq
XGBNB0.6310.2740.3820.2360.9730.4790.9710.0300.6230.1980.3010.1770.9800.4160.9710.022
XGBEQS0.2400.9310.3820.2360.9070.4630.9690.0300.1830.9320.3050.1800.9040.4040.9680.022
XGBILW0.2600.9130.4040.2530.9170.4820.9700.0300.2060.9060.3350.2010.9190.4300.9700.022
XGBDCE0.4310.6770.5270.3570.9630.5870.9710.0300.3980.5980.4780.3140.9710.5520.9710.022
LGBNB0.4800.2120.2800.1630.9650.3950.9600.0300.4800.1400.2140.1200.9770.3420.9640.022
LGBEQS0.2220.9450.3590.2190.8960.4430.9690.0300.1710.9420.2890.1690.8950.3890.9690.022
LGBILW0.2250.9430.3640.2220.8980.4470.9690.0300.1760.9390.2970.1740.8990.3960.9690.022
LGBDCE0.3910.7090.5040.3370.9570.5680.9680.0300.3510.6320.4520.2920.9660.5310.9670.022
RFNB0.6330.2180.3240.1930.9730.4340.9610.0300.6100.1680.2630.1520.9790.3850.9580.022
RFEQS0.2630.9190.4090.2570.9180.4860.9690.0300.1970.9250.3250.1940.9130.4210.9690.022
RFILW0.6240.1950.2980.1750.9720.4120.9530.0300.6040.1530.2440.1390.9790.3690.9500.022
LDANB0.3540.0990.1540.0840.9670.2840.9210.0300.3440.0640.1080.0570.9770.2360.9220.022
LDAEQS0.1260.9490.2230.1250.7960.3160.9320.0300.0950.9470.1730.0940.7950.2740.9330.022
ModelRd= 0.7Rd= 0.9
PRERECF1JACPJACNJACAUCFreqPRERECF1JACPJACNJACAUCFreq
XGBNB0.6020.1170.1960.1090.9860.3280.9700.0140.5550.0340.0640.0330.9930.1810.9650.007
XGBEQS0.1170.9320.2080.1160.8990.3230.9670.0140.0550.9320.1050.0550.8830.2210.9600.007
XGBILW0.1390.8980.2410.1370.9190.3540.9680.0140.0720.8740.1330.0710.9160.2550.9630.007
XGBDCE0.3530.4810.4070.2550.9800.5000.9700.0140.2660.2770.2710.1570.9890.3940.9650.007
LGBNB0.3960.1200.1820.1000.9850.3140.9630.0140.2260.0600.0930.0490.9910.2180.9530.007
LGBEQS0.1110.9410.1990.1100.8920.3140.9680.0140.0540.9370.1020.0540.8790.2170.9620.007
LGBILW0.1170.9350.2070.1160.8980.3220.9690.0140.0580.9290.1090.0580.8890.2270.9630.007
LGBDCE0.2880.5110.3660.2250.9740.4660.9630.0140.2300.2960.2590.1490.9880.3830.9590.007
RFNB0.5980.1110.1870.1030.9860.3190.9520.0140.5000.0350.0660.0340.9930.1840.9320.007
RFEQS0.1230.9310.2170.1220.9040.3310.9680.0140.0560.9360.1060.0560.8840.2220.9630.007
RFILW0.5800.1000.1700.0930.9860.3030.9410.0140.4950.0310.0590.0300.9930.1740.9160.007
LDANB0.3480.0350.0640.0330.9850.1810.9230.0140.3400.0110.0210.0110.9930.1030.9160.007
LDAEQS0.0620.9440.1160.0610.7940.2210.9340.0140.0310.9360.0600.0310.7870.1570.9300.007
Freq: landslide frequency. PRE: precision. F1:F1 score. JACP: JaccardP. JACN: JaccardN. JAC: Jaccard.
Table 5. Statistical analysis of the number of cells (unit: thousand) in the generated landslide susceptibility maps for every model.
Table 5. Statistical analysis of the number of cells (unit: thousand) in the generated landslide susceptibility maps for every model.
Rd = 0.0Rd = 0.3
ModelLandslide SusceptibilityFreqLandslide SusceptibilityFreq
0.0–0.20.2–0.40.4–0.60.6–0.80.8–1.00.0–0.20.2–0.40.4–0.60.6–0.80.8–1.0
XGBNB2037.572.447.933.914.50.0482094.264.332.213.32.20.030
XGBEQS1821.076.957.965.8184.80.0481840.178.557.564.8165.40.030
XGBILW1842.873.554.060.2175.90.0481874.672.052.457.2150.10.030
XGBDCE1969.058.846.757.074.80.0482027.653.143.346.336.00.030
LGBNB2023.487.053.536.36.00.0482087.674.331.111.12.10.030
LGBEQS1752.6105.772.878.7196.50.0481778.1104.473.475.8174.50.030
LGBILW1752.7104.973.478.6196.60.0481782.2105.472.674.4171.70.030
LGBDCE1929.970.558.274.773.00.0481997.862.453.762.330.10.030
RFNB2069.358.837.426.314.40.0482119.250.323.610.32.90.030
RFEQS1834.990.764.671.3144.80.0481829.5107.271.873.4124.40.030
RFILW2078.955.035.524.912.00.0482128.944.021.59.62.30.030
LDANB2015.5121.141.617.610.50.0482099.275.420.98.62.30.030
LDAEQS1448.9181.9120.7122.6332.20.0481466.4178.0125.7127.2309.00.030
Rd= 0.5Rd= 0.9
Model0.0–0.20.2–0.40.4–0.60.6–0.80.8–1.0Freq0.0–0.20.2–0.40.4–0.60.6–0.80.8–1.0Freq
XGBNB2125.952.519.57.11.20.0222193.410.71.80.40.00.007
XGBEQS1845.881.959.364.9154.40.0221800.596.071.178.6160.10.007
XGBILW1890.573.152.356.6133.70.0221902.778.457.862.9104.50.007
XGBDCE2064.449.536.534.021.90.0222153.828.313.17.53.50.007
LGBNB2125.756.614.75.83.60.0222193.48.31.40.82.30.007
LGBEQS1784.3112.675.474.9159.10.0221743.0131.089.495.2147.60.007
LGBILW1794.6109.474.373.5154.50.0221779.5122.082.391.5131.00.007
LGBDCE2036.457.647.145.619.70.0222141.534.916.39.73.90.007
RFNB2142.940.415.76.01.40.0222192.411.32.10.40.00.007
RFEQS1814.0120.779.376.8115.50.0221696.1180.2115.5103.6110.90.007
RFILW2150.434.814.45.61.10.0222194.19.91.80.40.00.007
LDANB2142.747.811.04.40.50.0222197.67.41.30.10.00.007
LDAEQS1466.6184.8129.7132.6292.50.0221436.9210.0148.3158.5252.60.007
Freq: landslide frequency.
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Zhang, S.; Wang, Y.; Wu, G. Earthquake-Induced Landslide Susceptibility Assessment Using a Novel Model Based on Gradient Boosting Machine Learning and Class Balancing Methods. Remote Sens. 2022, 14, 5945. https://doi.org/10.3390/rs14235945

AMA Style

Zhang S, Wang Y, Wu G. Earthquake-Induced Landslide Susceptibility Assessment Using a Novel Model Based on Gradient Boosting Machine Learning and Class Balancing Methods. Remote Sensing. 2022; 14(23):5945. https://doi.org/10.3390/rs14235945

Chicago/Turabian Style

Zhang, Shuhao, Yawei Wang, and Guang Wu. 2022. "Earthquake-Induced Landslide Susceptibility Assessment Using a Novel Model Based on Gradient Boosting Machine Learning and Class Balancing Methods" Remote Sensing 14, no. 23: 5945. https://doi.org/10.3390/rs14235945

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop