Susceptibility Analysis of the Mt. Umyeon Landslide Area Using a Physical Slope Model and Probabilistic Method

: Every year, many countries carry out landslide susceptibility analyses to establish and manage countermeasures and reduce the damage caused by landslides. Because increases in the areas of landslides lead to new landslides, there is a growing need for landslide prediction to reduce such damage. Among the various methods for landslide susceptibility analysis, statistical methods require information about the landslide occurrence point. Meanwhile, analysis based on physical slope models can estimate stability by considering the slope characteristics, which can be applied based on information about the locations of landslides. Therefore, in this study, a probabilistic method based on a physical slope model was developed to analyze landslide susceptibility. To this end, an inﬁnite slope model was used as the physical slope model, and Monte Carlo simulation was applied based on landslide inventory including landslide locations, elevation, slope gradient, speciﬁc catchment area (SCA), soil thickness, unit weight, cohesion, friction angle, hydraulic conductivity, and rainfall intensity; deterministic analysis was also performed for the comparison. The Mt. Umyeon area, a representative case for urban landslides in South Korea where large scale human damage occurred in 2011, was selected for a case study. The landslide prediction rate and receiver operating characteristic (ROC) curve were used to estimate the prediction accuracy so that we could compare our approach to the deterministic analysis. The landslide prediction rate of the deterministic analysis was 81.55%; in the case of the Monte Carlo simulation, when the failure probabilities were set to 1%, 5%, and 10%, the landslide prediction rates were 95.15%, 91.26%, and 90.29%, respectively, which were higher than the rate of the deterministic analysis. Finally, according to the area under the curve of the ROC curve, the prediction accuracy of the probabilistic model was 73.32%, likely due to the variability and uncertainty in the input variables.


Introduction
Landslides cause substantial economic and social losses, especially in urban areas where many people live. Landslides are destructive and represent the most frequent risk factors in mountainous areas; especially in urban areas where damage to forests and infrastructure, such as buildings and

Study Area
The study area, Mt. Umyeon, is located in Seoul, the capital of the Republic of Korea, within 126°59′-127°01′ E, 37°27′-37°28′ N ( Figure 2). Mt. Umyeon is located at the center of a densely populated area, with highways to the east, rivers and parks to the south, and major cultural facilities around the Mt. Umyeon area. The study area, which was selected based on the type of watershed affected by rainfall, is approximately 22.18 km 2 (width: 5.075 km, length: 4.37 km). The maximum height of Mt. Umyeon is 293 m; the southern slope has a large slope and a valley, while the northern slope is gentle. Mt. Umyeon is a relatively low mountainous region of gneiss formed by retardation and weathering. The terrain is vulnerable to landslides because the gneiss in the bedrock is distributed with severe weathering and many faults. In addition, the dark veins are partially infiltrated, and overall, weathering is severe and the outcrop is poor [38].

