Soil Depth Prediction Model Using Terrain Attributes in Gangwon-do, South Korea

Abstract

Soil depth is a crucial parameter in slope stability analysis in mountainous areas. The drilling survey is the most reliable method for determining soil depth, but it requires a high cost for the vast geographical area. Therefore, this study proposes a soil depth prediction model for mountainous areas that uses Terrain Attributes (TAs) from digital maps. Gangwon-Do, a predominantly mountainous region in South Korea, is selected as the study target area. The study area is classified by parent rock type into igneous rocks, metamorphic rocks, and sedimentary rocks. The correlation with TAs is analyzed through multi-collinearity using drilling data published in the Korea drilling information database. In addition, the most suitable combination of variables is selected through multi-collinearity analysis, and the regression model using STI, TWI, and SLOPE is found to be the most appropriate model (VIF < 10). The proposed model for soil depth shows significance at p < 0.001, and the correlation coefficient ($R^2$) is figured out for igneous rock (0.702), metamorphic rock (0.686), and sedimentary rock (0.693). In addition, the reliability of the proposed model was verified by using data from regions not included in the model development, and the correlation coefficients were igneous rock (0.867), metamorphic rock (0.801), and sedimentary rock (0.814). The model proposed is more suitable for Korean topography than the existing statistical models; it can help to increase the accuracy of slope stability analysis.

1. Introduction

Soil depth is defined as the vertical distance from the ground surface to bedrock [1,2,3,4]. The soil depth can be changed by erosion and weathering, and it is determined by various environmental variables, including slope, land use, curvature, parent material, weathering rate, climate, vegetation cover, upslope contributing area, and lithology [5,6]. In addition, the weathering depth varies depending on the rock type [7], and the distribution of soil depth will be different accordingly.

Soil depth plays an important role in landslide disaster prevention and management. The increasingly complex shallow landslide generation process integrates with physically based models such as SHALSTAB [8], SHETRAN [9], GEOtop FS+ [10], TRIGRS [11], and H-slider [12] to predict spatial patterns of landslide occurrences [12,13,14,15]. However, most of the analysis methods using landslide analysis models are carried out under the assumption that the soil depth is constant without consideration of the spatial distribution for soil depth [8,12,13,16,17,18]. Additionally, the flow of groundwater and the formation of wetting front, which have a significant effect on the stability of the slope, are also significantly affected by the soil depth and its spatial distribution [19,20].

Many researchers have conducted studies on how to predict soil depth, and these studies can be classified broadly into three groups: (1) physically based, (2) interpolation-based, and (3) regression methods.

The physically based model for soil depth prediction was developed as a landscape evolution model to solve the soil mass-balance equation over time [5]. Most landscape evolution models use linear sediment transport laws, which assume that the terrain slope and the sediment amount are equal. Although it develops into a simple analytical solution for mass balance, observational results report that landslides contribute to soil transport [5,21,22,23,24]. Recently, non-linear transport laws have been used to predict slope failure [23,25,26] and to describe the planar topography of steep slopes [27].

Interpolation-based soil depth prediction methods include IDW [28], Kriging [29,30], Co-Kriging [31], regression-kriging, and so on [3,32]. These days, many researchers have suggested the interpolation method to predict the spatial distribution of the soil depth, and this method often shows high accuracy compared to other methods. Penížek and Borůvka [33] reported that the co-kriging method showed higher accuracy than the ordinary kriging (OK), regression kriging (RK), and linear regression methods. However, although the interpolation method is convenient to use, it requires a large amount of data [5,34,35], and it cannot be said to show high accuracy in all cases because the deviation is significant due to the various environmental variables [36].

Regression methods correlate soil depth with TAs from digital elevation maps (DEMs) [33,37,38,39]. Goodman [40] reported that the coefficient of determination ($R^2$) of the slope angle, absolute height, and relative height with soil depth were 0.46, 0.36, and 0.37, respectively. Gessler et al. [41] also analyzed the statistical correlation between soil depth and plan curvature, compound topographic index, and upslope mean plan curvature. Additionally, Qiyoung et al. [4] suggested that the correlation coefficient between the regression model using slope and TWI as variables and soil depth was 64.1%, and Mehnatkesh et al. [2] reported that the correlation coefficient between slope, wetting index, catchment area, and STI with soil depth was 76%.

In this paper, the study area was classified into igneous rock, metamorphic rock, and sedimentary rock, which are representative rock type in Korea. The correlation between TAs and soil depth in the study area is analyzed. Additionally, Using the multiple linear regression method, we propose a soil depth prediction model for each rock type that is suitable for the geological condition in Korea.

