Stage Division of Landslide Deformation and Prediction of Critical Sliding Based on Inverse Logistic Function

Abstract

The cumulative displacement-time curve is the most common and direct method used to predict the deformation trends of landslides and divide the deformation stages. A new method based on the inverse logistic function considering inverse distance weighting (IDW) is proposed to predict the displacement of landslides, and the quantitative standards of dividing the deformation stages and determining the critical sliding time are put forward. The proposed method is applied in some landslide cases according to the displacement monitoring data and shows that the new method is effective. Moreover, long-term displacement predictions are applied in two landslides. Finally, summarized with the application in other landslide cases, the value of displacement acceleration, 0.9 mm/day², is suggested as the first early warning standard of sliding, and the fitting function of the acceleration rate with the volume or length of landslide can be considered the secondary critical threshold function of landslide failure.

1. Introduction

Judging the deformation stages and predicting the critical sliding time of a landslide is very significant to determine how and when the landslide will fail. The most common method of dividing the deformation stages and predicting the trend of landslides is generally based on the cumulative displacement-time curve [1,2]. From the perspective of the prediction time length, it can be divided into a long-term prediction (1 to 3 years) [3,4,5,6,7], medium-term prediction (3 to 12 months), short-term prediction (1 to 3 months) and critical sliding prediction (1 to 10 days) [8].

Saito [1,9] first proposed a landslide prediction model based on the three-stage creep theory of rock and soil materials and successfully predicted the Gaochangshan landslide in Japan. Since then, many different landslide prediction models have been developed [10,11,12,13,14,15,16,17,18]. Over the last decades, Wang and Nie [19] combined nonlinear regression analysis with the quadratic curve exponential smoothing method to predict the future displacement value of landslides. Crosta and Agliardi [20] discussed how to set velocity thresholds to achieve different levels for forecasting landslides. Li et al. [21] presented an application of a linear combination model with optimal weight in landslide displacement prediction. Bozzano et al. [22] described four years of continuous slope monitoring using an integrated platform, calibrating the empirical parameters of the Voight function to predict the landslide failure time. Lian et al. [23] established a set of predictors considering different environmental factors and proposed a switched method to select the appropriate individual predictor. Liao et al. [24] proposed and applied a step-like displacement prediction model based on a kernel extreme learning machine with grey wolf optimization (GWO-KELM). Krkač et al. [25] presented a methodology for the prediction of landslide movements using random forests. Meng et al. [26] used a vector autoregressive model to predict periodic displacement based on time series analysis. Landslide displacement can be divided into trend displacement, reflecting the long-term trend of landslides, which is the response of geologic structures, and periodic displacement, reflecting the volatility of landslides, which is mainly affected by external factors such as rainfall. Tang et al. [27] studied the evolution mechanism of the landslide considering the main influencing factors and providing an opinion for landslide prediction with consideration of the interior geology characteristics and external dynamic factors.

From the beginning of deformation to the final failure, landslide deformation has some inherent patterns. Terzaghi [28] first approached the correlation between the accelerating phase and landslide movements. Ter-Stepanian [29] and Tavenas and Leroueil [30] identified the presence of relevant creep deformations before failure. Furthermore, several authors discussed the possibility of defining alert threshold levels. Cruden and Masoumzadeh [31] proposed three velocity levels to correspond with three accelerating creep stages. Intrieri et al. [32] approached the task by studying the most critical periods of the entire dataset of the Torgiovannetto landslide. Moreover, according to the characteristics of the cumulative displacement and time curves of landslides, most scholars generally recognize the three-stage evolution mode of landslide deformation—the initial deformation stage [10,33,34,35,36], the constant velocity deformation stage, and the accelerated deformation stage. Some other scholars further subdivide the third stage according to the characteristics of acceleration variation. Xu et al. [33] proposed dividing the third stage into the initial acceleration, medium acceleration, and accelerated acceleration stages as shown in Figure 1.

The key to landslide prediction is to accurately estimate the dynamic variation and development trend of landslide stability. Generally, most scholars have used different types of functions for fitting, such as the Verhulst logistic function and Verhulst inverse function [37], which can be used to predict the displacement of landslides with time and divide the deformation stages. Meanwhile, random methods are the second most widely used type, such as the GM model [38,39,40,41,42,43], regression model [44], the BP network model [5,45,46], and the Markov Chain model [47,48,49]. These methods are applied to predict landslide displacement in a short time. The ARIMA model [50] and SVR model [51,52,53,54] are more suitable for short-term prediction. In addition, other scholars applied GIS technology to evaluate landslide susceptibility [55,56,57].

The prediction of landslides is mainly based on the relationship between the velocity change of displacement and time, but the research based on the acceleration and acceleration rate is limited. Considering the mechanism of landslide failure, the maximum acceleration is more consistent with the prediction criteria. At the same time, the landslide failure is mainly progressive failure, and the main feature of deformation is that the displacement acceleration is greater than zero and growing. However, scholars only deduce the formula of acceleration and acceleration maximum time for the Verhulst model, and no scholars have yet deduced and applied the formula of the acceleration and acceleration maximum time criterion for the inverse logic function curve, but they are not enough in the quantitative standard division of landslide evolution stage and the threshold of the critical sliding point.

