Rainfall-Induced Shallow Landslide Susceptibility Mapping at Two Adjacent Catchments Using Advanced Machine Learning Algorithms

Abstract

Landslides impact on human activities and socio-economic development, especially in mountainous areas. This study focuses on the comparison of the prediction capability of advanced machine learning techniques for the rainfall-induced shallow landslide susceptibility of Deokjeokri catchment and Karisanri catchment in South Korea. The influencing factors for landslides, i.e., topographic, hydrologic, soil, forest, and geologic factors, are prepared from various sources based on availability, and a multicollinearity test is also performed to select relevant causative factors. The landslide inventory maps of both catchments are obtained from historical information, aerial photographs and performed field surveys. In this study, Deokjeokri catchment is considered as a training area and Karisanri catchment as a testing area. The landslide inventories contain 748 landslide points in training and 219 points in testing areas. Three landslide susceptibility maps using machine learning models, i.e., Random Forest (RF), Extreme Gradient Boosting (XGBoost) and Deep Neural Network (DNN), are prepared and compared. The outcomes of the analyses are validated using the landslide inventory data. A receiver operating characteristic curve (ROC) method is used to verify the results of the models. The results of this study show that the training accuracy of RF is 0.756 and the testing accuracy is 0.703. Similarly, the training accuracy of XGBoost is 0.757 and testing accuracy is 0.74. The prediction of DNN revealed acceptable agreement between the susceptibility map and the existing landslides, with a training accuracy of 0.855 and testing accuracy of 0.802. The results showed that the DNN model achieved lower prediction error and higher accuracy results than other models for shallow landslide modeling in the study area.

1. Introduction

Landslides are one of the major natural hazards worldwide due to unfavorable geological causes, weathering patterns, shallow soil deposits and heavy rainfall. Landslides create dangerous menaces to the environment. The Korean peninsula is currently experiencing climate change effects, including annual temperature, precipitation and the rate of typhoon occurrence [1]. In the Korean Peninsula, one of the main causes of heavy downpour is typhoons. Every year, large-scale typhoons pass over South Korea, leaving a trail of devastation in their path. Shallow-seated landslides often result from heavy downpours accompanying typhoons and torrential rains during the summer season [2].

Landslide prediction is a crucial task because the inherent causes of landslides are complex due to mechanical and chemical weathering, nonhomogeneous soil deposits, the orientation of geological lineament, and the erosion of a slope. The external causes, including rainfall, earthquakes, the excavation of a slope, and the loading of a slope or its crest, are the most powerful agents that trigger the landslide [3,4,5,6,7].

Landslide susceptibility is the occurrence probability of a landslide under certain geo-environmental conditions, and it estimates “where” landslides are most likely to occur in the future [8,9]. In the last three decades, various researchers have proposed and applied GIS-based qualitative and quantitative methods [10]. The qualitative methods are subjective and depend on expert decisions [11,12,13], and quantitative methods are objective and founded on the geographical distribution of landslides and influencing factors (IFs). The weight of evidence, frequency ratio, relative effect models and logistic regression models have been frequently used in quantitative methods [14,15,16,17,18,19]. In recent years, machine learning techniques have been increasingly employed for resolving the non-linear spatial relationship of landslide occurrences, using either regression or classification methods [20,21,22,23,24,25].

Advances in computational technologies and the rapid development of advanced machine learning algorithms have amended the understanding and management of many geo-environmental problems. The main advantage of machine learning over traditional techniques is its capability to handle high-dimensional and complex non-linear datasets [26]. Various previous investigators have utilized machine learning in different fields of geosciences and hydrogeology, including spring potential mapping [27], gully erosion mapping [28], flood hazard mapping [29] and landslide studies [30,31]. The application of advanced machine learning algorithms to landslide susceptibility modeling has been very limited, despite the variety thereof and the potential advantages.

The main aim of this research is to evaluate landslide susceptibility based on machine learning techniques, i.e., Random Forest (RF), Extreme Gradient Boosting (XGBoost) and Deep Neural Network (DNN), in Deokjeokri and Karisanri catchments. The designed focus of this study is landslide susceptibility mapping, with training of the model done in Deokjeokri catchment and testing of the model done in Karisanri catchment, close to Deokjeokri catchment. Two adjacent catchments (Deokjeokri and Karisanri) were selected as study sites because the two sites have similar geological and geomorphological settings, and they had severe rainfall-induced landslide problems.