2. Methodology

2.1. Study Area and Terrain Attributes

The study area is selected as Gangwon-do in South Korea (Figure 1), with a total area of 16,827.93 km² and a population of 1,539,521 (population density of 90 people/km²). Gangwon-do has the most mountainous region and the lowest population density in Korea [42]. To estimate the terrain attributes of the study area, a digital map (scale 1:25,000) of data from the National Geographic Information Institute [43] is used. Terrain attributes are derived using ESRI’s ArcGIS (ArcMap 10.1). Definitions and results of each derived terrain attribute are shown in

Table 1. Definition of Terrain Attributes in the study area and its derived results.

Terrain Attributes	Define
	The SLOPE is defined as the angle between the tangent to the terrain surface and the horizontal plane and is a variable that determines the flow velocity by gravity [44]. $\beta =\arctan \left(\sqrt{p^2+q^2}\right)$$\mathrm{where}p=\frac{\partial z}{\partial x}, q=\frac{\partial z}{\partial y}$
	Topographic Wetness Index (TWI) is an index of the water content in the ground, and it can determine the spatial pattern of water flow [45,46]. TWI was developed by Beven and Kirkby [47] in TOPMODEL and is defined as a function of the upstream contributing area per unit area orthogonal to the slope and the flow direction. $\ln \left(\frac{\alpha }{\tan \beta }\right)$ Where α is the slope area of the unit grid length, β is the slope at a point on the surface.
	Sediment Transport Index (STI) was proposed by Moore and Burch [48], and is an indicator of the erosion and sedimentation process. It is defined by the non-linear equation of the specific catchment area and slope. $LS={\left(\frac{A_S}{22.13}\right)}^{0.6}{\left(\frac{\sin \beta }{0.0896}\right)}^{1.3}$ where$A_S$ is a specific catchment area, β is the slope at a point on the surface.
	Stream Power Index (SPI) is an index indicating the erosion risk of potential flows on a terrain surface. When the catchment area and slope increase, the SPI also increases as the flow velocity due to gravity increases [45,46] The SPI is defined by the following equation. $SPI=\ln \left(A_s\times \tan \beta \right)$
	Curvature has an important effect on landslides and is classified into three types: convex, plane, and concave. In general, it is considered to have a negative effect on the landslide due to ponding caused by more extreme rainfall penetration in the case of the concave curvature relative to that of the Convex [49].
	Specific catchment area (SCA) is defined as the value obtained by dividing the upslope area by the contour width and is a value commonly used in hydrology to analyze the flow of water on the slope [50]. $A_s=\underset{w\to 0}{\lim }\left(\frac{A}{w}\right)$

.

2.2. Statistics Analysis

2.2.1. Correlation Analysis

Public drilling data from the Geotechnical Information DB System is used to acquire soil depth data of the study area. This database includes drilling information for 332,889 holes in Korea [51]. In order to secure the reliability of variables which are overly sensitive to the streamline among TAs, data near the streamline were excluded. In addition, too old data (produced before 2000) or data, in which location information of digital maps and boring test data is inconsistent were excluded, since it can be considered an outlier.

Yoon [41] used the soil depth ratio (SR = soil depth/slope height) derived by the slope height for representative rock types in Korea for a total of 373 slopes. It was evaluated and reported that the average soil depth ratio (SR) was 0.09 for sedimentary rocks, 0.17 for metamorphic rocks, and 0.3 for igneous rocks. As such, the weathering depth will be different depending on the rock types, and the soil depth will also be different for this reason. Therefore, in this study, the Boring Test data in study area was classified into igneous, metamorphic, and sedimentary rocks, which are representative rock types in Korea, and the correlation between soil depth and TAs are analyzed. After the excluding outliners, a total of 297 sites were obtained. They were classified into 101 sites for igneous rock, 101sites for metamorphic rock, and 95 sites for sedimentary rock. Descriptive statistics for soil depth by representative rock types are shown in Figure 2. The used soil depth data range from 0.2 to 0.68 m, and the average depth ranges from 2.63 to 3.03 m. In particular, the sedimentary rock data is distributed relatively wide compared to the data for the other two types of rock.

Correlation analysis is performed between the soil depth of the study area and the calculated TAs. The correlation analysis results for each of the representative rock types are shown in Figure 3, Figure 4 and Figure 5 in the form of a correlation matrix. It shows the histograms of each parameter along the central diagonal, the scatter diagram for each pair of parameters in the lower triangle, and each Pearson correlation coefficient and significance level (p-value) is shown in the upper triangle. The results of the correlation analysis indicate that the correlation between soil depth and SLOPE is the highest in all three rock types, and the correlation between Curvature and SCA (Specific catchment area) is relatively low.