In brief, although scholars divide the cumulative displacement time history curve of landslides into three, four, or five stages, a unified understanding is formed on its typical characteristics. Based on this unified understanding, this paper intends to adopt a new method, which is based on the inverse logistic function with a consideration of inverse distance weight, to predict landslide displacements, divide the landslide deformation stages, and realize the prediction of the critical sliding time. Lastly, four failed landslides are taken as examples to check the proposed method, and two landslides in deformation are predicted.

2. Displacement Prediction Model of Landslides

2.1. The Inverse Logistic Function Model of the Displacement-Time Curve

As shown in Figure 1, the cumulative displacement-time curve of a landslide is a typical inverse logistic function. The general expression of the logistic curve function can be written as follows:

$$y=\frac{\alpha }{1+e^{-\beta \left(x-x_0\right)}}+\gamma$$

(1)

where x₀, α (α > 0), β (β > 0), and γ are the unknown parameters to be determined. If x and y are replaced by the cumulative displacement of a landslide, S(t) and t, respectively; then the inverse function of the logistic function can be obtained based on Equation (1),

$$S\left(t\right)=\frac{1}{\beta }\mathit{ln}\left(\frac{t-\gamma }{\alpha +\gamma -t}\right)+x_0$$

(2)

where t denotes the time length of monitoring or prediction, day, and S(t) denote the cumulative displacements of the landslide, mm. Equation (2) is the general expression of the cumulative displacement-time curve.

The displacement monitoring work of a landslide could start at any time. Therefore, displacement occurring before monitoring could not be considered. However, the initial displacement value of the cumulative displacement-time function is usually set to zero at the monitoring start time, t = 0.

$$S\left(0\right)=0$$

(3)

Then, Equation (3) can be derived from Equation (2),

$$x_0=\frac{1}{\beta }\mathit{ln}\frac{\alpha +\gamma }{-\gamma }$$

(4)

Then, Equation (2) can be transformed as follows:

$$S\left(t\right)=\frac{1}{\beta }\mathit{ln}\frac{\frac{\alpha +\gamma }{-\gamma }t+\alpha +\gamma }{\alpha +\gamma -t}$$

(5)

In order to simplify the above formula, we can assume $k=\frac{1}{\beta },a=\frac{\alpha +\gamma }{-\gamma },b=\alpha +\gamma$, then Equation (5) is simplified to as follows:

$$S\left(t\right)=k\mathit{ln}\frac{at+b}{b-t}$$

(6)

where k (k > 0), a, and b are all unknown parameters to be determined. Equation (6) is the general inverse logistic function expression of the cumulative displacement-time curve of a landslide. According to Equation (6), we know that the time, t, tends to the value of b when the displacement value of the landslide tends toward positive infinity. That is, the right trend line of Equation (6) is t = b > 0.

In nature, it is a common rule that the smaller the distance is, the closer the attribute. Accordingly, when the unknown value of some attribute of one objective needs to be predicted, the measured value near the predicted time has a greater impact on the predicted value than that far from the predicted time. That is, the point nearer to the predicted position has a larger weight than the farther point, so this method is called inverse distance weighting (IDW). The simplest form of the IDW method is the reciprocal of the distance, which makes the weight sometimes very large or even infinite and then causes an error prediction result. For this reason, a constant is added to the distance when using the reciprocal mode,

$${\omega }_i=\frac{1}{d_i+const}$$

(7)

where w_i denotes the weight and d_i denotes the distance of the point to the right.

Equation (7) assigns a large weight to the nearer points, and the weight decays quickly at a slightly farther distance. Although this is what we want, it makes the algorithm more sensitive to noise. The Gauss function is complex, but it overcomes the above shortcomings. Its form is as follows [58]:

$$\omega \left(t\right)={\omega }_0{e}^{-\frac{{(t-b)}^2}{{2\mathrm{c}}^2}}$$

(8)

where w(t) is the weight value at time t, w₀ is the basic weight, and c denotes a constant to be determined. For the cumulative displacement-time curve of a landslide, b denotes the time of the right trend line. Therefore, the weight range of the function is (${\omega }_0{e}^{-\frac{{(b)}^2}{{2\mathrm{c}}^2}}$, ${\omega }_0$). Here, let ${\omega }_0=1$ and $c=\frac{b}{\sqrt{2}}$. Then the Gauss inverse distance weight function is simplified as follows:

$$\omega \left(t\right)=e^{-{(1-t/b)}^2}$$

(9)

Equation (8) is the suggested formula to calculate the weight, which is in the range of [e⁻¹, 1). Furthermore, the optimal parameters of the objective function, Equation (5), can be obtained according to the least square method of IDW, which is expressed as follows:

$$min\mathsf{\Delta }={\displaystyle \sum }_{i=1}^n{\left\{\omega \left(t_i\right)\left[S^{\prime }\left(t_i\right)-S\left(t_i\right)\right]\right\}}^2$$

(10)

where S′(t_i) and S(t_i) are the fit and monitoring displacements at time t_i, respectively.

2.2. The Division Standard of Landslide Deformation Stages