2. Description of Study Area

Deokjeokri catchment and Karisanri catchment are located within the Inje region, Gangwon Province in the northeastern part of Korea (inset of Figure 1). Deokjeokri catchment occupies about 33.4 km², and Karisanri catchment about 22.2 km². These catchments are confined by moderate to steep mountains. The predominant lithological types of this area include banded gneisses and granites, as shown in Figure 1. The study area consists of residual soil that includes silt and clay. Shallow landslides are the major land degradation process in both catchments. The channels are narrowed and meandered. The major extrinsic cause of landslide in the study area was a heavy downpour. The landslide occurred under the indirect gist of two typhoons, Ewiniar and Bilis, in July 2006. The total rainfall from 14–16 July 2006 was 402 mm, and the maximum hourly rainfall intensity, which occurred on 15 July 2006, was 62 mm/hr. The extreme rainfall that occurred on 14–16 July resulted in 17 deaths and 12 people unaccounted for in Inje County.

3. The Collected Dataset and Methods

3.1. Landslide Inventory

A landslide inventory map records the geographic location, size, date and type of landslide [32,33,34,35,36]. The precise detection of landslide emplacements is crucial for the prediction and assessment of the landslide susceptibility models. A new generation of satellite images, like Word view and Geo-eye, aided a quick reproduction of landslide inventory maps [37,38]. Recently, the application of interferometry techniques to radar images has been extensively used in inventory mapping [39]. However, due to the limitation of resources, in this study, we (1) examined aerial photographs (provided in the web portal site; www.map.kakao.com, www.map.naver.com) to prepare a landslide inventory map. Both web portals give an excellent bird’s eye view with a 25 cm resolution. We also (2) performed fieldwork. The field observations showed that the predominant mass movements are shallow slides that sometimes ensued as debris flows during periods of heavy rainfalls, and occurred on natural and man-made slopes with highly weathered bedrock.

As shown in Figure 2, a total of 748 landslides were identified, and afterward digitized in ArcGIS 10.3 environment, for further analysis in Deokjeokri catchment (training area), and 219 landslides were identified at Karisanri catchment (testing area). Some glimpses of landslide observed during field visits are presented in Figure 3. In this study, a single point feature from the center of each landslide’s initiation area was extracted and rasterized into 10 × 10 m pixel size. The aerial extent of the initiation of each landslide was measured to be less than 100 m². This means it can represent the landslide initiation area based on the least grid size of the digital elevation model (DEM). A Gaussian Kernel density function of the sampling strategy (an unbiased sample technique) was used [40] for non-landslide pixels. A double number of non-landslide pixels (i.e., 1496 from Deokjeokri and 438 from Karisanri) was sampled in each catchment.

3.2. Landslide Influencing Factor (IF)

The main process for susceptibility modeling is the gathering and designing of a spatial database in a GIS environment [41]. The landslide influencing factors (IFs) considered in this study are the intrinsic factors such as topographic, hydrologic, forest, soil and geological factors, which are collected from available resources and fieldworks.

Table 1. Data type, scale and producer.

Spatial Data		Scale	Resolution	Producer
Landslide inventory		Point data		Satellite image, aerial photographs, field survey, ArcGIS 10.2
Topographic factor	Aspect	1:5000	10 × 10 m	ArcGIS 10.2
	Elevation
	Slope
	Internal relief
	Curvature
Hydrologic factor	Drainage proximity (h)	1:5000		ArcGIS 10.2
	Drainage proximity (v)
	Drainage density
	Stream Power Index		10 × 10 m
	Sediment Transport Index
	Topographic Wetness Index
Soil factor	Soil type	1:25,000	10 × 10 m	Korea Forest Service (KFS)
Forest factor	Forest type	1:25,000	10 × 10 m	Korea Forest Service (KFS)
Geological factor	Geology	1:50,000	10 × 10 m	Korean Institute of Geoscience and Mineral Resources (KIGAM)

represents the spatial database obtained from different sources. In this study, based on reviewing previous research [22,23,42,43] and the availability of data, 14 IFs were utilized to identify the hillslope features of landslide occurrence. Among the 14 IFs, topographic (except aspect map) and hydrologic factors are continuous, but forest, soil and geology are categorical factors. The selection of an appropriate pixel is necessary to achieve high precision in landslide susceptibility mapping [44,45]. Tarolli and Tarboton (2006) found that a 10 × 10 m resolution is adequate to achieve better landslide susceptibility prediction performance, and thus a 10 × 10 m resolution DEM was created from 5 m interval contours. Similarly, all collected IFs (thematic layers) were rasterized into 10 × 10 m pixel size.