2.2.2. Multi-Collinearity Analysis

When there is a high correlation between the independent variables in the multi-regression calculation, the problems of multi-collinearity may occur, which may decrease reliability of the regression coefficient estimates. Multi-collinearity can be evaluated through a variance inflation factor (VIF) such as that given in Equation (1).

$$\mathrm{Tolerance}=\left(1-{\mathrm{R}}^2\right),\mathrm{VIF}={\mathrm{Tolerance}}^{-1}$$

(1)

The larger the VIF, the smaller the tolerance, and the more likely it is multicollinear. If the VIF is greater than 5, it is judged that the possibility of multi-collinearity is high. In general, if the VIF is greater than 10, it is determined that multi-collinearity is problematic [52].

Multi-collinearity analysis is performed on six cases consisting of linear regression expressions of six TAs, which are classified into three types of representative rock (

Table 2. The result of multi-collinearity analysis.

Case No.	Parameter	Igneous Rock		Metamorphic Rock		Sedimentary Rock
		Tolerance	VIF	Tolerance	VIF	Tolerance	VIF
1	SLOPE	0.012	82.070	0.010	96.059	0.270	3.706
	TWI	0.002	591.425	0.002	523.360	0.021	46.671
	STI *	0.001	1166.374	0.001	1089.255	0.027	37.516
	SPI	0.001	1836.403	0.000	2307.441	0.061	16.270
	CURV	0.731	1.368	0.700	1.429	0.937	1.068
	SCA **	0.015	66.595	0.013	76.471	0.032	30.870
2	SLOPE	0.121	8.287	0.142	7.048	0.270	3.699
	TWI	0.016	63.426	0.018	56.247	0.027	37.110
	STI *	0.062	16.075	0.068	14.708	0.055	18.134
	CURV	0.738	1.355	0.740	1.352	0.941	1.063
	SCA **	0.015	66.008	0.014	72.581	0.032	30.814
3	SLOPE	0.141	7.075	0.167	5.975	0.271	3.695
	TWI	0.266	3.763	0.230	4.351	0.538	1.859
	STI *	0.239	4.181	0.221	4.516	0.384	2.604
	CURV	0.779	1.284	0.749	1.335	0.961	1.040
4	SLOPE	0.141	7.074	0.172	5.809	0.272	3.678
	TWI	0.282	3.545	0.254	3.939	0.545	1.833
	STI *	0.244	4.107	0.221	4.515	0.389	2.573
5	SLOPE	0.577	1.733	0.639	1.566	0.556	1.798
	TWI	0.577	1.733	0.639	1.566	0.556	1.798
6	SLOPE	1.000	1.000	1.000	1.000	1.000	1.000

* ln(STI), ** ln(SCA).

). Table 2 shows the results of the analysis by sequentially removing variables that cause multi-collinearity problems. As a result of the analysis, there is a multi-collinearity problem in Case No. 1 and 2, and Case No. 3–6 have a VIF of 4 or less. Thus, there is no multi-collinearity problem.

2.2.3. Multi-Linear Regression Analysis

Multiple regression analysis is one of the statistical methods used to establish and model the relationship between two or more independent variables for one dependent variable. If the data are linearly related, the multi-linear regression (MLR) of the observation $\text{'}y\text{'}$ for ‘$n$’ independent variables ‘$x$’ is as follows:

$$y=\underset{\_}{x}·\underset{\_}{\beta }+\epsilon ={\beta }_0+{\beta }_1{x}_1\cdots {\beta }_n{x}_n+\epsilon$$

(2)

where ${\beta }_0$ represents the intercept, and ${\beta }_1\cdots {\beta }_n$ represent the regression coefficients. $\epsilon$ is an error that is not explained by the model and is assumed to follow a normal distribution with a mean of zero and a variance of ${\sigma }^2$. Estimates of regression coefficients can be obtained by $\widehat{\beta }={\left({\underset{\_}{x}}^T\underset{\_}{x}\right)}^{-1}{\underset{\_}{x}}^{T}y$.

The selection of appropriate independent variables for multi-regression analysis is one of the critical issues. Generally, the methods of selecting independent variables include forward selection, backward elimination, stepwise method, and so on. However, the regression equation does not necessarily include the selected variables. If multi-collinearity exists between the independent variables, another model can be selected without a separate analysis.