Equation (6) can be obtained from the landslide cumulative displacement-time monitoring curve, and then the function expressions of velocity, acceleration, acceleration rate with time can be obtained by the first derivative, second derivative and third derivative of Equation (5). Their formulas are expressed as follows,

$$\{\begin{array}{l}V\left(t\right)=\frac{\mathrm{kb}\left(a+1\right)}{\left(at+b\right)\left(b-t\right)}\\ a_c\left(t\right)=\frac{-\mathrm{kb}\left(a+1\right)\left(ab-b-2at\right)}{{\left(at+b\right)}^2{\left(b-t\right)}^2}\\ a_c^{’}\left(t\right)=\frac{2\mathrm{kb}\left(a+1\right)}{{\left(at+b\right)}^3{\left(b-t\right)}^3}\left[3a^2{t}^2-3\mathrm{ab}\left(\mathrm{a}-1\right)t+\left(a^2-a+1\right)b^2\right]\end{array}$$

(11)

where V(t), a_c(t) and a_c′(t) are in mm/day, mm/day² and mm/day³, respectively.

In mathematics, curvature is the rate of change in the direction of a curve with respect to the distance along the curve, and its value equals the reciprocal of the radius of the circle that most closely conforms to the curve at a given point. Therefore, the larger the curvature is, the greater the change in direction of the curve. In a general curve function, y = f(x), the curvature, K, can be calculated by the following formula:

$$K=\frac{\left|y^{″}\right|}{{\left(1+{y^{\prime }}^2\right)}^{3/2}}$$

(12)

where $y^{\prime }$ and $y^{″}$ are the first and second derivatives, respectively, of y with respect to x.

(1)

The curvature function of the cumulative displacement-time curve is

$$K_1=\frac{\left|a_c\left(t\right)\right|}{{\left(1+V^2\left(t\right)\right)}^{3/2}}$$

(13)

There are two maximum curvature extreme points in the complete ideal displacement-time curve of landslides. The corresponding time of the first maximum curvature extreme point ($K_1^{\max 1}$) is the cut-off point of the initial deformation stage and the uniform deformation stage, and the second ($K_1^{\max 2}$) corresponds to the cut-off point of the uniform deformation stage and the accelerated deformation stage. However, in reality, the monitoring work may begin at any moment, which causes the monitored cumulative displacement-time curve to miss the first maximum curvature extreme point. If so, in most cases, the start moment of the monitor is in the uniform deformation stage of the landslide; then, only the second maximum curvature is encountered. Sometimes, the curve even misses two maximum curvature extreme points, which means the landslide has reached the acceleration deformation stage.

(2)

The curvature function of the velocity-time curve is

$$K_2=\frac{\left|a_c^{\prime }\left(t\right)\right|}{{\left(1+a_c^2\left(t\right)\right)}^{3/2}}$$

(14)

Generally, there is only one maximum curvature extreme point ($K_1^{\max 2}$) for the velocity-time curve of Equation (13), which corresponds to the cut-off point of the accelerated deformation stage and the critical sliding stage.

A whole cumulative displacement-time curve of landslides could be divided into four different stages: the initial deformation stage (the time before the time of $K_1^{\max 1}$), the uniform deformation stage (the time between $K_1^{\max 1}$ and $K_1^{\max 2}$), the acceleration deformation stage (the time between $K_1^{\max 2}$ and $K_2^{\max }$), and the critical sliding stage (the time after $K_2^{\max }$), as shown in Figure 2.

The above proposed method needs to be verified to determine whether it is available and accurate when used to predict the displacement of landslides. Here, four famous failed landslides in China, including the Xintan landslide, Wolongsi New landslide, Huangci landslide, and Saleshan landslide, are taken as typical examples for prediction by our proposed method. Finally, two ongoing landslides are predicted.

3. Verification with Typical Landslides

In order to certify the accuracy of the proposed method in this paper, four cases, Xintan Landslide, Wolongsi landslide, Huangci landslide and Saleshan’s new landslide, were applied to divide the deformation stages and predict the sliding time. Furthermore, the prediction results were compared with other methods. Six other landslide cases were applied, and two ongoing landslides were predicted. Figure 3 shows the specific location of all the landslide cases.

3.1. Case 1: Xintan Landslide

The Xintan landslide is an accumulation landslide that failed on 12 June 1985. The area of the landslide is approximately 0.68 × 10⁶ m², and its volume reaches 30.0 × 10⁶ m³, most of which slid into the Yangtze River. The landslide is approximately 2000 m in length from north to south and 450 m in width. The longitudinal average gradient of the landslide is 23 degrees. The height difference of the landslide is approximately 800 m. Permian limestone and Devonian quartz sandstone are exposed in the area of Guangjiaya at the back of Xintan landslide, while the landslide is located above the Silurian sandstone and sand shale, which are easily weathered, broken, and impermeable. Xintan landslide is the revival of ancient landslide deposits. The main reason is that there are often-produced collapse deposits from the Guangjiaya, which push the ancient landslide deposits moving forward. Furthermore, the surface land reclamation is digging and planting in disorder, and the surface drainage is not smooth.

Since 1978, the upper area of the Xintan landslide has experienced obvious deformation, and two displacement monitoring points (No. A3 and B3) were set up on the landslide in January 1978. The cumulative displacement data of point No. B3 on the Xintan landslide is shown in

Table 1. The measured displacements of the Xintan landslide at point B3 (mm) [59].

	Month
Year	Jan	Feb	Mar	Apr	May	June	July	Aug	Sept	Oct	Nov	Dec
1978	13.3	16.7	20.5	24.8	28.8	43.4	60.5	66.5	75.0	77.6	77.8	87.8
1979	89.4	93.3	98.2	99.7	104.4	115.5	126.4	149.9	302.0	367.0	404.5	412.7
1980	427.9	437.8	442.4	451.4	456.4	471.4	520.4	588.4	688.4	690.1	721.8	723.9
1981	731.4	739.1	747.5	762.7	762.7	784.0	788.2	791.1	817.1	819.1	824.1	832.1
1982	837.1	840.6	848.1	926.1	926.1	1096.3	1117.8	1254.0	1504.8	1747.8	1823.6	1895.1
1983	1995.8	2024.5	2056.4	2087.6	2087.6	2148.5	2277.9	2430.0	2582.1	2734.2	2886.3	3038.1
1984	3190.5	3342.6	3494.7	3646.8	3646.8	3951.0	4188.7	4521.7	4673.5	5195.6	5467.3	5705.0
1985	5863.3	6012.7	6177.2	6304.5	6304.5	-	-	-	-	-	-	-

, which is used to predict through our proposed method, as shown in Figure 4 and Figure 5.

It is found that the annual cumulative displacement from 1978 to 1981 is hundreds of millimeters. However, from 1982 to 1985, cumulative displacement increased sharply to thousands of millimeters. The Xintan landslide continued to deform several years at an annual displacement of 1000 mm. Sun [60] considered that the Xintan landslide was in the initial creep stage before 1979, in the constant creep stage from August 1979 to July 1982, in the accelerated creep stage from July 1982 to 15 May 1985, entered the stage of catastrophic failure between 15 May 1985 and 12 June 1985, and then arrived in the critical sliding stage. The Rockfall and Landslide Research Institute of Hubei Province [61] thought that the Xintan landslide was in the potential deformation and creep stage before 1982, the revival stage from March 1982 to May 1983, the rapid deformation stage from June 1983 to June 1985, and the sliding stage from 9 to 12 June 1985. Furthermore,

Table 2. The prediction results of the Xintan landslide based on different monitoring data.

Parameters		Deformation Stages
		The Start Point of the Acceleration Deformation Stage	The Start Point of the Critical Sliding Stage
~1 April 1985	time	12 December 1980	18 June 1985
k = 2935	V(t)	1.07	53.3
a = −0.51	a(t)	8.75 × 10−4	1.01
b = 2810	a′(t)	1.24 × 10−6	0.110

lists the deformation stage division and prediction results of the Xintan landslide by the proposed method using the monitoring data during different durations of time.

According to the inverse logic curve, the prediction time of the landslide is 18 June 1985, while the actual landslide is 12 June 1985. The prediction results of the landslide are compared with those of other different methods as shown in

Table 3. Prediction results of Xintan landslide.

Model	Prediction Physical Parameter	Prediction Time	Δ/d *
Inverse Logic Function Model	Displacement acceleration	18 June 1985	6
Original synergetic model [59]	Displacement acceleration	16 October 1987	856
Improved synergetic model [59]	Displacement acceleration	14 March 1985	−90
Displacement GM (1,1) Model [62]	Displacement velocity	13 December 1981	−1277
Velocity GM (1,1) Model [62]	Displacement velocity	1 July 1982	−1077
Original Verhulst model [63]	Displacement velocity	24 July 1987	772
Improved Verhulst model [63]	Displacement acceleration	3 February 1986	236
Improved Verhulst model [63]	Displacement accelerated acceleration	7 February 1985	125

* where Δ is the absolute error, which refers to the difference between the predicted landslide occurrence time t_y and the actual landslide occurrence time t_s, namely, Δ is t_y minus t_s, and Δ unit is day. When the actual occurrence time of landslide is earlier than the prediction time, Δ > 0. When the actual occurrence time of landslide is later than the prediction time, Δ < 0.

. It shows that the prediction accuracy of our method is much better than any other one.

3.2. Case 2: New Wolongsi Landslide

The New Wolongsi landslide is a loess landslide. The thickness of the landslide mass is more than 50 m to as much as 90 m, it is 645 m long and 650 m wide, and the volume scale is more than 20.0 × 10⁶ m³. It is a super-deep and super-large landslide. Damaged several times in the past, the failure range of each time is basically the same. It is the fourth terrace of the Weihe River at the back edge of the landslide. The upper part of the terrace is the loess of 97 m thickness, the lower part is river alluvium of 54 m thickness, and the lower part is Neogene mild clay and gravel. The gravel usually contains water. The rear sliding surface is arc-shaped with a steep slope, the middle front part is nearly horizontal, and the front edge is a reverse slope that overlaps the floodplain terrace. Since the crack was first discovered in March 1971, the No. 5 fissure on the New Wolongsi landslide was monitored with manual measurement for 66 days, and slip damage occurred on 5 May.

Taking the cumulative displacement-time monitoring data, which are shown in

Table 4. The measured displacements of the New Wolongsi landslide [59].

Month	March
Date	15	16	17	18	19	20	21	22	23	24	25	26	27	28	29	30	31
displacement (mm)	1.0	1.5	1.7	2.5	3.2	4.0	4.4	5.1	5.9	6.3	7.0	7.3	7.8	8.2	8.4	8.7	9.0
Month	April
Date	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17
displacement (mm)	9.2	9.4	10.0	10.1	10.3	10.4	10.5	10.8	11.1	12.0	13.0	13.6	14.0	15.0	16.1	16.4	17.2
Month	April													May
Date	18	19	20	21	22	23	24	25	26	27	28	29	30	1	2	3	4	5
displacement (mm)	17.0	18.0	19.0	19.0	20.0	23.0	24.0	25.2	26.0	27.0	28.2	30.0	31.0	32.0	33.0	42.0	47.0	61.0

, of the No. 5 fracture as an example, the proposed method is used for prediction, and the results are shown in Figure 6 and Figure 7 and

Table 5. The prediction results of the New Wolongsi landslide based on different monitoring data.

Parameters		Deformation Stages
		The Start Point of the Acceleration Deformation Stage	The Start Point of the Critical Sliding Stage
~3 May 1971	time	27 April 1971	3 May 1971
k = 16.00	V(t)	1.42	4.03
a = 0.10	a(t)	0.109	1.00
b = 54.0	a′(t)	0.017	0.507

.

The prediction time of the landslide by the inverse logic curve is on 3 May 1971, which is basically consistent with the actual landslide time on 5 May 1971. The comparison with other methods is as shown in

Table 6. Prediction results of Wolongsi’s new landslide.

Model	Prediction Physical Parameter	Prediction Time	Δ/d
Inverse Logic Function Model	Displacement acceleration	3 May 1971	−2
Original synergetic model [59]	Displacement acceleration	26 May 1971	21
Improved synergetic model [59]	Displacement acceleration	5 May 1971	0
Improved Verhulst model [63]	Displacement acceleration	9 May 1971	4
Improved Verhulst model [63]	Displacement accelerated acceleration	28 April 1971	−7
Original Verhulst model [63]	Displacement velocity	24 April 1971	19

. It shows that most of the methods, including our proposed method, have good prediction accuracy.

3.3. Case 3: Huangci Landslide

The Huangci landslide is not entirely a new landslide; it has been damaged several times in its history. The landslide is close to the edge of the fourth terrace, and the settlement cracks, which evolved to the back-edge cracks of Huangci landslide caused by irrigation, have appeared for more than ten years. The sliding surface of the west bedding plane has also existed for several years. Huangci landslide is nearly 500 m wide at the rear and 300 m wide on the front; it is 370 m long in the north and south directions, and the volume of landslide is nearly 6.0 × 10⁴ m³. The thickness of the loess layer at Huangci landslide is 42.67 m, the thickness of light-yellow loess (Q₃) at the upper part is 36.40 m, and the thickness of brown-red/brown-yellow loess (Q₂) at the lower part is 6.27 m. The soft plastic layer with a thickness of 9 m is 31.73–40.90 m below the ground, and the plastic layer is 40.90–42.67 m below the ground. The lower part of the loess layer is an 8 m-thick gravel layer. Below the gravel is Cretaceous mudstone with argillaceous sandstone. The bedrock belongs to monoclinic structure; it inclines to Huangci village in front of the landslide.

A large amount of irrigation water is the main inducing factor of the Huangci landslide. The initial monitoring time of the Huangci landslide was on 1 August 1994, and the monitoring time interval was 15 days. The ground displacement monitoring data of No. A6 by GPS are shown in

Table 7. The measured displacements of the Huangci landslide at point A6 [59].

Date	Year/Month	1994/8			1994/9		1994/10		1994/11		1994/12		1995/1
	Date	1	16	31	15	30	15	30	14	29	14	29	13	28
displacement (mm)		1.0	3.3	6.1	8.3	10.5	15.3	18.1	21.6	27.5	33.6	39.0	48.7	64.3

. At 2:30 a.m. on 30 January, 1995, the Huangci landslide slipped for 90 min before it completely stopped. The leading edge of the landslide went straight to Huangci Village. Our proposed method is also used on the Huangci landslide, and the results are shown in Figure 8 and Figure 9 and

Table 8. The prediction results of the Huangci landslide based on different monitoring data.

Parameters		Deformation Stages
		The Start Point of the Acceleration Deformation Stage	The Start Point of the Critical Sliding Stage
~28 January 1995	time	29 January 1995	11 February 1995
k = 25.9	V(t)	1.44	5.19
a = 0.05	a(t)	0.080	1.04
b = 216	a′(t)	8.92 × 10−3	0.416

.

According to the inverse logic curve, the prediction time of the landslide is February 11, 1995, and the actual landslide time is January 30, 1995. The comparison with other methods is shown in

Table 9. Prediction results of Huangci landslide.

Model	Prediction Physical Parameter	Prediction Time	Δ/d
Inverse Logic Function Model	Displacement acceleration	11 February 1995	−13
Original synergetic model [59]	Displacement acceleration	26 March 1995	55
Improved synergetic model [59]	Displacement acceleration	29 January 1995	−1
Verhulst Model [62]	Displacement acceleration	25 December 1994	−66
Displacement GM (1,1) Model [62]	Displacement velocity	17 December 1994	−74
Velocity GM (1,1) Model [62]	Displacement velocity	25 December 1994	−66

. It shows that the improved synergetic model has the best prediction accuracy, followed by our proposed method, and the other four methods have a larger prediction error.

3.4. Case 4: Saleshan Landslide

The Saleshan landslide occurred in Dongxiang County, Gansu Province. Its failure first happened on 7 March 1983. The Saleshan landslide is mainly composed of red-layered clay rock of the Linxia Formation. The red layer of the Linxia Formation is covered with loess, in which some vertical discontinuities were developed. The sliding bed of the Saleshan landslide is mainly composed of mudstone.

On 25 March 1986, the Saleshan landslide failed again. Its volume was approximately 2.40 × 10⁶ m³, and the sliding distance was approximately 250 m. Displacement monitoring was carried out after it occurred the first time, and the displacement monitoring data are shown in

Table 10. The measured displacements of the Saleshan landslide [59].

Year	1984							1985												1986
Month	6	7	8	9	10	11	12	1	2	3	4	5	6	7	8	9	10	11	12	1	2
Displacement (mm)	4	10	0	1.5	2.5	5.5	7.2	5.2	5.2	4.8	2.5	4.8	7	12	15	17	27	30	41	50	75

. The prediction results using our proposed method are shown in Figure 10 and Figure 11 and

Table 11. Prediction results of the Saleshan landslide based on different monitoring data.

Parameters		Deformation Stages
		The Start Point of the Acceleration Deformation Stage	The Start Point of the Critical Sliding Stage
~15 February 1986	time	10 February 1986	23 February 1986
k = 25.4	V(t)	1.36	5.02
a = −0.60	a(t)	0.078	1.02
b = 654	a′(t)	0.042	2.15

.

The prediction time of the landslide by the inverse logic curve is 23 February 1986, while the actual time of the landslide is 25 March 1986. The comparison with other methods is shown in

Table 12. Prediction results of Saleshan’s new landslide.

Model	Prediction Physical Parameter	Prediction Time	Δ/d
Inverse Logic Function Model	Displacement acceleration	23 February 1986	−30
Original synergetic model [58]	Displacement acceleration	22 May 1986	58
Improved synergetic model [58]	Displacement acceleration	5 March 1986	−20
Original Verhulst model [62]	Displacement velocity	4 July 1986	101
Improved Verhulst model [62]	Displacement acceleration	15 April 1986	21
Improved Verhulst model [62]	Displacement accelerated acceleration	16 February 1986	−37

. It shows that the improved synergetic and Verhulst models have the best prediction accuracies, followed by our proposed method, and the other three methods have larger prediction errors.

The factors affecting the prediction accuracy are as follows: (1) The more monitoring data that is accumulated, the higher the prediction accuracy will be. (2) The monitoring points could be representative. The deformation of different parts of the landslide is not synchronous. Therefore, the choice of the monitoring points applied to predict should be able to reflect the overall deformation trend of the landslide. For example, for the pull-type landslide, it is more appropriate to set monitoring points at the front edge of the landslide. For the thrust-type landslide, it is more appropriate to set monitoring points at the back edge of landslide.

4. Landslide Predictions

The long-term displacement predictions of the Liangshuijing landslide and Gapa landslide were carried out by the proposed method, and the deformation evolution stage of each landslide was studied and divided; then, the critical sliding time and the change of each parameter were predicted.

4.1. Case 1: Long Displacement Prediction of Liangshuijing Landslide

Liangshuijing landslide is located on the right bank of the Yangtze River with a slope angle of 30 to 35 degree. The height of the leading edge of the landslide is approximately 100 m, the height of the trailing edge is approximately 319.5 m, the relative height difference is approximately 219.5 m, the plane longitudinal length is approximately 434 m, the lateral width is 358 m, the area is approximately 11.82 × 10⁴ m², the average thickness of the slide body is approximately 34.5 m, and the total volume is approximately 4.08 × 10⁶ m³. The professional monitoring was just carried out after the landslide experienced a large deformation in March 2009.

Manual horizontal displacement monitoring began on 5 April 2009 and stopped on 22 April 2009. Automatic horizontal displacement monitoring began on 20 April 2009, and the surface displacement monitoring point ZJC22 on the right side of the central part of the landslide was selected for prediction, as shown in

Table 13. The measured displacements of Liangshuijing landslide.

Time yy/mm/dd	Displacement /mm	Time yy/mm/dd	Displacement /mm	Time yy/mm/dd	Displacement /mm
2009/04/05	0.00	2009/06/14	215.90	2009/08/23	269.87
2009/04/10	18.83	2009/06/19	229.71	2009/08/28	266.11
2009/04/15	59.00	2009/06/24	242.26	2009/09/02	273.64
2009/04/20	96.65	2009/06/29	241.00	2009/09/07	273.64
2009/04/25	105.44	2009/07/04	243.52	2009/09/12	266.11
2009/04/30	112.97	2009/07/09	243.52	2009/09/17	269.87
2009/05/05	115.48	2009/07/14	252.30	2009/09/22	273.64
2009/05/10	133.05	2009/07/19	251.05	2009/09/27	272.39
2009/05/15	134.31	2009/07/24	259.83	2009/10/02	278.66
2009/05/20	156.90	2009/07/29	259.83	2009/10/07	273.64
2009/05/25	166.95	2009/08/03	254.81	2009/10/12	278.66
2009/05/30	184.52	2009/08/08	258.58	2009/10/17	277.41
2009/06/04	197.07	2009/08/13	264.85	2009/10/22	281.17
2009/06/09	207.11	2009/08/18	264.85	2009/10/27	287.45

. The monitoring data were maintained until October 2009.

The cumulative variation curve of horizontal displacement at the monitoring points is shown in Figure 12, and the monitoring results are shown in

Table 14. The prediction results of the Liangshuijing landslide based on different monitoring data.

Parameters		Deformation Stages
		The Start Point of the Uniform Deformation Stage	The Start Point of the Acceleration Deformation Stage	The Start Point of the Critical Sliding Stage
~15 November 2009	time	30 May 2009	3 September 2086	26 October 2086
k = 89.51	V(t)	1.42	1.43	9.14
a = 3520.0	a(t)	−0.023	0.023	0.932
b = 28,337.8	a′(t)	7.14 × 10−4	2.17 × 10−3	0.571

.

4.2. Case 2: Long-Term Displacement Prediction of the Gapa Landslide

The Gapa landslide is located on the right bank of the Yalong River. This area is approximately 11 km upstream of the Jinping first-class hydropower station, and the volume of the Gapa landslide is approximately 13.0 × 10⁶ m³, The landslide has a length of 980 m, a width of 320 to 400 m, an area of 0.28 km².

There are seven surface displacement monitoring points, denoted as GNSS1–GNSS6, and GNSS9 near the landslide area. GNSS1 and GNSS2 are located in the number one sliding body, which is the most likely failure area (see Figure 13). Therefore, GNSS1 and GNSS2 displacement monitoring data are selected for prediction, and the monitoring results are shown in

Table 15. The prediction results of the Gapa landslide based on different monitoring data.

Parameters		Deformation Stages
		The Start Point of The Uniform Deformation Stage	The Start Point of the Acceleration Deformation Stage	The Start Point of the Critical Sliding Stage
GNSS1	time	10 February 2017	25 September 2034	26 July 2037
k = 1877	V(t)	1.98	1.98	40.20
a = 11.41	a(t)	−1.57 × 10−3	1.57×10−3	0.85
b = 7907	a′(t)	2.99 × 10−6	8.56×10−6	0.10
GNSS2	time	6 June 2016	19 May 2029	16 October 2031
k = 1706	V(t)	2.15	2.15	42.9
a = 7.40	a(t)	−1.95 × 10−3	1.96 × 10−3	1.07
b = 5793	a′(t)	4.37 × 10−6	1.20 × 10−5	0.14

.

The prediction results show that the deformation is not synchronized in different parts of the Gapa landslide. In the number one sliding body, the former part of the Gapa landslide is destroyed before the trailing part, which illustrates that the Gapa landslide is a retrogressive type of landslide. Furthermore, the Gapa landslide will not fail in the next 10 years.

5. Warning Threshold of Critical Sliding

In order to find out one or some parameters as an effective tool to predict the critical sliding point of landslides, the accumulated displacements, velocities and accelerations of above four landslides and some other landslide cases from references at the moment of the critical sliding are listed in

Table 16. Velocity and acceleration characters at the critical sliding moment of different landslides.

Landslides	Parameters
	Degree of Slope Surface (°)	Degree of Sliding Surface (°)	Depth (m)	Different Elevator (m)	Width (m)	Length (m)	Volume (106 × m3)	V (mm/day)	a (mm/day2)	a′ (mm/day3)
Xintan landslide	25–28	25–28	40–50	800	270–700	2000	30	53.3	1.01	0.110
Wolongsi landslide	0–25	0–45	50–90	200	645	650	20	4.65	1.00	0.507
Huangci landslide	15–35	15–20	40	100	300–500	370	6	5.19	1.04	0.416
Saleshan landslide	15–45	25–35	30–40	425	275	250	2.4	5.02	1.02	2.150
Liangshuijing landslide	30–35	25–50	34.5	221.5	358	434	4.08	9.14	0.93	0.571
Baishi Landslide [64]	40	38	25	270	260	300	2.0	4.00	1.07	0.569
Tianhuangping [65]	35–50	40.3	40	70	60	100	0.062	6.45	1.02	0.329
Xikou Landslide [47]	65	37	24	108	80–90	80	0.168	8.89	0.99	0.220
Jianshanbeibang landslide [66]	/	/	5–15	240	240	105	0.252	15.26	1.02	0.710
Jimingshi Landslide [67]	45	30	15	230	150	250–300	0.6	4.04	1.00	1.49
Shuiwenzhan Landslide	28–33 20–21	28–33	64–120 (80)	500	570	1000	15	28.75	0.85	0.138
Gapa landslide	30–34	23–34	60	470	360	980	13	42.90	1.07	0.140

.

According to Table 16, Figure 14 is obtained. Obviously, it is found that the values of acceleration at the critical sliding moment are very close to all landslides; their values are at a range of between 0.85 and 1.07 mm/day². Therefore, the acceleration can be thought as a key index to predict the critical sliding point, and the threshold of displacement acceleration judging that the critical sliding point is 0.90 mm/day².

Furthermore, the relationship of the acceleration rate with the geometric parameters of landslides has skewed distribution characteristics, as shown in Figure 14, while Gamma function is one of the commonly used skewed distribution functions. Gamma function is shown as the following formula:

$$\{\begin{array}{l}f\left(x,\alpha ,\beta \right)=\frac{1}{{\beta }^{\alpha }\Gamma \left(\alpha \right)}{x}^{\alpha -1}{e}^{-\frac{x}{\beta }}, x>0\\ \Gamma \left(\alpha \right)={{\displaystyle \int }}_0^{+\infty }{e}^{-t}{t}^{\alpha -1}dt\end{array}$$

(15)

where α and β are shape parameter and scale parameter. Then, we can obtain the shape and scale parameters of Gamma function, which fits the relationship of the acceleration rate (mm/day³), a_c’, with the volume of landslide (× 10⁶ m³), Vol, and the length of landslide (m), L. Based on the data in Table 16, the critical threshold function can be fitted as shown in Equation (15).

$$\{\begin{array}{l}{a}_c^{’}\left(L\right)=900f\left(L,5,8\right)+\frac{1-e^{-L}}{10\left(1+e^{-L}\right)}\\ a_c^{’}\left(Vol\right)=8f\left(Vol,3,1\right)+\frac{1-e^{-Vol}}{10\left(1+e^{-Vol}\right)}\end{array}$$

(16)

In the above formula, the function, f, is Gamma function as Equation (15), and the second part on the right side is an additive term, which is for improving fitting precision.

Sometimes, the geometric parameters are evaluated; they are not so accurate, and will reduce the prediction, effect shown in Equation (16). The displacement acceleration is almost a constant when a landslide enters a sliding state. Therefore, the value of displacement acceleration, 0.9 mm/day², is suggested as the first early warning standard of sliding, and Equation (16) can be used as the secondary critical threshold function of landslide failure.

6. Conclusions

The prediction model and early warning criteria are the key to landslide prediction, most of which are judged by the characteristics of the cumulative displacement. However, due to the different sizes of landslides and the complexity, nonhomogeneity, and uncertainty of their attributes, the cumulative displacement before failure is different by a few centimeters to dozens of meters, and the deformation velocity of landslides at the moment of critical sliding can also vary from a few millimeters per day to dozens of centimeters per day. Therefore, the general threshold values of displacement and velocity are difficult to predict for different landslides. Moreover, there is no quantitative standard for dividing the landslide deformation stages.

Based on the analysis and study of the variations in the cumulative displacement, deformation rate, displacement-time curvature, and velocity-time curvature of typical landslides, a suitable method is built to predict the displacement of landslides that can calculate the critical sliding time of landslides. This proposed method can predict the landslide displacement, divide the landslide deformation stages, and gain the threshold values of the critical sliding point by displacement acceleration and acceleration rate. Some meaningful conclusions and understandings are obtained:

(1)

The cumulative displacement-time curves of landslides are in good agreement with the inverse logistic function. A general expression of the cumulative displacement-time function of landslides based on the logistic inverse function is proposed. Furthermore, the least squares formula of the inverse logistic function prediction with the IDW method is suggested to fit the inverse logistic function based on the displacement monitoring data of landslides.

(2)

Based on the prediction model of the inverse logistic function with the IDW method, a new standard to divide the landslide deformation stages is proposed. This method was applied in some typical landslides and verified as effective and accurate. Lastly, this method was also used to predict deformation in two landslides—the Gapa and Liangshuijing landslides—and the critical sliding times were calculated. It should be noted that these two landslides remain in the initial deformation stage or just enter into the uniform deformation stage according to the monitoring displacement data. The prediction results are suggested to be continuously renewed when more monitoring displacement data are obtained.

(3)

The prediction accuracy is mainly decided by the cumulative monitoring time and the distance far from the critical sliding time, and the critical sliding prediction time will be more accurate when the monitoring data are closer to the failure times of the landslides; otherwise, the prediction time may be farther from the real failure time. Therefore, we should predict the critical sliding time continually with the increase in monitoring data, which will improve the prediction accuracy.

(4)

The displacement acceleration is recommended as a key index to predict the critical sliding of landslide, and its threshold value is suggested as 0.90 mm/day². Furthermore, the supplementary index of the critical sliding moment is the displacement acceleration rate, and the critical threshold function is suggested as Equation (15).