3.2.1. Topographic Factors

The aspect does not show a direct relation with landslide occurrence, however, the aspect of the terrain is related to parameters such as sunlight exposure and precipitation [46,47]. The aspect was divided into north, northeast, east, southeast, south, southwest, west, northwest and flat, as shown in Figure 4a. The elevation is another important element in landslide studies. One study shows that a higher altitude tends to cause landslides [46]. The elevation of the study area ranges from 195.1 to 1513.67 m asl, as presented in Figure 4b. The slope is a very powerful driver of the landslide. In fact, the slope affects the surface and subsurface flow, and thus the soil moisture content, the soil formation and the likelihood of soil erosion [48]. In the study area, the slope ranges for 0° to 62.3° as shown in Figure 4c. Internal relief indicates the topographic breakage—the potential energy available for mass movements [49]. Internal relief was prepared by calculating the local height difference within nine pixels (unit area), and in the study area it ranges from 0.15 to 72.1 m, as shown in Figure 4d. Terrain curvature represents the morphology of the hillslope (Figure 4e), and it controls the slope erosion process by the convergence or divergence of water flow [50].

3.2.2. Hydrologic Factors

Hydrological factors are vital factors in rainfall-induced slope instability. Two types of drainage proximity, i.e., horizontal and vertical Ifs, were used in this study. The drainage proximity (h) (Figure 5a) was obtained from the Euclidean function in the GIS environment. A new topo-hydrological factor, i.e., drainage proximity (v), was used. It gives the height above the nearest drainage network (Figure 5b). This index considers the height differences along flow paths. Drainage density is the length of the drainage within a pixel (Figure 5c). The presence of higher drainage density means less percolation but faster surface flow [51].

The SPI indicates the erosive processes, as caused by runoff. As the specific catchment area and slope gradient increase, the quantity of water contributed by the upslope and the surface runoff also increases [52]. Equation (1) defines the SPI:

$$SPI=A_s\times tan\beta$$

(1)

where A_s is the specific catchment area, and β is the slope gradient. The spatial distribution of SPI is presented in Figure 5d.

The STI reflects the erosive power of the overflow drain and the process of erosion and deposition. To calculate STI, two components of the slope which are responsible for soil loss, i.e., length (L) and steepness (S), were considered, as suggested by Moore and Burch (1986) [53], given in Equation (2),

$$STI={\left(\frac{A_s}{22.13}\right)}^{0.6}{\left(\frac{sin\beta }{0.0896}\right)}^{1.3}$$

(2)

where A_s is the specific catchment area, and β is the slope gradient. The spatial distribution of STI is shown in Figure 5e.

TWI is related to soil conditions and surface runoff [54]. TWI is the tendency of water to aggregate in the catchment and the propensity of gravitational powers to move downslope [55]. TWI is defined as

$$TWI=log(\frac{\alpha }{tan\beta })$$

(3)

where α is the cumulative upslope area, and β is the slope angle. The distribution of TWI in the study area is shown in Figure 5f.

3.2.3. Forest and Soil Factors

Rickli et al. (2002) [56] addressed the effect of forests in hillslope reinforcement. The crucial role of trees and forests is to forestall mass movements by fortifying and drying soils. The catchments mainly consist of oak, Japanese larch, Japanese red pine, Japanese pine and Korean pine forest, as shown in Figure 6a. The non-forest area consists of agricultural land.

Soil properties affect degrees of drainage and erosion, which influence the landslide occurrence. The particle size and pore distribution affect water movement and the holding of infiltrated water [57]. Finer soils can hold higher volumes of water than coarse-textured soils [58]. In the study area, the predominant soil types are sandy loam soils (Figure 6b). However, downstream areas were covered by silty clay loam and clay loam. The northern and southern elevated areas of the catchments are rocky and covered by a very thin layer of soil.

3.2.4. Geologic Factor

