Application of GWO-ELM Model to Prediction of Caojiatuo Landslide Displacement in the Three Gorge Reservoir Area

: In order to establish an e ﬀ ective early warning system for landslide disasters, accurate landslide displacement prediction is the core. In this paper, a typical step-wise-characterized landslide (Caojiatuo landslide) in the Three Gorges Reservoir (TGR) area is selected, and a displacement prediction model of Extreme Learning Machine with Gray Wolf Optimization (GWO-ELM model) is proposed. By analyzing the monitoring data of landslide displacement, the time series of landslide displacement is decomposed into trend displacement and periodic displacement by using the moving average method. First, the trend displacement is ﬁtted by the cubic polynomial with a robust weighted least square method. Then, combining with the internal evolution rule and the external inﬂuencing factors, it is concluded that the main external trigger factors of the periodic displacement are the changes of precipitation and water level in the reservoir area. Gray relational degree (GRG) analysis method is used to screen out the main inﬂuencing factors of landslide periodic displacement. With these factors as input items, the GWO-ELM model is used to predict the periodic displacement of the landslide. The outcomes are compared with the nonoptimized ELM model. The results show that, combined with the advantages of the GWO algorithm, such as few adjusting parameters and strong global search ability, the GWO-ELM model can e ﬀ ectively learn the change characteristics of data and has a better and relatively stable prediction accuracy.


Introduction
Landslide is one of the worst geological disasters in many areas of the world, with the characteristics of wide distribution, strong paroxysm, high frequency and great harmfulness. Landslide disaster can not only cause a large number of casualties and property losses, but also bring huge damage to resources, environment and ecology [1]. The Three Gorges Reservoir area in China is prone to In this paper, a typical step-wise-characterized landslide called Caojiatuo landslide in the TGR area is taken as an example. Based on the time series analysis principle, the displacement is decomposed into two components: the trend displacement and periodic displacement [27]. The trend displacement is extracted by using the moving average method, and predicted by using the cubic polynomial function of the robust weighted least-square method. In this paper, the GWO is innovatively introduced into the ELM, and a landslide displacement prediction model based on GWO-ELM is proposed to predict the periodic displacement. Root-mean-square error (RMSE) and mean absolute percent error (MAPE) are applied as evaluation criteria to analyze the predicted results. This new model helps to build a more accurate landslide warning system.

Time Series Model
Landslide displacement is a nonstationary time series varying with time [5]. The study shows that the accumulative displacement of the landslide is caused by the interaction of the internal geological conditions (lithology, geological structure, landform, etc.), external environmental factors (precipitation, reservoir water level, etc.) and other random factors (uncertainty) [28]. The displacement of the landslide under the control of its own internal basic geological conditions is expressed as an approximately monotone increasing function varying with time, that is, the trend displacement. The displacement under the influence of external environmental factors is shown as an approximately periodic function which changes with seasonal rainfall and reservoir water level, namely the periodic displacement. Based on the above features, the time series decomposition theory is adopted to decompose the accumulative displacement of the landslide: where X(t) represents the monitored displacement; α(t) is the trend displacement; β(t) is the periodic displacement; γ(t) is the random displacement. The random displacement is mainly caused by the random influencing factors (wind load, vehicle load and vibration load, etc.). Due to the limitation of the current monitoring technology, it is difficult to obtain the data of such factors, and the influence of them is relatively small. Generally, it is not considered in the research [29]. Therefore, Equation (1) can be simplified as follows: Generally speaking, the trend displacement can be extracted by the moving average method, which can smooth the short-term fluctuation of time series and reflect the long-term trend. The specific calculation method is to assume that the monitoring value of landslide displacement at time t is s(t), and the trend displacement at time t is β(t). The calculation formula is as follows: where t = k, k + 1, . . . , n; s(t − n + 1) is the monitoring value of landslide displacement at time t − n + 1; n is the moving cycle. Based on the time series decomposition theory, after the extraction of the trend displacement is completed, the periodic displacement α(t) is obtained by removing the trend displacement β(t) from the accumulative displacement X(t).