Helsel and Hirsch [53] described the existing selection method and its shortcomings. They stated that there are many advantages of choosing a suitable model by evaluating the combination of statistics used in model selection and independent variables used in the model.

In this study, the adjusted coefficient of determination ($R_{adj}^2$) and RMSE are used to determine the number of variables to be used in the soil depth prediction model.

The coefficient of determination is a statistic that indicates the model’s suitability for a given dataset. Therefore, the large value of the coefficient of determination means that the proportion of changes that the independent variable can explain is large. However, the coefficient of determination may be increased by some independent variables having high explanatory power, even if some of the independent variables in the model are not very descriptive. There is a correction coefficient as a statistic that complements these shortcomings of the adjusted coefficient of determination ($R_{adj}^2$) and is defined as follows:

$$R_{adj}^2=1-\left[\left(1-R^2\right)\left(n-1\right)/\left(n-k-1\right)\right]$$

(3)

where n is the number of points in the data sample, k is the number of independent variables included in the regression equation and is the number of variables in the model excluding the constant. These $R_{adj}^2$can be used to compare the suitability of the model [54].

Additionally, it is the most widely used method through RMSE (root mean squared residual) for error analysis of regression models, and the definition is given as Equation (4).

$$RMSE=\sqrt{{\displaystyle \sum }_{i=1}^n{\left({\left(Y_m\right)}_i-{\left(Y_p\right)}_i\right)}^2/n}$$

(4)

where $Y_m$ is the measurement value, $Y_p$ is the predicted value, $n$ is the number of datapoints.

To determine a suitable model, $R^2$, $R_{adj}^2$, RMSE for each case are calculated as

Table 3. The results of the multi-collinearity analysis.

Case No.	Igneous Rock			Metamorphic Rock			Sedimentary Rock
	${\mathit{R}}^2$	${\mathit{R}}_{\mathit{a}\mathit{d}\mathit{j}}^2$	RMSE	${\mathit{R}}^2$	${\mathit{R}}_{\mathit{a}\mathit{d}\mathit{j}}^2$	RMSE	${\mathit{R}}^2$	${\mathit{R}}_{\mathit{a}\mathit{d}\mathit{j}}^2$	RMSE
1	0.707	0.688	0.863	0.697	0.677	0.773	0.856	0.846	0.672
2	0.704	0.689	0.866	0.691	0.674	0.781	0.698	0.681	0.973
3	0.704	0.691	0.868	0.690	0.678	0.781	0.695	0.682	0.977
4	0.702	0.693	0.870	0.686	0.676	0.876	0.693	0.683	0.981
5	0.665	0.658	0.922	0.666	0.659	0.864	0.621	0.613	1.09
6	0.607	0.603	0.999	0.609	0.605	0.878	0.604	0.600	1.114

Case 1: SLOPE, TWI, STI*, SPI, CURV, SCA**, Case 2: SLOPE, TWI, STI*, CURV, SCA**, Case 3: SLOPE, TWI, STI*, CURV, Case 4: SLOPE, TWI, STI*, Case 5: SLOPE, TWI, Case 6: SLOPE.

, and the results are shown in Figure 6

(1)

As the number of variables used in the model decreases, $R_{adj}^2$ decreases: Especially, a large decrease in Cases No. 5 and 6.

(2)

Since Cases No. 1 and 2 have multi-collinearity problems, suitable models are Cases No. 3 and 4.

(3)

Case No. 4, compared to Case No. 3, has a relatively large error.

(4)

Applied Model in this study is Case No. 4, since it uses fewer variables than in Case No. 3.

3. Predicted Soil Depth Model

Previous studies to predict soil depth using TAs have been attempted continuously, and as a result, statistical models for soil depth prediction have been proposed. Representative statistical models for soil depth prediction are summarized in

Table 4. Previously proposed statistical soil depth prediction models.

Proposer	Model	Range of Soil Thicknesses (cm)	Number of Datapoints
Penizek [33]	$88.8+1.020\mathrm{AS}+0.057\mathrm{AL}-2.491\mathrm{S}$	40–160	553
Gessler [41]	$-57.95+12.83{\mathrm{CV}}_{\mathrm{PL}}+21.46\mathrm{TWI}$	0–200	30
Qiyong [4]	$97.883-1.457\mathrm{S}+2.688\mathrm{TWI}-0.073\mathrm{AL}$	3.1–198.4	137
Han [55]	$17.918\mathrm{TWI}+0.550$(Cell size = 10m)	30–150	79
Mehnatkesh [2]	$122.13-0.11\mathrm{S}+0.012\mathrm{TWI}+0.012\mathrm{CA}-0.23\mathrm{STI}$	30–150	100