Landslides are significantly determined by the lithological characteristics of the hillslope because the individual lithological unit relates to different weathering degrees [59,60]. As described earlier, the lithology of this area includes banded gneisses and granites (Figure 1). Prominent lineaments and dense joint spacing with weathering and erosions were observed in both catchments, which create weak zones. The drainage system follows the lineaments. The bulky granite peaks are one of the typical characteristics of the Inje area.

3.3. Modeling Approaches

As illustrated in Figure 7, the fundamental working phases of the adopted methodology in this study comprise four main steps: (1) the preparation of landslide inventory and related database, (2) the preparation of training and testing dataset (in this research Deokjeokri catchment was selected as a training area and Karisanri catchment as a testing area), (3) the landslide IFs were tested and selected using multicollinearity test, and (4) the preparation and comparison of the performance of landslide susceptibility maps using RF, XGBoost and DNN models.

The number of IFs can be minimized via multicollinearity testing, reducing high-dimensional data. For the susceptibility analysis, the IFs were examined through a multicollinearity test, which is useful for selecting applicable IFs. The variance inflation (VIF) and tolerance (TOL) are two indicators of the level of multicollinearity [61,62]. A VIF value more prominent than or equivalent to 5 and a TOL value under 0.2 demonstrates a genuine multicollinearity issue [63,64]:

$$\mathrm{TOL}=1-R^2$$

(4)

$$\mathrm{VIF}=\frac{1}{\mathrm{TOL}}$$

(5)

where the coefficient of determination (R²) is the proportion of the variance in the target variable [65,66,67].

The quality of the predictions was examined using the receiver operating characteristic curve (ROC) technique. One of the popular methods for assessing model accuracy is using a receiver operating characteristic curve (ROC) in the landslide susceptibility assessment. The ROC is a useful method for describing the quality of probabilistic prediction systems [68]. The area under the curve can be utilized as a measurement to evaluate the general execution of the model [69], assuming that the bigger the area is, the better the presentation of the model becomes. A rough guide for grouping the exactness is the customary scholastic point framework. The following rankings were considered for accuracy test [70]: 0.90–1 (excellent), 0.80–0.90 (good), 0.70–0.80 (fair), 0.60–0.70 (poor) and 0.50–0.60 (fail).

In order to compare the reliability of the obtained landslide susceptibility maps, a difference in the three landslide susceptibility methods was evaluated in a 2 × 2 contingency table. The numbers of correctly predicted landslides and non-landslides are indicated by TP (True positives) and TN (True negatives), and the numbers of incorrectly predicted landslides and non-landslides are FN (False negatives) and FP (False positives), respectively. The aim of two-class prediction methods is to separate true cases from false ones. Based on the four descriptors, several further measures can be calculated. Sensitivity, also called true positive rate (TPR), and specificity (true negative rate, TNR) show the ratios of the landslide and non-landslide cases correctly identified by the landslide susceptibility models. Accuracy (ACC) refers to the correct classification of index. The mathematical basis of these parameters [71] is given as in equations (6)–(8):

$$Sensitivity=\frac{TP}{TP+FN}$$

(6)

$$Specificity=\frac{TN}{FP+TN}$$

(7)

$$ACC=\frac{TP+TN}{TP+FP+FN+TN}$$

(8)

The machine learning methods that were applied in the present study are briefly described below.

3.3.1. Random Forest (RF)

Random Forest (RF), introduced by Breiman in 2001 [72], has been demonstrated to be an incredible method in classification and regression. It has the ability to deal with high-dimension, continuous and categorical data. RF does not require an assumption about the statistical distribution of the data, and it can be considered robust with respect to changes in the composition of the dataset [73]. These properties are very useful in applications using variables that have mutual nonlinear interactions. In this work, we used the “ModelMap” package [74] in the R-statistical environment. The “ModelMap” package allows an interface between several existing R-packages to automate and construct the maps.

The training algorithm for RF applies the general technique of bootstrap aggregating or bagging to tree learners. The training set comprises a simple and fast way of learning a function that translates the dependent variable (target) X₀ and independent IFs data (X₁, X₂,…..,X_n) to output Y_i, where (X₁, X₂,…..,X_n) can be a mix of categorical and numerical variables. X₀ can be binary for classification, or numeric for regression. After training query point x, each tree independently predicts the average of the Y_i, as shown in Equations (9) and (10)

$$f_n^j(x)=\frac{1}{N^e(A_n(x))}{\displaystyle \sum _{\begin{array}{l}{Y}_i\in A_n(x)\\ I_i=e\end{array}}{Y}_i}$$