Extreme Learning Machine Model
ELM is a training method based on the Single-hidden Layer Feedforward Neural Network (SLFNs) [29]. The output weight of the network can be obtained by one-step calculation and analysis.
Water 2020, 12, 1860 4 of 22 In this process, the only optimal solution can be obtained and the learning speed is very high just by setting the threshold value of the neurons in the hidden layer. For N arbitrary and different samples ( , y i2 , . . . , y im ] T ∈ R m , and there are L neuron nodes in the hidden layer. Then the output of the standard feedforward neural network, the excitation function of which is g(x) can be expressed as: where a i = [a i1 , a i2 , . . . , a in ] T is the input weight from the input layer to the nodal point of the ith hidden layer; b i is the threshold value of the ith neuron in the hidden layer; β i = [β i1 , β i2 , . . . , β im ] T is the output weight connecting the ith node in the hidden layer; a i × x i is the inner product of vector a i and x i ; the excitation function g(x) can be selected as "Sigmoid", "RBF", or "Sine", etc.
If the standard feedforward neural network with L hidden layer neurons and excitation function g(x) can approach these N samples with zero error, then there exist a i , b i , c i making the following equation true: Equation (5) can be simplified to: where H is the output matrix of the neural network. In ELM, the output weight and threshold are randomly given, and the hidden layer matrix H namely becomes a definite matrix, which makes the training of feedforward neural network transform to a problem of finding the least square solution of the output weight matrix; the output weight matrix β is: where H + represents the Moore-Penrose generalized inverse matrix of the hidden layer output matrix H, which can be obtained analytically by orthogonal projection or singular decomposition method, etc.

Gray Wolf Optimization
Gray Wolf Optimization (GWO) algorithm is a new swarm intelligent optimization algorithm proposed by MIRJALILI [25]. GWO simulates the process of the gray wolf hunting behavior such as wolf tracking, encircling, chasing, and attacking prey, so as to achieve the goal of optimization. GWO has the advantages of having a simple principle, few adjustment parameters, strong global search ability and so forth.
The optimization process of GWO is as follows: a group of gray wolves is randomly generated in the search space. According to the fitness level, the wolves numbered α, β, δ in the gray wolf group evaluate and locate the position of the prey. The rest of the individuals take this as the standard, and calculate the distance between themselves and the prey to complete the prey capture and realize the optimization process [25]. Definition 1. Social Hierarchy. The GWO model has a strict social class, which can be divided into α, β, δ, ω wolves on the basis of the social class. The social class in the algorithm is reflected in the fitness level.

Definition 2.
Encircling Prey. In the process of wolf hunting, they need to s encircle the prey and determine the location of the prey. The mathematical equations of encirclement behavior are as follows: Water 2020, 12, 1860 5 of 22 where w is the current iteration number; X p (w) is the position vector of the prey in generation w; X(w) is the position vector of individual gray wolf in generation w; X(w + 1) is the position vector of individual gray wolf in generation w + 1; C is coefficient vector; r 1 , r 2 belong to random vector, with the value range of [0, 1]; µ is a vector of convergence; the value of a declines linearly from 2 to 0 during the iteration.
Definition 3. Hunting. The hunting process of wolves on prey is represented by the constant replacement of hunting location information. Specifically, in the iterative process, the algorithm saves the current locations of the best three wolves (α, β, δ), and updates the location of wolves in other search units (ω) according to their information to obtain the optimal solution. Hunting behavior can be expressed as: where k = α, β, γ; i = 1, 2, 3; X p (w + 1) is the position vector of the prey in w + 1 generation.

Proposed Prediction Model Based on GWO-ELM
In this paper, based on the principle of time series, the accumulative displacement is decomposed into the trend and periodic displacement. For trend displacement, the robust weighted least squares cubic polynomial fitting projections is adopted. For periodic displacement, the GWO is introduced to the ELM model for training. The framework system of the proposed extreme learning machine model, namely GWO-ELM model, is shown in Figure 1.