AS: Aspect, AL: Altitude, S: Slope,${\mathrm{CV}}_{\mathrm{PL}}$: Plane Curvature.

.

The soil depth prediction model derived in this study is shown in

Table 5. The result of soil depth prediction models in this study.

Rock Type	Model	p-Value	Number of Data
Igneous	$0.626-0.536\mathrm{STI}+0.509\mathrm{TWI}-0.021\mathrm{SLOPE}$	<0.001	101
Metamorphic	$4.442-0.702\mathrm{STI}+0.228\mathrm{TWI}-0.034\mathrm{SLOPE}$	<0.001	101
Sedimentary	$4.289-0.667\mathrm{STI}+0.242\mathrm{TWI}-0.041\mathrm{SLOPE}$	<0.001	95

, and the soil depth prediction model results for the study area are shown in Figure 7a. In addition, the scatter plot for each representative rock type is shown in Figure 7b–d. $R_{adj}^2$ for representative rock types are 0.698 (Igneous rock), 0.676 (Metamorphic rock), 0.683 (Sedimentary rock), and the RMSE for representative rock types are 0.870 (Igneous rock), 0.876 (Metamorphic rock), and 0.981 (Sedimentary rock).

For the verification of the developed regression equation, the measurement data and regression equation are compared with the data for the three regions not included in the development of the regression model. The soil depth distribution of the target area is prepared (

Table 6. Results of soil depth prediction models at the verification site.

Rock Type	R	${\mathit{R}}^2$	${\mathit{R}}_{\mathit{a}\mathit{d}\mathit{j}}^2$	RMSE	p-Value	Number of Datapoints
Igneous	0.931	0.867	0.859	0.724	<0.001	18
Metamorphic	0.895	0.801	0.794	1.104	<0.001	30
Sedimentary	0.902	0.814	0.807	0.288	<0.001	29

). The numbers of data in igneous, metamorphic, and sedimentary rock areas for verification are 18, 30, and 29 sites, respectively, and the determinants between regression and instrumentation data are 0.867, 0.801, and 0.814, respectively. Additionally, $R_{adj}^2$ for representative rock types are 0.859 (Igneous rock), 0.794 (Metamorphic rock), 0.807 (Sedimentary rock), and RMSE for representative rock types are 0.724 (Igneous rock), 1.104 (Metamorphic rock), 0.288 (Sedimentary rock).

The distribution of soil depth for the entire area to be verified (Figure 8a,c,e) and the model verification results for each rock distribution (Figure 8b,d,f) are shown in the below figures. As a result of the verification, the $R_{adj}^2$is from 0.794 to 0.859, which is judged to be an appropriately predictable model. However, RMSE is relatively large in metamorphic rock regions. It is considered that this is because the soil depth data of the metamorphic rock verification in the target area are too concentrated near the valley, and a significant amount of the data includes values measured inside the valley.

4. Conclusions

This paper proposes a statistical model that can predict soil depth, an essential factor in landslide disasters. Since erosion and weathering characteristics depend on the rock type, a soil depth prediction model for each rock type was developed: igneous rock, metamorphic rock, and sedimentary rock.

-

The regression model uses open data provided by the Geotechnical Information DB System; a total of 297 sites were obtained. They were classified into 101 sites for igneous rock, 101 for metamorphic rock, and 95 for sedimentary rock.

-

As a result of analyzing the correlation between the six TAs obtained from the numerical map and soil depth, the variables with the highest correlation are SLOPE, and curvature and SCA are found to have relatively low correlation. In addition, a model using three variables (SLOPE, STI, TWI) is determined from $R_{adj}^2$, and RMSE values for multi-collinearity analysis and the combination of six cases for variables.

-

For the models of igneous rock, metamorphic rock, and sedimentary rock, the $R_{adj}^2$ values are 0.698, 0.676, and 0.683, respectively, and the RMSE values are 0.870, 0.876, and 0.981. Additionally, Verification sites used data from 18 igneous rock sites, 37 metamorphic rock sites and 30 sedimentary rock sites. The $R_{adj}^2$ values are 0.859, 0.794, 0.807, and the RMSEs are 0.724, 1.104, 0.288 for igneous, metamorphic, and sedimentary rocks, respectively.

Therefore, it is expected that the three types of soil depth prediction models proposed will have better applicability and utilization for soil in Korea. However, the data used in this paper are open data provided by the Geotechnical Information DB System and are rather less in number to represent the wide study area. In addition, it should be noted that the relationship between TAs and soil depth used for model development was developed, assuming a linear correlation.