(9)

and the forest averages the predictions of each tree

$$f_n^{(M)}(x)=\frac{1}{M}{\displaystyle \sum _{j=1}^M{f}_n^j(x)}$$

(10)

where $A_n(x)$ denotes the leaf containing x, and $N^e(A_n(x))$ denotes the number of estimation points it contains.

3.3.2. Extreme Gradient Boosting (XGBoost)

The XGBoost is based on the Gradient Boosting concept. It uses a more regularized model formalization to control over-fitting, which gives it better performance. It was developed by Chen and Guestrin (2016) [75]. In XGBoost, training is done using an additive strategy. Given samples i with independent IFs data (X₁, X₂, …, X_n), and a tree ensemble model uses the K additive function to predict the outputs Y.

$$Y_i=\varphi (x_i)={\displaystyle \sum _{k=1}^K{f}_k(x_i),f_k\in }F$$

(11)

where F is the set of all possible trees. The $f_k$ is a function at each of the k steps of the descriptor values in $x_i$ to a certain output. The differentiable loss function that measures the difference between the predictions and the target is simply a mean-square error. After each step of boosting, the algorithm scales the newly added weights. In this study, an open source code “xgboost” for R-environment was utilized [76].

3.3.3. Deep Neural Network (DNN)

DNN is inspired by the information processing and communication pattern in the biological nervous system [77,78]. In DNN, layer representations are learned via models called a neural network. From a technical point of view, DNN algorithms are multilayered neural networks. The main body of DNN consists of an input layer that receives the data and an output layer that gives the prediction, and the core part is the hidden layer that processes the data. With the rapid development of DNN, a great number of libraries and packages has been developed to set a DNN with minimal effort. In this work, “H2O” packages [79] are designed for the R-statistical environment, and run seamlessly on both ‘CPU’ and ‘GPU’ devices.

Recently, DNN has been the subject of most of the advancements in artificial intelligence. A DNN is a multi-layered model for learning complex non-linear relationship between input data. For a given dependent variable i and known model parameters θ, the Poisson distribution of the response variable can be predicted, as given in Equation (12)

$${\lambda }_{m,\theta }\left(x\right)=\exp \left({\gamma }_i^T{a}_m^{N_{h,.}}\left(x,\theta \right)\right)$$

(12)

where $a_m^{N_{h,.}}\left(x,\theta \right)$ builds the shared part of the model. In this study, the ReLU activation function was chosen because it displayed empirically less optimization problems.

4. Results and Analysis

4.1. Selection of IFs

The 14 landslide IFs were analyzed for multicollinearity using TOL and VIF parameters. The results show that the highest value of VIF is 10.47, and the lowest value of TOL is 0.10, as presented in

Table 2. Multicollinearity test for the 14 landslide IFs.

Statistic	R²	Tolerance	VIF
Aspect	0.05	0.95	1.05
Elevation	0.40	0.60	1.67
Slope	0.70	0.30	3.33
Internal relief	0.90	0.10	10.47
Curvature	0.49	0.51	1.95
Drain prox. (h)	0.45	0.55	1.82
Drain prox. (v)	0.55	0.45	2.23
Drainage density	0.52	0.48	2.10
SPI	0.71	0.29	3.45
STI	0.78	0.22	4.55
TWI	0.74	0.26	3.85
Forest	0.25	0.75	1.33
Soil	0.02	0.98	1.02
Geology	0.09	0.91	1.09

Remark: Bold numbers are failed in multicollinearity test.

, which shows that IR has a multicollinearity problem among the 14 landslide IFs. So, IR was eliminated for further analysis.

4.2. Application of RF in Landslide Susceptibility Mapping

Three training parameters need to be defined in the RF model, as follows: ntree, the number of bootstrap samples for the original data (the default value is 500); mtry, the number of different predictors tested at each node, i.e., as many as the IFs; node size, the minimum size of the terminal nodes of the trees [72]. The elements not included in the ntree bootstrap sample are referred to as out-of-bag data (OOB) for that bootstrap sample. OOB is used to estimate an error rate called the OOB error rate (Figure 8a). The result suggests that when the training of the data was applied to building the model, the error rates stay constant after 100 trees.