Geologic Aspects
The Caojiatuo landslide [30] is in the TGR area of China. After the construction of the Three Gorges Dam, the hydrogeological conditions in this region has been dramatically changed and the TGR area has become a world-famous landslide-prone area. The Caojiatuo landslide, as one of the typical landslides, locates on the left bank of the Yangtze river, Wushan county, Hubei province ( Figure 2). The plane shape of the landslide is close to a square approximately 600 m long and 500 m wide. The landslide occupies a mass with an area approximately 30 × 10 4 m 2 and a volume approximately 1500 × 10 4 m 3 ( Figure 3). The sliding surface was the interface between the loose accreted body and the underlying bedrock, with the length of approximately 600 m, the width of approximately 500 m, dip angle of 30 • , and dip direction of 150 • (Figure 4). The thickness of the contact zone was approximately 2 m, mainly comprising clay with a small amount of debris.
The landslide body is composed of colluvial and residual silty clay, which was yellow-brown and easily softened by water [2]. The material also contains approximately 35% rubble, which reduces from the surface to the bottom of the landslide mass, and consists of gray-green silty siltstone of size 5-10 cm. The bedrock mainly consisted of silty siltstone, gray-green marlite, and a small amount of magenta mudstone of the Badong Formation in the Triassic.

Monitoring Data Analysis
The water level of the TGR reached 135 m in 2003, and the Caojiatuo landslide was first deformed and tensile cracks appeared along the trailing edge. With the development of the landslide, the cracks increased year by year. In order to observe the surface deformation of the landslide and avoid catastrophic geological disasters, since 2007, the landslide was divided into two regions: zone A and zone B and nine GPS sensors ( Figure 3) were installed to monitor the stability of the landslide per month, especially during the rainy season. Figure 5 shows the relationship between the monthly displacement of the 9 GPS observation points of Caojiatuo landslide and the reservoir water level and precipitation during 2007-2013. It can be seen that the displacement of block A (GPS-2, GPS-3, GPS-5, GPS-6) is a typical step-wise evolution process. The displacement of the landslide surface increased rapidly between May and July every year, and gradually approached a stable state in the following period. The displacement of block B (GPS-1, GPS-4, GPS-7, GPS-8, GPS-9) shows a different variation feature. The deformation of the monitoring points in block B was relatively small, and the fluctuation with the changes of the reservoir water level and precipitation was not obvious. Therefore, the deformation of the landslide mass of block B is defined as the trail-type deformation.  The accumulative displacement in block A was obviously affected by seasonal precipitation and reservoir water level, which is conducive to the study of factors affecting landslide deformation. The displacement values monitored by the GPS-6 and GPS-3 observation points in block A are relatively large. The monitoring data are representative and can effectively reflect the movement process of the whole landslide. Therefore, the data of GPS-6 and GPS-3 are selected for analysis in this paper. According to the curves of the monitored accumulative displacement, monthly precipitation and reservoir water level in Figure 5, the data can be divided into three stages: (1) Stage 1 (7 January-9 May): The displacement growth was slow in this stage. According to the curve of reservoir water level, the reservoir began to fill water every year in September. Around January 2009, the maximum water level of the reservoir reached 175 m for the first time. With the increase of the reservoir water level, the displacement of the monitoring points in the current month did not change significantly, indicating that the influence of the rising reservoir water level on the displacement was weak. The water level of the reservoir reached the lowest value around May each year, and the displacement of the corresponding monitoring point increased to a certain extent after a month. The reason for the phenomenon was that the change of the water level inside the landslide lagged behind the water level. According to the precipitation curve, when seasonal heavy rainfall occurred in the landslide area (mainly from May to September every year), the landslide displacement continued to increase, which reflected the correlation between landslide displacement and rainfall. Therefore, in Stage 1, the displacement of monitoring points was mainly caused by rainfall, and the influence of reservoir water level was small.
(2) Stage 2 (9 May-12 January): The displacement increases sharply at the early time, while increased slowly in the subsequent period. The water level of the reservoir fluctuated between 145 m and 175 m. The landslide displacement corresponding to the peak water level of the reservoir remained relatively stable, which proved once again that the rising water level of the reservoir had a weak effect Water 2020, 12, 1860 9 of 22 on the landslide displacement. Around May 2009, the area of the landslide was subjected to continuous heavy rainfall, and the water level of the reservoir reached the lowest value. The largest step-wise deformation in the landslide occurred, which could be considered as a result of the combined effect of reservoir water decline and high precipitation.
(3) Stage 3 (12 January-13 December): In this stage, displacement increased rapidly. The reservoir water level had the same fluctuation rule, but the seasonal rainfall generally intensified. When the water level of the reservoir reached the peak point, the displacement of the monitoring point did not change. However, under the combined effect of the peak precipitation and the lowest reservoir water level in May each year, the landslide displacement increased rapidly and the step-wise characteristics were strengthened.