Study Area
The study area, Mt. Umyeon, is located in Seoul, the capital of the Republic of Korea, within 126 • 59 -127 • 01 E, 37 • 27 -37 • 28 N (Figure 2). Mt. Umyeon is located at the center of a densely populated area, with highways to the east, rivers and parks to the south, and major cultural facilities around the Mt. Umyeon area. The study area, which was selected based on the type of watershed affected by rainfall, is approximately 22.18 km 2 (width: 5.075 km, length: 4.37 km). The maximum height of Mt. Umyeon is 293 m; the southern slope has a large slope and a valley, while the northern slope is gentle. Mt. Umyeon is a relatively low mountainous region of gneiss formed by retardation and weathering. The terrain is vulnerable to landslides because the gneiss in the bedrock is distributed with severe weathering and many faults. In addition, the dark veins are partially infiltrated, and overall, weathering is severe and the outcrop is poor [38]. In this region, a number of landslides occurred due to heavy rains, with a cumulative rainfall of 587.5 mm for 3 days from 26 to 28 July, 2011 [37]. About 150 large and small landslides occurred [38], and the area of debris flow was very wide compared to about 11 square kilometers selected as the radius of the study area ( Figure 3). Landslides were presumed to result from heavy rains over a period of about 1 h following the weakening of ground due to previous heavy rains of 230.0-266.5 mm from about 15 h before the landslide [38]; the estimated time of the landslide was 09:00 on 27 July, 2011. A number of landslides in the form of debris flow have been reported and landslide occurrence locations were collected based on field investigation and visual analysis of aerial photographs before and after the landslides [38] by points [3,5,13].  Table 1 shows the spatial datasets used in this study. To build the input dataset, all input data were constructed in a raster format with a 5 × 5 m grid structure using the ArcGIS 10.3 software. For the analysis of landslide susceptibilities, relevant factors were selected through literature review of previous studies [41][42][43][44][45][46]. First, a 1: 5000 topographic map was obtained, from which we collated data such as elevation, slope, and specific catchment area. A digital elevation model was constructed by extracting a contour vector layer, including elevation attributes, from a digital In this region, a number of landslides occurred due to heavy rains, with a cumulative rainfall of 587.5 mm for 3 days from 26 to 28 July, 2011 [37]. About 150 large and small landslides occurred [38], and the area of debris flow was very wide compared to about 11 square kilometers selected as the radius of the study area ( Figure 3). Landslides were presumed to result from heavy rains over a period of about 1 h following the weakening of ground due to previous heavy rains of 230.0-266.5 mm from about 15 h before the landslide [38]; the estimated time of the landslide was 09:00 on 27 July, 2011. A number of landslides in the form of debris flow have been reported and landslide occurrence locations were collected based on field investigation and visual analysis of aerial photographs before and after the landslides [38] by points [3,5,13]. In this region, a number of landslides occurred due to heavy rains, with a cumulative rainfall of 587.5 mm for 3 days from 26 to 28 July, 2011 [37]. About 150 large and small landslides occurred [38], and the area of debris flow was very wide compared to about 11 square kilometers selected as the radius of the study area ( Figure 3). Landslides were presumed to result from heavy rains over a period of about 1 h following the weakening of ground due to previous heavy rains of 230.0-266.5 mm from about 15 h before the landslide [38]; the estimated time of the landslide was 09:00 on 27 July, 2011. A number of landslides in the form of debris flow have been reported and landslide occurrence locations were collected based on field investigation and visual analysis of aerial photographs before and after the landslides [38] by points [3,5,13].  Table 1 shows the spatial datasets used in this study. To build the input dataset, all input data were constructed in a raster format with a 5 × 5 m grid structure using the ArcGIS 10.3 software. For the analysis of landslide susceptibilities, relevant factors were selected through literature review of previous studies [41][42][43][44][45][46]. First, a 1: 5000 topographic map was obtained, from which we collated data such as elevation, slope, and specific catchment area. A digital elevation model was constructed by extracting a contour vector layer, including elevation attributes, from a digital  Table 1 shows the spatial datasets used in this study. To build the input dataset, all input data were constructed in a raster format with a 5 × 5 m grid structure using the ArcGIS 10.3 software. For the analysis of landslide susceptibilities, relevant factors were selected through literature review of previous studies [41][42][43][44][45][46]. First, a 1: 5000 topographic map was obtained, from which we collated data such as elevation, slope, and specific catchment area. A digital elevation model was constructed by extracting a contour vector layer, including elevation attributes, from a digital topographic map Remote Sens. 2020, 12, 2663 5 of 17 ( Figure 4a). Finally, the slope and specific catchment area (SCA) were constructed based on the digital elevation model (Figure 4b,c). Second, the z-model was used to build the depth data of the fracture surface. The z-model is a submarine model that reflects the topographic characteristics of the slope by calculating its depth with respect to the elevation [48]. According to the results of ground drilling and seismic surveys by the Korean Society of Civil Engineers [37], the distribution of the soil thickness in the Mt. Umyeon area is 1.18-10.13 m. Therefore, the minimum and maximum depths of the soil in the z-model were set to 1 and 10 m, respectively, and soil thickness was calculated ( Figure 4d) [37].