The error rate was found to be 12%; therefore 88% of the model building result was accurate, which makes this a reasonably good model. Figure 8b depicts the variable importance by means of the mean decrease accuracy in percentage. The result shows that drainage proximity (h) has the highest importance followed by slope and forest.

4.3. Application of XGBost in Landslide Susceptibility Mapping

Basic training using XGBoost is the most critical part of the modeling. It contains several hyper-parameters to be set [80]. We hold the trees constant at the default of 1000 rounds. The learning rate was selected as 0.1 to avoid overfitting the model, and the maximum depth of the tree was chosen as 3. Figure 9a shows the log loss during the modeling process. After tuning the model, it was found that the log loss of the training process is 0.48 and the log loss for the testing data is 0.56 in 90 iterations, respectively. The model has the capability to check the important variables using the bias of all IFs. Figure 9b illustrates that drainage proximity (h) is the most important IF, followed by slope and curvature during the XGBoost modeling.

4.4. Application of DNN in Landslide Susceptibility Mapping

There is no universal rule for the determination of an optimum neural net structure. In this study, three hidden layers and 13 neurons in each hidden layer were used. The input IFs were normalized and then a “ReLU” activation function was used to run the model. The loss function is an important output of the model, which is used to measure the inconsistency between the dependent variable and the predicted outcome. Figure 10a presents the loss function during the modeling process. It has been found that as the number of epochs increases, the loss values for the training and test data decrease. After 10 epochs, the curves were stable. The relative importance of the model is presented in Figure 10b. Drainage proximity (h) is the most important IF, followed by slope and forest.

4.5. Evaluation Measures

Three landslide susceptibility maps obtained from RF, XGBoost and DNN were classified into five classes, namely very low, low, moderate, high and very high, using the Jenks natural breaks classification method, as presented in Figure 11a–c, respectively.

The outcomes show a common pattern, ensuring that with an increase in the index value, the occurrence of the probability of the landslide also increases, i.e., a higher index represents a higher susceptibility to landslide occurrence. The highly susceptible area is more confined near to the drainage and steep areas. The flat and gentle slopes represent the low susceptible areas.

The ROC evaluation for the three methods in the landslide susceptibility mapping shows that the area under the curve (AUC) values of the training data are 0.756 (fair), 0.757 (fair) and 0.855 (good) for RF, XGBoost and DNN, respectively (Figure 12a). Similarly, the AUC values of testing data are 0.703 (fair), 0.74 (fair) and 0.802 (good) for RF, XGBoost and DNN, respectively (Figure 12b).

The comparison of the three methods shows that the prediction capabilities of these models are different. Among the three models, DNN shows the highest prediction capability, followed by the XGBoost and then the RF model. The measurement of the landslide susceptibility models was validated utilizing the statistical evaluation parameters. The results during the testing stage for the implemented models are presented in

Table 3. Performance results of implemented models in the testing phase.

Data Set	Model	TP	FP	FN	TN	Specificity (%)	Sensitivity (%)	ACC (%)
Testing	RF	155	145	64	293	66.89	70.77	68.19
	XGBoost	159	106	60	332	75.79	72.6	74.73
	DNN	184	72	35	366	83.56	84.01	83.71

, which clearly shows that the execution of all models produced very good outcomes. During the testing phase, the correctly captured presence of landslides (TP) and absence of landslides (TN) is higher in the DNN model than in the RF and XGBoost models.

The DNN model exhibited the highest values of specificity (83.56%) and sensitivity (84.01%). The ACC of DNN was 83.71% (testing); on other hand, this value was 74.73% (testing) for XGBoost and 68.19% (testing) for the RF model, respectively. The performance results show that the prediction of the probable occurrence of a landslide by the DNN model surpasses the remaining models by a significant degree, followed by XGBoost and RF.

The consistent performance and model uncertainty were estimated using a four-fold plot. The quadruple plot is a visual representation of a confusion matrix that summarizes the number of TP, TN, FP and FN. Figure 13 shows the aggregated data of landslide presence or absence, i.e., the pixel frequencies, numerically in the corners of the plot. Thus, there were 219 landslides in the test catchment, of which the RF model correctly classified 155 (70.78%), compared with 438 non-landslides, of which 293 (66.89%) were correctly classified. In the case of XGBoost, 72.60% of the landslides were correctly classified, and 75.80% of the non-landslides were identified as no landslides pixels. Similarly, the DNN model identified 84.02% true landslides and 83.56% no landslides.