Prediction of Trend Displacement
According to the above discussion, the displacements of observation points GPS-6 and GPS-3 are selected as the study object. The trend displacement of the Caojiatuo landslide is extracted by using the moving average method (Equations (3) and (15)). The data from January 2008 to January 2013 are used as the training set to fit the equation, and the data from January 2013 to December 2013 are used to verify the accuracy of the fitting equation. The reservoir water level scheduling and GPS data in the TGR area all take a period of one month, that is, the movement period is set as n = 12, which represents the 12 months of each year. Based on the time series addition principle, the trend displacement is removed from the monitoring accumulative displacement of the landslide to obtain the periodic displacement. The extraction results of the trend and period displacement of the observation points GPS-6 and GPS-3 are shown in Figures 6 and 7, respectively.    The trend displacement is controlled by the geological conditions of the landslide and represents the main trend of landslide deformation. The trend displacement of GPS-6 and GPS-3 data are divided into four and three parts respectively, and the cubic polynomial based on the robust weighted least square method is used to fit the trend displacement on the time axis (Equation (15)). The calculated results of each coefficient of the fitting function are shown in Table 1. The fitting results are shown in Figures 8 and 9. The trend displacement of the fitting equations is close to the measured value.
where f (t) means the predictive values of the trend displacement at the time t; and a, b, c, d are the coefficients, where a cannot be equal to zero.

The External Influencing Factors
The prediction of periodic displacement is the key of landslide displacement prediction. In this model, the factors affecting the periodic displacement are taken as the input sequence and the periodic displacement as the output sequence. The selection of influencing factors directly affects the accuracy of model training. According to the above data analysis, the reservoir water level and precipitation are the two main factors affecting the periodic displacement.
Precipitation is one of the main external triggers of landslide deformation in the TGR area. On the one hand, the saturated soil is formed by water infiltration into the landslide, which increases the bulk density of the soil and reduces the stability of the landslide. On the other hand, the rain will induce the chemical reaction with the substance in the landslide, resulting in argillization and softening, which will further induce the occurrence of the landslide. Figure 10 shows the relationship between the periodic displacement of GPS-6 and GPS-3 and precipitation in the last one or two months. The precipitation in the last one or two months has a certain influence on the periodic displacement of the observation point. Both of them show fluctuation change, and the peak of periodic displacement generally lags behind the peak of precipitation. Therefore, the one-month and two-month cumulative antecedent rainfall are selected as the influencing factors of precipitation on the periodic displacement.  The periodic fluctuation of reservoir water level is another major factor leading to the step-wise deformation of the landslide. The change of reservoir water level will affect the distribution of groundwater dynamic field and the seepage field inside the landslide mass, which induce the deformation and damage of the landslide [31]. Figure 11 shows the relationship between the periodic fluctuation of the reservoir water level and the periodic displacement at the observation points GPS-6 and GPS-3. It can be clearly seen that the periodic fluctuation of water level has a certain correlation with the periodic displacement. In particular, when the water level drops sharply, the landslide deformation will be seriously affected. The effect of reservoir water level is usually a slow process, and there will be a so-called "hysteresis effect". Therefore, reservoir level change in one-month and two-month period as well as the average elevation of reservoir level in the current month are selected as the influencing factors of reservoir water level on periodic displacement. The same landslide affected by the same external triggers (precipitation, reservoir water level, etc.) may show completely different evolution characteristics in different stages. In the initial stable state, the landslide may remain stable even under the action of a strong external force. On the contrary, if the landslide is unstable at the beginning, or even in a highly unstable state, a relatively weak external force will eventually lead to great damage [32]. Therefore, this paper proposes to take the evolution state as a supplement to the factors of precipitation and reservoir water level to characterize the influence of other periodic factors on the landslide displacement. The landslide displacements over the past one, two and three months are selected respectively to represent the current evolution state of the landslide [33] as another three influencing factors.
Based on the above analysis, eight influencing factors related to periodic displacement are obtained. The correlation between periodic displacement and external triggering factors is analyzed based on the gray relational grade (GRG) theory [34][35][36][37][38][39][40][41][42][43][44][45][46][47]. According to the GRG theory, when the discrimination coefficient is 0.5, if the relational degree value corresponding to the influence factor is greater than or equal to 0.8, it means that the two are closely related. Table 2 lists the input terms of periodic displacement and their corresponding GRG values. It can be seen that the GRG values of the eight influencing factors selected in this paper are all greater than 0.8, indicating that there is a large correlation between the periodic displacement and the selected input items. Figure 12 intuitively reflects the GRG relational grade between the accumulative displacement and the evolution stage of the landslide in recent months. It can be seen that the GRG values corresponding to the landslide displacement of the one, two and three months selected as the input terms are all greater than 0.8, indicating a great correlation with the periodic displacement. Through the GRG theory, the data shown in the table and figure verify the rationality of the parameter selection in this paper.