Spatial Datasets
Third, unit weight, cohesion, friction angle, and hydraulic conductivity were obtained through field surveys with indoor experiments. Direct shear tests were conducted on four sampling points (SPs) on the obtained samples of cohesion and friction angle, and borehole shear tests and in-situ permeability tests were conducted at eight drilling points (DPs) in the study area ( Figure 5) [38]. Table 2 summarizes the results obtained from this process. In order to interpolate and analyze geological characteristics of the entire study area based on the obtained data [32,41], kriging spatial interpolation analysis from ArcGIS, a spatial processing interpolation technique that reflects the correlation between the distance from the surrounding value and the value located around it, was performed to construct data in raster format ( Figure 6).   In sequence, automatic weather system (AWS) data were used to construct rainfall intensity data. The Korea Meteorological Administration (KMA) operates over 510 AWS sites to monitor the atmospheric conditions near the ground in real time. To construct rainfall data in the study area, rainfall data was measured at 1-h intervals at seven AWS stations near Mt. Umyeon (Table 3). The rainfall intensity over the 16 h from 19:00 26 July to 09:00 27 July, 2011 was obtained to calculate the rainfall intensity at each AWS. Kriging interpolation was also applied to construct rainfall intensity data for the entire study area (Figure 7).  In sequence, automatic weather system (AWS) data were used to construct rainfall intensity data. The Korea Meteorological Administration (KMA) operates over 510 AWS sites to monitor the atmospheric conditions near the ground in real time. To construct rainfall data in the study area, rainfall data was measured at 1-h intervals at seven AWS stations near Mt. Umyeon (Table 3). The rainfall intensity over the 16 h from 19:00 26 July to 09:00 27 July, 2011 was obtained to calculate the rainfall intensity at each AWS. Kriging interpolation was also applied to construct rainfall intensity data for the entire study area (Figure 7).  Finally, a 2011 landslide occurrence map was prepared to validate the analysis. A total of 103 landslide locations were extracted by superimposing a 1: 5000 digital topography onto an aerial photograph with a spatial resolution of 50 cm before and after the landslide (Figure 8). Finally, a 2011 landslide occurrence map was prepared to validate the analysis. A total of 103 landslide locations were extracted by superimposing a 1: 5000 digital topography onto an aerial photograph with a spatial resolution of 50 cm before and after the landslide (Figure 8).

Physically Based Model
We used a physical slope model, an infinite slope, in which the depth of the ground is shorter than the length of the landslide. The infinite slope model is the most suitable model for GIS-based landslide analysis and is widely used to analyze landslides caused by rainfall [49][50][51][52]. The stability of the infinite slope model is expressed by the factor of safety (FS), which is the ratio of shear stress to shear strength. The shear stress is the stress that acts on the fracture surface, taking into account the weight and influence of the slope material, whereas the shear strength is based on the action of the pore water pressure on the vertical stress acting perpendicularly to the fracture surface.

Physically Based Model
We used a physical slope model, an infinite slope, in which the depth of the ground is shorter than the length of the landslide. The infinite slope model is the most suitable model for GIS-based landslide analysis and is widely used to analyze landslides caused by rainfall [49][50][51][52]. The stability of the infinite slope model is expressed by the factor of safety (FS), which is the ratio of shear stress to shear strength. The shear stress is the stress that acts on the fracture surface, taking into account the weight and influence of the slope material, whereas the shear strength is based on the action of the pore water pressure on the vertical stress acting perpendicularly to the fracture surface. Therefore, an FS value of less than 1 indicates a destructive state; and is expressed as [53]:  [36] and wet index of the infinite slope stability model, to calculate the groundwater level (z w ). SHALSTAB is a hydraulic model that considers only the ground flow while ignoring the outflow of the surface in the hydraulic model proposed by [54], which considers the flow and surface runoff in shallow ground. It is expressed as follows: where q is rainfall [m/day], a is the watershed area [m 2 ], b is the width of the contour line [m], and T is transmissivity [m 2 /day]. Furthermore, w is the wetness index, which is the ratio of the groundwater level to the depth of the slope (z w /D), the relative depth of the actual groundwater level with respect to the slope. The wetness index w is between 0 and 1; its maximum value is 1 because the groundwater level does not exceed the depth of the slope. The groundwater level is expressed as:

Monte Carlo Simulation
In order to calculate the safety factor through a mathematical model, input data such as shear strength (cohesion force and friction angle) of the ground or groundwater level are required. However, these data are usually obtained through limited field surveys or indoor experiments, which is absolutely insufficient compared to the size of a wide area. Uncertainty intervenes in these data [55]. In this study, the probabilistic analysis technique was applied to quantitatively consider uncertainty and reflect it in the analysis. Such methods can be used to estimate the probability that the value of a state function satisfies a threshold by assuming that the input variable is randomly selected from a specific distribution. The probabilistic analysis method replaces the safety factor to evaluate the risk of slope using the probability of failure and is evaluated as the most effective method to quantify uncertainty among the various techniques proposed so far.
Probabilistic methods include MC simulation, first-order second moment (FOSM), the point estimate method (PEM), etc. FOSM and PEM are approximate methods, so their accuracy is relatively low compared to MC simulation [31,56]. MC simulation can represent the variability of input variables by randomly generating input variables, and is suitable for analyzing one or more random variables [57]. In this study, analysis was performed using Monte Carlo simulation among the probabilistic analysis techniques. The MC simulation method first determines the model of the state function. The state function used in landslide susceptibility analysis depends on FS; and we calculate the probability of failure, which is the probability that FS is less than 1. Second, the probability distribution of the input variables is calculated. Because information about the probability distribution of input variables is not available, this is considered to be a normal distribution with an appropriate average and standard deviation [28,32,33,58]. Third, N random values were extracted from the distribution of input variables, and N values of FS were calculated by substituting the values of the randomly generated input variables. Finally, the probability of failure, i.e., the proportion of the N FS values that are less than 1, was calculated as follows:

Applications and Validation
In our MC simulation, the cohesion and frictional angle were considered as random variables. According to previous studies, the probability distributions of geotechnical characteristics follow a normal or lognormal distribution, so the probability distribution of the input variables was assumed to be normal in this study [59][60][61]. Since the mean and standard deviation are required to define the distribution of input variables, the average was calculated based on the constructed data, and the standard deviation was calculated from the coefficient of variation. In previous studies, the ranges of the coefficient of variation of the internal friction angle and cohesion were 10-20% and 25-30%, respectively [56,62,63]. Because the entire study area was composed of gneiss and thus the geological characteristics were relatively similar, the minimum value of each coefficient of variation was assumed to be 10% for the internal friction angle and 25% for the cohesive coefficient, and the number of repetitions was set to 100,000.
The probability of failure is established by probabilistic methods, in contrast to the deterministic analysis in which we interpret an area to be unstable when FS is less than 1 [64]. Therefore, based on the results of previous studies, landslide-susceptible areas were classified based on failure probabilities of 1%, 5%, and 10% [64]. We then calculated the landslide prediction rate, which is the ratio of the number of landslides in a landslide-susceptible area. Furthermore, we carried out a deterministic analysis of the same dataset and calculated the FS by taking a simple average of the data. Finally, we compared the results from the MC simulation to those obtained using deterministic techniques.
Finally, the area under the curve (AUC) of the receiver operating characteristic (ROC) was calculated for validation purposes. The x-axis of the ROC graph is the ratio of the expected landslide area where it shows high susceptibility of landslides, and the y-axis is the probability of landslides. The AUC is calculated from the area under the ROC graph and has values between 0 and 1. The closer it is to 1 (100%), the higher the accuracy [45].

Result
The results of the MC simulation are summarized in Table 4 and shown in Figure 9. When a 1% probability of failure was considered to indicate a susceptible area, 57.75% and 42.25% of the study area was found to be unstable and stable, respectively, and the landslide prediction rate was 95.15%. When the susceptible area was set based on a probability of failure of 5%, the unstable area was 54.23%, the stable area was 45.77%, and the landslide prediction rate was 91.26%. When we defined the susceptible areas based on a probability of failure of 10%, the unstable area was 52.02%, the stable area was 47.98%, and the landslide prediction rate was 90.29%. By contrast, according to the deterministic analysis method, the unstable area was 42.72%, the stable area was 57.28%, and the prediction rate, which in this case indicated the proportion of the predicted landslides that occurred, was 81.55%. Additionally, the AUC calculation of the ROC graph of the MC simulation was 0.7332 (73.32%) (Figure 10). Remote Sens. 2020, 12, x FOR PEER REVIEW 13 of 19   Based on the results of both analytical methods, we conclude that the stable area occupies more than 40% of the study area. According to the results of the MC simulation, as the reference probability of failure decreases, the proportion of unstable regions increases with increasing landslide prediction rate. Our MC simulation predicted relatively high proportions of unstable areas and landslide prediction rates compared to the deterministic analysis. An AUC value of 0.7 or more can be interpreted as indicating good predictive performance [65,66]. The AUC obtained from the MC simulation (0.7332) was greater than 0.7, which could explain the high landslide prediction rate. The good predictive performance of the MC simulation is attributed to the variability of the input variables.  Based on the results of both analytical methods, we conclude that the stable area occupies more than 40% of the study area. According to the results of the MC simulation, as the reference probability of failure decreases, the proportion of unstable regions increases with increasing landslide prediction rate. Our MC simulation predicted relatively high proportions of unstable areas and landslide prediction rates compared to the deterministic analysis. An AUC value of 0.7 or more can be interpreted as indicating good predictive performance [65,66]. The AUC obtained from the MC simulation (0.7332) was greater than 0.7, which could explain the high landslide prediction rate. The good predictive performance of the MC simulation is attributed to the variability of the input variables.