5. Discussion

The most effective means of reducing casualties and losses ensuing from landslides is landslide risk management; therefore, high-quality landslide susceptibility maps are an important tool [81] to predict a probable occurrence of the landslide. Landslide susceptibility mapping has been studied more via qualitative and quantitative methods. Previous studies focused on improving the performance of landslide prediction using various machine learning methods. The sample size used for modeling was not sufficient for these methods. Therefore, we used the Gaussian Kernel density function to select a double number of pixels without landslides. The purpose of the study was to apply and compare different machine learning models to obtain accurate landslide susceptibility maps. Therefore, in the present study, three classification algorithms (RF, XGBoost and DNN) were utilized and compared for landslide susceptibility mapping in the Inje area. The accuracy assessment of the results is the next step in comparing several models of landslide prediction. Without accuracy assessment, the results could not be interpreted, the methods could not be provided support or input, and they would have no scientific meaning. For this purpose, four metrics (sensitivity, specificity, ACC and AUC of ROC) were used in the evaluation of the model. Basically, the RF model is one of the more popular and widely used models to determine landslide susceptibility [73,82,83], while on other hand, XGBoost and DNN have been explored in limited studies for the spatial prediction of the occurrence of landslides, which showed promising results in the study area. An important step in landslide susceptibility modeling is determining the most significant IFs to decrease over-fitting. From multicollinearity analysis, IR was identified as a variable having collinearity issues; therefore IF was discarded in the modeling process. With regard to important IFs, three models revealed that horizontal drainage proximity (h) is the most important IF among 13 IFs. This finding justifies the statement given by Nobre et al. (2011) [84] that the terrain nearer to drainage has a strong correlation with the soil–water saturation regime, and it is a significant controlling factor for landslide occurrences, especially in the soil slopes. Gökceoglu and Aksoy (1996) [85] also state that a stream adversely affects the hillslope because most of the landslides occur in proximity to the stream. Similarly, the slope is second most important factor in these implemented models, as the slope is the most powerful agent of landslides [86]. On the other hand, SPI or TWI are less significant factors, compared to other IFs.

Although the applied machine learning models are “black box models”, applying machine learning is a field in which interest has been growing in recent years. The prescient outcomes do not distinguish the reasons for the landslide, however they do show the relation between landslides and terrain features. Nevertheless, such results may yield understanding of landslide events, and of which IFs affect landslides in the study area. To get more definitive results, more experiments should be done to validate the model with time series landslide inventories.

Previous studies illustrated that the RF method can adequately depict the probable occurrence of landslides and spatial data [82,83,87,88]. There are several pieces of research on landslide susceptibility that have been done using a gradient boosting algorithm [89,90,91]. It was found that the Boosted model performs better than the RF model. Recently, DNN has become popular among researchers [92,93]. Wang et al. (2017) [94] stated that the DNN framework achieved higher or comparable prediction accuracy, and this study also supports the previous findings [95,96,97]. Thus, DNN can give better predictions and can improve the understanding of landslide-susceptible areas.

6. Conclusions

Landslide susceptibility mapping is an important tool to predict the probability of occurrences of landslides in areas with hilly terrain. Therefore, high-quality landslide prediction models are of paramount importance. Two adjacent catchments were selected for the application and comparison of three machine learning models, i.e., RF, XGBoost and DNN. Deokjeokri catchment was selected as the training area and Karisanri catchment as the testing area. A landslide inventory consisting of 748 landslides in Deokjeokri catchment and 219 landslides in Karisanri catchment was prepared by using GIS. In total, 14 thematic layers were collected as landslide-influencing factors. Among these, IR failed in a multicollinearity test. From the comparison in this study, the three models have different accuracy levels. The results of DNN (training AUC = 0.855, testing AUC = 0.802) were significantly better than those of the XGBoost (training AUC = 0.757, testing AUC = 0.74) and RF (training AUC = 0.756, testing AUC = 0.703) models. The performance results of the testing phase also indicated that DNN has higher accuracy (83.71%) than XGBoost (74.73%) and RF (68.19%). The three models revealed that drainage proximity (h) is the most important IF. The DNN model is a very promising alternative for shallow landslide modeling in the study area. It is obvious that a high-quality landslide susceptibility map can decrease costs and time. Therefore, the resulting model can be used for spatial planning in the study area to reduce the risk of landslides and develop environment management.