Prediction by GWO-ELM Model
After selecting appropriate influencing factors as input items, GWO-ELM model is used to predict the periodic displacement. The specific process of prediction is as follows: (1) Data preprocessing. Normalize various influencing factors and periodic displacements to the interval of [0, l] according to different dimensions. In this way the influence of data dimension could be eliminated.
(2) Model parameter initialization. Based on the analysis in this paper, select the one-month and two-month cumulative antecedent rainfall, reservoir level change in one-month and two-month period as, the average elevation of reservoir level, the landslide displacements over the past one, two and three months as the influencing factors. The number of nodes in the input layer is therefore set to 8. The input item of the model is the predicted value of periodic displacement. The number of node in the output layer is therefore set to 1. The other model parameter settings are shown in Table 3. Table 3. Model parameter initialization settings. (3) Model parameter optimization. Apply the GWO algorithm to search the parameters of the ultimate learning machine model, and then obtain the optimal connection weight matrix between the input layer and the hidden layer, as well as the threshold value of the neuron in the hidden layer. The number of wolf groups in GWO-ELM is set to 30 and the number of iterations is 220. The optimal penalty factor and kernel function parameters can be obtained through the GWO algorithm. In order to demonstrate the superiority of GWO-ELM model, the nonoptimized ELM model is introduced in this paper for comparison under the condition of the same parameter setting. The parameter settings for the two types of machine learning models are shown in Table 3.

Number of Neurons
(4) Training model and prediction. Construct the prediction model by using the parameters obtained from optimization. The model is used to predict the periodic displacement of the landslide. Finally, the prediction results of the periodic displacement of each model are obtained.