Discussion
In this study, we conducted a landslide susceptibility analysis of the Mt. Umyeon area, a representative example of urban landslide damage. An interpolation method was used to construct data for the entire study area because we could only obtain subsurface data and geotechnical characteristics for a few points. In addition, it was infeasible to obtain sufficient data to infer the probability distribution due to the cost and constraints of acquiring the actual property values, so we assumed the probability distribution to be normal. Since it was difficult to obtain a sufficient amount of field survey data compared to the scope of the study area, but it was difficult to apply the unsaturated soil theory, which requires many experimental values and well reflects the behavior of the groundwater level of the actual ground; saturated soil theory with relatively simple interpretation was applied in this study. In addition, several hydraulic parameters were assumed based on existing studies. Therefore, a more accurate groundwater model can be constructed with sufficient data on hydraulic parameters obtained through experiments. Future studies should increase the number of data acquisition points and experiments at each point to compensate for these limitations.
The results of the landslide susceptibility analysis of the Mt. Umyeon area can be summarized as follows. The MC simulation results showed that the landslide prediction rate (90.25, 91.26, and 95.15%) was significantly higher than the landslide prediction rate (81.55%) of the deterministic method. Comparing the results of the deterministic analysis method and the MC simulation, which is a probabilistic analysis method, by using the coincidence ratio with the location coordinate of the landslide, the landslide prediction rate of 8.7% to 13.6% is higher than the deterministic analysis method. Through these results, the probability analysis method can be applied considering the variability of the rainfall intensity, so accuracy can be improved by considering the uncertainty inherent in the rainfall intensity. MC simulation result also yielded a high AUC value (0.7332). We attributed this improvement in results to the fact that MC simulation can include the variability and uncertainty of the input variables. Moreover, the proportion of unstable areas around Mt. Umyeon exceeded 40%. This is thought to be due to the weathered gneiss distributed throughout the area and the high rainfall intensity.

Conclusions
Following the Mt. Umyeon landslide, a representative landslide case in the Seoul metropolitan city, South Korea, it has become necessary to analyze landslide susceptibility. The purpose of this study is to analyze the susceptibility of landslide disasters considering uncertainty before the occurrence of widespread landslides. Analysis based on physical slope models can be used to evaluate landslide susceptibility in areas without prior landslide occurrence. A spatial database of Mt. Umyeon, the study area, was constructed and analyzed by applying it to the infinite slope model that is similar to the characteristics of landslide occurrence in the study area. In addition, landslide susceptibility was applied and analyzed based on a physically based model and MC simulation, a probabilistic analysis technique considering uncertainty. An infinite slope model was used as the physically based model. In the GIS environment, it was possible to analyze landslide susceptibility in the study area, and by using an infinite slope model and using probabilistic techniques, it was possible to evaluate the landslide susceptibility as a quantitative indicator of probability of failure. Furthermore, to evaluate the accuracy of our landslide predictions, we applied a deterministic method and compared the results. Finally, the accuracy of the landslide prediction was calculated based on the AUC.
This study confirms that it is possible to evaluate high-accuracy landslide susceptibility without prior information of landslide occurrences by combining a physical slope model with probabilistic method. By varying the reference probability of failure from 1% to 10% in the MC simulation, it was possible to adjust the safety level as needed. This means that the reference failure probability can be varied according to the purpose of analysis. For example, a susceptibility map with a high standard probability of failure should be used to determine the locations of disaster prevention structures to minimize costs. Conversely, if the danger zone is temporarily set to minimize human damage, the susceptibility map should be based on the minimum probability of failure. In this way, the same dataset and probabilistic technique can be used for different purposes.
To prevent a repeat of the damage incurred by the Mt. Umyeon landslide, it is necessary to carry out landslide susceptibility studies of areas where landslides have not occurred. In particular, prior landslide susceptibility analysis should be carried out in areas with high population densities to minimize large-scale damage. The methodology presented herein can be used to prepare measures to reduce