Analysis of Prediction Results
Based on the data of Caojiatuo landslide, this paper selects the data of GPS-6 and GPS-3 observation points from January 2008 to January 2013 as the training and validation set of the model. The remaining data from January 2013 to December 2013 are used as a test set to measure the accuracy of the model. The training results of GPS-6 and GPS-3 corresponding to the two models are shown in Figures 13 and 14, respectively. It can be seen that both ELM and GWO-ELM model can learn the fluctuation characteristics of data well. It is proven that the training process of the model is correct by the good coincidence between the peak points and the low points on the curves.
After the correct training process, the displacements of the test set are predicted. Figures 15 and 16 show the predicted results of periodic displacement of GPS-6 and GPS-3 corresponding to the two models, respectively. Tables 4 and 5 correspond to the values of prediction accuracy and error in Figures 15 and 16, respectively. The figures and charts indicate that both models can well predict the periodic displacement value of the landslide, and the time corresponding to the extreme displacement is consistent with the observed one. Moreover, the results of the GWO-ELM prediction model are significantly better than that of the nonoptimized ELM model. After the optimization of ELM model by GWO algorithm, the prediction accuracy has been improved significantly. This is because the ELM model applies random assignment in the weighting matrix between the input layer and hidden layer and threshold of hidden layer neurons. A large number of variable parameters make the parameter variables difficult to control and the result is generally not stable. After introducing and taking the advantages of the GWO algorithm (less needed adjusting parameters, strong global search ability), the optimum connection weights matrix and threshold can be effectively extracted. Therefore, the model can achieve better prediction accuracy and data stability.  because the ELM model applies random assignment in the weighting matrix between the input layer and hidden layer and threshold of hidden layer neurons. A large number of variable parameters make the parameter variables difficult to control and the result is generally not stable. After introducing and taking the advantages of the GWO algorithm (less needed adjusting parameters, strong global search ability), the optimum connection weights matrix and threshold can be effectively extracted. Therefore, the model can achieve better prediction accuracy and data stability.   because the ELM model applies random assignment in the weighting matrix between the input layer and hidden layer and threshold of hidden layer neurons. A large number of variable parameters make the parameter variables difficult to control and the result is generally not stable. After introducing and taking the advantages of the GWO algorithm (less needed adjusting parameters, strong global search ability), the optimum connection weights matrix and threshold can be effectively extracted. Therefore, the model can achieve better prediction accuracy and data stability.      Tables 4 and 5, the prediction accuracy improved significantly. At the same time, it should be pointed out that the relative errors of the predicted values of the two models at the points near the bottom of the curve are larger than those at other positions. The absolute error has a similar pattern, but it is not as obvious as the relative error. This is related to the algorithm of ELM. When the predicted point is located near the low value point on displacement curve, the information of the previous time node has obvious influence on the next time node. This will result in a large error near this point. The GWO algorithm can control the error and improve the accuracy of prediction.

Prediction of the Accumulative Displacement
The predicted cumulative displacement can be obtained by superposing the predicted trend displacement with the predicted periodic displacement. The monitoring data of the GPS-6 and GPS-3 observation points in the Caojiatuo landslide and the predicted results of the GWO-ELM model are shown in Figures 17 and 18, respectively. The predicted trend displacement and periodic displacement are highly accurate, so the predicted values of accumulative displacement are also highly consistent with the measured ones. shown in Figures 17 and 18, respectively. The predicted trend displacement and periodic displacement are highly accurate, so the predicted values of accumulative displacement are also highly consistent with the measured ones.

Conclusions
The displacement of the Caojiatuo landslide is characterized by step-wise deformation. Based on the analysis of the monitoring data, it is concluded that the precipitation and water level fluctuation in the reservoir area are the key external factors affecting the landslide deformation. External triggering factors cause changes in various hydraulic states of the landslide mass, resulting in sudden and large displacement. Over time, the effects of such factors tend to be stable. Thus, it shows the step-wise-characterized deformation.
Based on the time series theory, the step-wise displacement of the Caojiatuo landslide in the TGR area is divided into trend and periodic displacement by using the moving average method, considering the influence of the internal basic geological conditions (lithology, geological structure, landform and geomorphology, etc.) and the external environmental factors (precipitation, reservoir water level, etc.). The trend and periodic displacement are predicted by polynomial fitting method and GWO-ELM model respectively. These processes have clear physical meaning and is proved to be an effective analytical method.
The GWO-ELM model based on time series theory proposed in this paper can accurately predict the displacement of the Caojiatuo landslide at different stages, such as the creep stage and step deformation stage. It shows a better prediction ability than the nonoptimized ELM model. GWO-ELM model has great application potential in landslide displacement prediction in landslide prone areas such as TGA, and can be used to establish an effective early warning system for geological disasters.