# Application of Grey Deformation Prediction Model Optimized by Double Coefficient for Tailings Dam

Authors: Cunji Chu; Gangnian Xu
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

As a nonuniform structure, tailings dam undergoes complex and chaotic nonlinear changes, under the joint influence of multiple dynamic or nondynamic factors. These changes make it difficult to predict the deformation of tailings dam accurately with a numerical model. To solve the problem, this paper proposes a grey deformation prediction model optimized by double coefficient (GDPM-DC). Unlike a single grey prediction model, the GDPM-DC does not mutate significantly but adapts well to specific scenarios. Besides, the model can smoothen and stabilize the original data and thus achieve accurate prediction of the deformation of tailings dam. The main results are as follows. The GDPM-DC made more accurate predictions than the traditional grey model (1, 1) (GM (1, 1)), the grey model based on logarithmic transformation (GM-LT), and the grey model optimized by weight coefficient (GM-WC). It significantly improved the overall prediction accuracy of the vertical and transverse deformations of the dam and controlled the relative error of the predicted seepage pressure to 2.79%–3.71%. Moreover, the model could forecast the trend component and random fluctuation component of seepage pressure effectively, fit the increasing trend in stages 1–3 and the decreasing trend in stages 3–9, and realize the quantitative prediction of deformation law for the operating tailings dam. The research results provide a meaningful reference for the instability analysis and safety management of tailings dam.


---

## Body

## 1. Introduction

The number, capacity, and height of tailings ponds in China, a large mineral producer, are among the highest in the world. The boom of mining produces a huge amount of waste rocks and tailings, which must be properly stored and treated [1]. Tailings are considered a major hazard source with high potential energy. The failure of the tailings dam often causes heavy casualties, huge property losses, and irreparable environmental pollution [2]. In recent years, the incidence of tailings dam accidents generally shows a downward trend, thanks to the economic development and technological innovation. However, major accidents become even more frequent. The statistics of 1910–2010 show that 55% of the major dam break accidents occurred after 1990, and 74% of the dam break accidents after 2000 belong to major or particularly major accidents [3, 4].According to the surveys on the causes of tailings dam failure, the problems of seepage stability (dam displacement and change of saturation line) account for 1/3 of all accidents. Other vulnerability factors may also bring the excessive accumulation of pore pressure of dam, which in turn reduces the effective stress and shear strength of tailings [5, 6]. The factors affecting the seepage stability of the tailings dam are often coupled with each other. It is difficult to accurately describe the nonlinear coupling by traditional calculation methods [7]. Finite-element modeling also faces great limitations in this respect because it requires complex preparations: selection of professional calculation software, establishment of a three-dimensional (3D) model, grid meshing, and parameter setting [8]. Besides, the safety evaluation of the tailings dam is largely qualitative. The most popular approach is the safe checklist. The theoretical evaluation is rarely consistent with the actual situation. To better process nonlinear information, it is necessary to develop an accurate calculation method to analyze the existing data, quantify the features of dam deformation, and simplify the risk evaluation of the tailings dam [9].In recent years, various prediction models have been examined and applied widely to the stability analysis of the tailings dam. The typical models include hybrid computational intelligence approach, backpropagation (BP) neural network, and grey system model [10, 11]. Among them, the hybrid computational intelligence approach has been widely used in landslide displacement prediction of dam due to its advantages of speed and precision. The computational intelligence-based models include but are not limited to extreme learning machine (ELM), support vector regression (SVR), and AI model [12]. Ma et al. applied artificial intelligence (AI) to the modeling of dam landslide movement processes because of its advantages of providing excellent generalization ability and accurately describing complex and nonlinear behavior, and the mutual information (MI)-based measures are proposed for input variable selection (IVS) and incorporated into optimized support vector regression (SVR) for the displacement prediction of seepage-driven landslides [13]. Zhang et al. proposed a hybrid displacement prediction model based on time series theory and various intelligent algorithms to study the effect of frequency components and verified the effectiveness of the hybrid displacement prediction method based on CEEMD reconstruction and DTW-ACO-SVR model [14]. Then, in order to consider the seasonal variations of the rainfall and reservoir level, an optimal combination of methods in which the frequency components of decomposed inducing factors are considered for predicting dam landslide displacement was proposed. Combining CEEMD-LCSS reconstruction and ABC-SVR yields a robust model with remarkable prediction performance [15]. Therefore, the hybrid computational intelligence approach is mostly a hybrid prediction model based on SVR method, which predicts the overall displacement development trend in the long-term period in the future by identifying key variables, extracting high-frequency inducing factors such as the trend and periodic components, and constructing the density distribution of the landslide displacement prediction. Although it shows good performance in the face of small sample data, it usually needs to analyze and process the monitoring data in a long period of time. Compared with other models, it is a prominent problem in specific applications. The application of backpropagation (BP) neural network in tailings dam mostly focuses on the depth prediction of saturation line. Lang et al. used the combination of principal component analysis method and BP network to predict saturation line and used the extracted principal component analysis method to analyze the correlation between input variables but ignored the problem that BP network is easy to fall into local optimal solution [16]. Dai et al. introduced the Pearson correlation coefficient and variable combination method to eliminate the data redundancy between input variables and built a saturation line prediction model based on principal component analysis (PCA) and long-term and short-term memory (LSTM) neural network [17]. Zheng et al. combined with the advantages of particle swarm optimization in network optimization proposed to dynamically optimize particle swarm parameters, adjust the weight and threshold of BP network, and establish a groundwater level prediction model based on IPSO-BP neural network [18]. Zhang et al. used the orthogonal experimental design method and the good nonlinear mapping function of neural network and established a prediction model for rapid calculation of safety factor of tailings pond based on improved neural network [19]. The backpropagation (BP) neural network method usually takes the back slope distance, dry beach length, water level elevation, seepage flow parameters, etc. as the input of the network. Although the influence of several influencing factors on the saturation line is considered, it depends on the sample size and training times, and the operation is cumbersome in the practical application of the project. At present, the application of the grey system model in tailings dam is still in the development stage. It mainly predicts the vertical and transverse deformation of the dam. Wei et al. analyzed the influence of the length of the settlement monitoring data column, the unequal interval of the monitoring time, and the accumulation mode of the data on the prediction accuracy and proposed a method of dynamically predicting the settlement deformation process by using the grey equal dimension prediction model [20]. Ren et al. optimized the original data, initial value, and background value, respectively, and proposed an improved GM (1, 1) prediction model, which significantly improved the prediction accuracy [21]. He et al. used Markov model to deal with the randomness of time series and overcome the limitations of high volatility series and constructed a new dimensional unbiased grey Markov combination prediction model [22].Dam deformation is a complex and chaotic nonlinear process. It is important to select a prediction model that can forecast the deformation of the tailings dam properly, with a high accuracy. The grey prediction model has the characteristics of small sample size, short period, and verifiability. It can be modeled by continuous differential equation [23]. Compared with the hybrid computational intelligence approach and BP neural network prediction method, the grey prediction method does not need to consider the interaction of multiple factors, and the prediction model can be established only by extracting the deformation monitoring data of a certain item. However, this model has a strict requirement on the curve form of the original data. If the data are highly oscillating or abnormal, grey prediction will have a poor applicability and a low accuracy [24]. In order to solve this problem, it is necessary to establish a prediction model that can simultaneously predict the vertical deformation, transverse deformation, and saturation line of tailings dam. Therefore, the weight coefficient is optimized and combined with the Markov chain theory to propose the grey deformation prediction model optimized by double coefficient (GDPM-DC). By studying the transition probability between states, the next state of dam deformation is predicted, making up for the limitations brought by the data volatility and the traditional GM (1, 1) prediction model. Based on the monitoring data of a tailings dam in Xinjiang, this paper proposes the GDPM-DC and uses the model to establish the relationship between dam deformation and grey theory. The GDPM-DC has the advantages of less sample data, more convenient operation, and higher short-term prediction accuracy. It can overcome the defects of local optimization and overall poor prediction sequence and solve the problem that it is difficult for tailings dam management units to judge the stability of tailings dam in time during operation, so as to explore the superiority and applicability of the model in practical projects.


## 2. Preliminaries

According to the theory on grey model (1, 1) (GM (1, 1)) [25, 26], each period of dam deformation was regarded as the original time sequence. The vertical deformation, transverse deformation, and seepage pressure of the dam were measured at the monitoring points and viewed as the original equidistant sequence. The measured data were sorted into a regular data sequence and used to establish a grey prediction model of tailings dam deformation. The model was applied to solve and correct the subsequent prediction results.


### 2.1. GM (1, 1)

GM (1, 1) is composed of a first-order differential equation with only one variable. The discrete sequence needs to be preprocessed before being handled by this model. The original sequence of the monitoring data of the tailings dam can be expressed as(1)X0k=X01,X02,…,X0n.Using the 1-accumulating generation operator (1-AGO), the original sequence can be converted into sequenceX1k:(2)X1k=X11,X12,…,X1n.X1k can also be expressed as(3)X1k=∑i=1kX0i=X0k+X1k−1,k=1,2,…,n.Through accumulative calculation, the original irregular sequence can be transformed into an incremental sequence. The accumulative sequenceX1k is a continuous function of t, and the whitening equation of GM (1, 1) can be established as(4)dX1tdt+aX1t=u,where a is the development grey number and b is the grey action quantity. The integral of formula (4) on k−1,k can be converted into(5)X1k−X1k−1+a∫k−1kX1tdt=u,k=2,3,…,n.In GM (1, 1),a∫k−1kX1tdt can be approximated as a/2X1k+X1k−1 [27]. Thus, the background value Z1k of GM (1, 1) can be generated:(6)Z1k=12X1k+X1k−1.After discretizing the whitening equation, formulas (5) and (6) are compared to obtain the grey differential equation:(7)X0k+aZ1k=u.By the least squares method,a and u can be calculated by [28](8)a^=a,uT=BTB−1BTYn,where data matrix B and data column Yn can be expressed as(9)B=−Z121−Z131⋮⋮−Z1n1,Yn=X02,X03,…,X0nT.Under the initial condition ofX11=X01, the response function of the grey differential equation can be obtained as(10)X^1k=X01−uae−ak−1+ua.At this time, theX^1k value predicted by GM (1,1) is a 1-AGO sequence.Then, the 1-inverse accumulating generation operator (IAGO) is employed to obtain the restored valueX^0k:(11)X^0k=X^1k−X^1k−1,k=2,3,…,n.However, GM (1, 1) requires the original sequence to be highly smooth and has a good effect in predicting the monotonically changing and exponentially increasing sequences. After the original sequence goes through 1-AGO, the resulting sequence will exhibit an exponential growth trend, but the smoothness improvement is rather limited.Currently, many scholars convert the original sequence via logarithmic transformation [29–31], that is, make Y0k=lnX0k and then calculate formulas (2)–(11). The transformation can enhance the smoothness of the original sequence and improve the prediction accuracy on the sequence with a moderate growth trend. Nevertheless, the grey model based on logarithmic transformation (GM-LT) needs to restore the new sequence generated by function transformation. The restoration may suppress the prediction accuracy.


### 2.2. Optimization and Selection of Weight Coefficient

Formula (9) shows that the prediction accuracy of GM (1,1) is closely related to parameters a and u, which depend on the background value Z1k of the model. Therefore, the prediction accuracy can be improved by optimizing the background value of the model [32].In GM (1, 1),Z1k=1/2X1k+X1k−1. The background value, as the area of a trapezoid, can be expressed as(12)Z1k=a∫k−1kX1tdt.If sequenceX1k changes gently, formulas (6) and (12) will have a small error. If the sequence changes and fluctuates violently, the model will face a high lag error and deviation. The common ways to optimize the background value of the model include interval area reconstruction, optimal generation coefficient reconstruction, and integral reconstruction [33].Based on the integral mean value theorem, this paper introduces the dynamic weight coefficientα of the background value. When α∈0,1, the background value can be expressed as(13)Z1k=αX1k+1−αX1k−1.GM (1, 1) assumes that no sudden variable exists in the interval fromX1k−1 to X1k, or such a variable is so small as to be negligible. Thus, the current value and the subsequent value are of equal weights. In this case, the value of α in formula (13) is 0.5. However, the error between the actual and modeled background values cannot be ignored. The value of α also varies with calculation models. Therefore, α needs to be optimized to reduce the source of the error.This paper tries to minimize the mean relative error (MRE) between the measured valueX0k and the predicted value X^0k:(14)MRE=1n−1∑k=2nΔk=1n−1X0k−X^0kX0k.The value ofα is increased at a step length of Δα=0.01 in the interval α∈0,1 until α=1. The MRE at each α value is calculated on MATLAB. Through iterative calculation, the α value corresponding to the minimize MRE is selected as the optimal weight coefficient α∗. Next, α∗ is substituted into formula (13) to obtain a new value. On this basis, the grey model is optimized by weight coefficient and used to make predictions.


## 2.1. GM (1, 1)

GM (1, 1) is composed of a first-order differential equation with only one variable. The discrete sequence needs to be preprocessed before being handled by this model. The original sequence of the monitoring data of the tailings dam can be expressed as(1)X0k=X01,X02,…,X0n.Using the 1-accumulating generation operator (1-AGO), the original sequence can be converted into sequenceX1k:(2)X1k=X11,X12,…,X1n.X1k can also be expressed as(3)X1k=∑i=1kX0i=X0k+X1k−1,k=1,2,…,n.Through accumulative calculation, the original irregular sequence can be transformed into an incremental sequence. The accumulative sequenceX1k is a continuous function of t, and the whitening equation of GM (1, 1) can be established as(4)dX1tdt+aX1t=u,where a is the development grey number and b is the grey action quantity. The integral of formula (4) on k−1,k can be converted into(5)X1k−X1k−1+a∫k−1kX1tdt=u,k=2,3,…,n.In GM (1, 1),a∫k−1kX1tdt can be approximated as a/2X1k+X1k−1 [27]. Thus, the background value Z1k of GM (1, 1) can be generated:(6)Z1k=12X1k+X1k−1.After discretizing the whitening equation, formulas (5) and (6) are compared to obtain the grey differential equation:(7)X0k+aZ1k=u.By the least squares method,a and u can be calculated by [28](8)a^=a,uT=BTB−1BTYn,where data matrix B and data column Yn can be expressed as(9)B=−Z121−Z131⋮⋮−Z1n1,Yn=X02,X03,…,X0nT.Under the initial condition ofX11=X01, the response function of the grey differential equation can be obtained as(10)X^1k=X01−uae−ak−1+ua.At this time, theX^1k value predicted by GM (1,1) is a 1-AGO sequence.Then, the 1-inverse accumulating generation operator (IAGO) is employed to obtain the restored valueX^0k:(11)X^0k=X^1k−X^1k−1,k=2,3,…,n.However, GM (1, 1) requires the original sequence to be highly smooth and has a good effect in predicting the monotonically changing and exponentially increasing sequences. After the original sequence goes through 1-AGO, the resulting sequence will exhibit an exponential growth trend, but the smoothness improvement is rather limited.Currently, many scholars convert the original sequence via logarithmic transformation [29–31], that is, make Y0k=lnX0k and then calculate formulas (2)–(11). The transformation can enhance the smoothness of the original sequence and improve the prediction accuracy on the sequence with a moderate growth trend. Nevertheless, the grey model based on logarithmic transformation (GM-LT) needs to restore the new sequence generated by function transformation. The restoration may suppress the prediction accuracy.


## 2.2. Optimization and Selection of Weight Coefficient

Formula (9) shows that the prediction accuracy of GM (1,1) is closely related to parameters a and u, which depend on the background value Z1k of the model. Therefore, the prediction accuracy can be improved by optimizing the background value of the model [32].In GM (1, 1),Z1k=1/2X1k+X1k−1. The background value, as the area of a trapezoid, can be expressed as(12)Z1k=a∫k−1kX1tdt.If sequenceX1k changes gently, formulas (6) and (12) will have a small error. If the sequence changes and fluctuates violently, the model will face a high lag error and deviation. The common ways to optimize the background value of the model include interval area reconstruction, optimal generation coefficient reconstruction, and integral reconstruction [33].Based on the integral mean value theorem, this paper introduces the dynamic weight coefficientα of the background value. When α∈0,1, the background value can be expressed as(13)Z1k=αX1k+1−αX1k−1.GM (1, 1) assumes that no sudden variable exists in the interval fromX1k−1 to X1k, or such a variable is so small as to be negligible. Thus, the current value and the subsequent value are of equal weights. In this case, the value of α in formula (13) is 0.5. However, the error between the actual and modeled background values cannot be ignored. The value of α also varies with calculation models. Therefore, α needs to be optimized to reduce the source of the error.This paper tries to minimize the mean relative error (MRE) between the measured valueX0k and the predicted value X^0k:(14)MRE=1n−1∑k=2nΔk=1n−1X0k−X^0kX0k.The value ofα is increased at a step length of Δα=0.01 in the interval α∈0,1 until α=1. The MRE at each α value is calculated on MATLAB. Through iterative calculation, the α value corresponding to the minimize MRE is selected as the optimal weight coefficient α∗. Next, α∗ is substituted into formula (13) to obtain a new value. On this basis, the grey model is optimized by weight coefficient and used to make predictions.


## 3. Model Construction

The Markov prediction theory has two preconditions: nonaftereffect and stability. To meet the two preconditions, the Markov model is adopted to process the data sequence with GM (1, 1). The model can predict a highly volatile random sequence [34].


### 3.1. Establishment of State Transition Matrix

#### 3.1.1. Step 1: Segmentation of Original Sequence

The original sequence is often unsmooth, with an uncertain change trend. To overcome these defects, the original sequence is modeled in a piecewise manner on the basis of optimizing the weight coefficient. Firstly,X0k=X01,X02,…,X0n is divided into m subsequences of different number of samples [35]:(15)Xi0=Xi0k|i=1,2,…,m;k=i,i+1,…,n,where m≤n−4. The subsequences of different number of samples help to enhance the correlation between mutation data. Taking each subsequence as the original sequence, the optimal weight coefficient α∗i is solved for each subsequence. The predicted value of each subsequence can be restored by repeating formulas (7)∼(12). Then, the mean of the item corresponding to the predicted value of each subsequence can be obtained. The combined prediction sequence can be used to construct the residual model:(16)X^0k=X^1011,X^102+X^2022,L,X^10n+X^20n⋯X^m0nm.


#### 3.1.2. Step 2: Calculation of Residual Sequence

The basic idea of the Markov model is to derive the state transition matrix from the original data sequence and predict the future trend based on the matrix. The core concepts are “state” and “state transition.”To obtain the residual sequence, it is necessary to solve the difference between the original sequence and the prediction sequence:(17)ε0k=X0k−X^0k,where k=1,2,3,⋯,n. The relative error εk between the predicted value of GM (1, 1) and the measured value can be calculated by(18)εk=X0k−X^0kX0k=ε0kX0k.The relative error can be written asεk=ε1,ε2,⋯,εn.


#### 3.1.3. Step 3: Definition of States

The state transition of Markov chain focuses on setting up the state transition matrix and classifying the states. During state classification, it is necessary to compute the meanμ and mean squared error (MSE) σ of the relative error sequence εk and then divide the sequence into four regions: minε,μ−0.5σ,μ−0.5σ,μ, μ,μ+0.5σ, and μ+0.5σ,maxε. Each region constitutes a state interval Si. If εk∈a1i,a2i, the state interval of εk is Si, where i=1,2,3,4; a1i and a2i are the upper and lower bounds of the state, respectively.


#### 3.1.4. Step 4: Construction of State Transition Matrix

After gathering the statistics onSi, the number of data in each state interval is assumed as Mi, and the number of transitions from Si to Sj through K steps is assumed as Nijk. Then, the state transition probability can be calculated by(19)Pijk=NijkMi.Then, thePijk of each state is arranged to form the state transition matrix [36]:(20)Pk=p11kp12k⋯p1nkp21kp22k⋯p2nk⋮⋮⋮pn1kpn2k⋯pnnk,where ∑j=1nPijk=1.


### 3.2. Calculation of Autocorrelation Coefficient

The autocorrelation coefficient is applied and optimized based on the Markov prediction model. The principle is to assign a weight to each predicted value and find the weighted sum of the prediction probabilities in the same state. During the calculation, the first is to perform the autocorrelation analysis ofεk. At this time, εk can be regarded as the index value xk with K periods, and the autocorrelation coefficient with a lag of K is defined as(21)rK=∑k=1n−Kxk−x¯xk+K−x¯∑k=1nxk−x¯2,where rK is the K-th order autocorrelation coefficient, K∈S; xk is the index value of the k-th period; x¯ is the mean of the index value; and n is the length of the index sequence.The above autocorrelation analysis shows that the greater theK-th order autocorrelation coefficient, the more stable the change of the system state illustrated by the state transition matrix and the larger the weight to be assigned [37]. Then, the autocorrelation coefficient of each order is normalized to obtain the weight wK:(22)wK=rK∑K=1SrK.The autocorrelation coefficient can weaken or eliminate the lag effect of the prediction sequence and reflect the relationship of several subsequent sequences, in addition to that of the current sequence.


### 3.3. Optimization of Residual Model

The next period data ofX0k are predicted, using autocorrelation coefficient of the prediction sequence. This predicted value can be called extrapolated value Y^0k+1. First, choose the K periods that are closest to the extrapolated value; the state corresponding to the relative error of each period is taken as the initial state of that period; the number of transfer steps is set as 1,2,⋯,K in the order of the proximity to the extrapolated value. From the state transition matrix corresponding to a specific number of steps, the row vector Pik=pi1k,pi2k,⋯,pink corresponding to the initial states of the K nearest periods is organized into a new probability matrix as follows:(23)P=p111p121⋯p1n1p212p222⋯p2n2⋮⋮⋮pi1jpi2j⋯pinj.The prediction probabilities in the same state are weighted and summed up to obtain the transition probability of relative error under different step lengths:(24)Pi=∑i=1nwK×Pik.According to the definition of maximum likelihood, the state corresponding to the maximumPi is the state Si of period k+1, i.e., maxPi,i∈S. Then, the relative error of this period falls between the upper and lower bounds of Si. So, the final extrapolated value can be obtained as(25)Y^0k+1=X^0k+1×1−a1i+a2i2.Figure1 details the modeling process of the GDPM-DC.Figure 1
Modeling process of the GDPM-DC.Following the above modeling process, this paper considers two coefficients: the weight coefficient and autocorrelation coefficient. Based on the selection of different coefficients, four models were established to separately predict dam deformation (Table1).Table 1
Four calculation models.
NumberModelWeightcoefficientAutocorrelationcoefficient1GM (1,1)0.5Unreferenced2GM-LT0.5Unreferenced3GM-WCIterative calculationUnreferenced4GDPM-DCIterative calculationAutocorrelation calculationNote. GM-WC is short for the grey model optimized by weight coefficient.The GDPM-DC is a combination and optimization based on GM (1, 1). Firstly, the dynamic weight coefficient is referenced, and the GM-WC is established to improve the smoothness of the sequence. The feasibility of model optimization is illustrated by improving the forecast accuracy. Then, the model is combined with Markov prediction theory to further optimize the accuracy by calculating the autocorrelation coefficient. However, when selecting the double coefficients, the most complex work is to analyze the regular or irregular monitoring data. The two models are reasonably combined by segment modeling of subsequence in order to further improve the prediction accuracy and enhance the random agility.


## 3.1. Establishment of State Transition Matrix

### 3.1.1. Step 1: Segmentation of Original Sequence

The original sequence is often unsmooth, with an uncertain change trend. To overcome these defects, the original sequence is modeled in a piecewise manner on the basis of optimizing the weight coefficient. Firstly,X0k=X01,X02,…,X0n is divided into m subsequences of different number of samples [35]:(15)Xi0=Xi0k|i=1,2,…,m;k=i,i+1,…,n,where m≤n−4. The subsequences of different number of samples help to enhance the correlation between mutation data. Taking each subsequence as the original sequence, the optimal weight coefficient α∗i is solved for each subsequence. The predicted value of each subsequence can be restored by repeating formulas (7)∼(12). Then, the mean of the item corresponding to the predicted value of each subsequence can be obtained. The combined prediction sequence can be used to construct the residual model:(16)X^0k=X^1011,X^102+X^2022,L,X^10n+X^20n⋯X^m0nm.


### 3.1.2. Step 2: Calculation of Residual Sequence

The basic idea of the Markov model is to derive the state transition matrix from the original data sequence and predict the future trend based on the matrix. The core concepts are “state” and “state transition.”To obtain the residual sequence, it is necessary to solve the difference between the original sequence and the prediction sequence:(17)ε0k=X0k−X^0k,where k=1,2,3,⋯,n. The relative error εk between the predicted value of GM (1, 1) and the measured value can be calculated by(18)εk=X0k−X^0kX0k=ε0kX0k.The relative error can be written asεk=ε1,ε2,⋯,εn.


### 3.1.3. Step 3: Definition of States

The state transition of Markov chain focuses on setting up the state transition matrix and classifying the states. During state classification, it is necessary to compute the meanμ and mean squared error (MSE) σ of the relative error sequence εk and then divide the sequence into four regions: minε,μ−0.5σ,μ−0.5σ,μ, μ,μ+0.5σ, and μ+0.5σ,maxε. Each region constitutes a state interval Si. If εk∈a1i,a2i, the state interval of εk is Si, where i=1,2,3,4; a1i and a2i are the upper and lower bounds of the state, respectively.


### 3.1.4. Step 4: Construction of State Transition Matrix

After gathering the statistics onSi, the number of data in each state interval is assumed as Mi, and the number of transitions from Si to Sj through K steps is assumed as Nijk. Then, the state transition probability can be calculated by(19)Pijk=NijkMi.Then, thePijk of each state is arranged to form the state transition matrix [36]:(20)Pk=p11kp12k⋯p1nkp21kp22k⋯p2nk⋮⋮⋮pn1kpn2k⋯pnnk,where ∑j=1nPijk=1.


## 3.1.1. Step 1: Segmentation of Original Sequence

The original sequence is often unsmooth, with an uncertain change trend. To overcome these defects, the original sequence is modeled in a piecewise manner on the basis of optimizing the weight coefficient. Firstly,X0k=X01,X02,…,X0n is divided into m subsequences of different number of samples [35]:(15)Xi0=Xi0k|i=1,2,…,m;k=i,i+1,…,n,where m≤n−4. The subsequences of different number of samples help to enhance the correlation between mutation data. Taking each subsequence as the original sequence, the optimal weight coefficient α∗i is solved for each subsequence. The predicted value of each subsequence can be restored by repeating formulas (7)∼(12). Then, the mean of the item corresponding to the predicted value of each subsequence can be obtained. The combined prediction sequence can be used to construct the residual model:(16)X^0k=X^1011,X^102+X^2022,L,X^10n+X^20n⋯X^m0nm.


## 3.1.2. Step 2: Calculation of Residual Sequence

The basic idea of the Markov model is to derive the state transition matrix from the original data sequence and predict the future trend based on the matrix. The core concepts are “state” and “state transition.”To obtain the residual sequence, it is necessary to solve the difference between the original sequence and the prediction sequence:(17)ε0k=X0k−X^0k,where k=1,2,3,⋯,n. The relative error εk between the predicted value of GM (1, 1) and the measured value can be calculated by(18)εk=X0k−X^0kX0k=ε0kX0k.The relative error can be written asεk=ε1,ε2,⋯,εn.


## 3.1.3. Step 3: Definition of States

The state transition of Markov chain focuses on setting up the state transition matrix and classifying the states. During state classification, it is necessary to compute the meanμ and mean squared error (MSE) σ of the relative error sequence εk and then divide the sequence into four regions: minε,μ−0.5σ,μ−0.5σ,μ, μ,μ+0.5σ, and μ+0.5σ,maxε. Each region constitutes a state interval Si. If εk∈a1i,a2i, the state interval of εk is Si, where i=1,2,3,4; a1i and a2i are the upper and lower bounds of the state, respectively.


## 3.1.4. Step 4: Construction of State Transition Matrix

After gathering the statistics onSi, the number of data in each state interval is assumed as Mi, and the number of transitions from Si to Sj through K steps is assumed as Nijk. Then, the state transition probability can be calculated by(19)Pijk=NijkMi.Then, thePijk of each state is arranged to form the state transition matrix [36]:(20)Pk=p11kp12k⋯p1nkp21kp22k⋯p2nk⋮⋮⋮pn1kpn2k⋯pnnk,where ∑j=1nPijk=1.


## 3.2. Calculation of Autocorrelation Coefficient

The autocorrelation coefficient is applied and optimized based on the Markov prediction model. The principle is to assign a weight to each predicted value and find the weighted sum of the prediction probabilities in the same state. During the calculation, the first is to perform the autocorrelation analysis ofεk. At this time, εk can be regarded as the index value xk with K periods, and the autocorrelation coefficient with a lag of K is defined as(21)rK=∑k=1n−Kxk−x¯xk+K−x¯∑k=1nxk−x¯2,where rK is the K-th order autocorrelation coefficient, K∈S; xk is the index value of the k-th period; x¯ is the mean of the index value; and n is the length of the index sequence.The above autocorrelation analysis shows that the greater theK-th order autocorrelation coefficient, the more stable the change of the system state illustrated by the state transition matrix and the larger the weight to be assigned [37]. Then, the autocorrelation coefficient of each order is normalized to obtain the weight wK:(22)wK=rK∑K=1SrK.The autocorrelation coefficient can weaken or eliminate the lag effect of the prediction sequence and reflect the relationship of several subsequent sequences, in addition to that of the current sequence.


## 3.3. Optimization of Residual Model

The next period data ofX0k are predicted, using autocorrelation coefficient of the prediction sequence. This predicted value can be called extrapolated value Y^0k+1. First, choose the K periods that are closest to the extrapolated value; the state corresponding to the relative error of each period is taken as the initial state of that period; the number of transfer steps is set as 1,2,⋯,K in the order of the proximity to the extrapolated value. From the state transition matrix corresponding to a specific number of steps, the row vector Pik=pi1k,pi2k,⋯,pink corresponding to the initial states of the K nearest periods is organized into a new probability matrix as follows:(23)P=p111p121⋯p1n1p212p222⋯p2n2⋮⋮⋮pi1jpi2j⋯pinj.The prediction probabilities in the same state are weighted and summed up to obtain the transition probability of relative error under different step lengths:(24)Pi=∑i=1nwK×Pik.According to the definition of maximum likelihood, the state corresponding to the maximumPi is the state Si of period k+1, i.e., maxPi,i∈S. Then, the relative error of this period falls between the upper and lower bounds of Si. So, the final extrapolated value can be obtained as(25)Y^0k+1=X^0k+1×1−a1i+a2i2.Figure1 details the modeling process of the GDPM-DC.Figure 1
Modeling process of the GDPM-DC.Following the above modeling process, this paper considers two coefficients: the weight coefficient and autocorrelation coefficient. Based on the selection of different coefficients, four models were established to separately predict dam deformation (Table1).Table 1
Four calculation models.
NumberModelWeightcoefficientAutocorrelationcoefficient1GM (1,1)0.5Unreferenced2GM-LT0.5Unreferenced3GM-WCIterative calculationUnreferenced4GDPM-DCIterative calculationAutocorrelation calculationNote. GM-WC is short for the grey model optimized by weight coefficient.The GDPM-DC is a combination and optimization based on GM (1, 1). Firstly, the dynamic weight coefficient is referenced, and the GM-WC is established to improve the smoothness of the sequence. The feasibility of model optimization is illustrated by improving the forecast accuracy. Then, the model is combined with Markov prediction theory to further optimize the accuracy by calculating the autocorrelation coefficient. However, when selecting the double coefficients, the most complex work is to analyze the regular or irregular monitoring data. The two models are reasonably combined by segment modeling of subsequence in order to further improve the prediction accuracy and enhance the random agility.


## 4. Case Analysis

### 4.1. Deformation Prediction of GM-WC

During construction and operation, the tailings dam undergoes vertical and transverse displacement and variation of seepage pressure. The collective result of these processes is the deformation of the tailings dam. Taking a tailings dam in Xinjiang as an example, this paper carries out automatic monitoring of the dam deformation. The monitoring data were processed for prediction, which better reflect the dam stability under seepage. Figure2 shows the monitoring plan for the tailings dam.Figure 2
Monitoring plan for the tailings dam.Six monitoring points were mainly arranged at the abutment and middle of the tailings dam. The long-term monitoring was realized using instruments like level gauge, inclinometer, and osmometer. The monitoring data were compiled into the original sequence of the prediction model (Table2).Table 2
Monitoring data at points 01–06.
Monitoring contentMonitoring pointMonitoring period/day4812162024283236Monitoring stage123456789Vertical deformation/mmSZY-017.411.211.913.214.014.615.015.616.5SZY-026.59.810.411.612.312.813.213.714.5SZY-036.59.39.911.011.712.212.513.013.8SZY-047.811.812.513.914.815.415.816.417.4SZY-056.810.310.912.212.913.413.914.415.3SZY-065.98.89.410.411.111.511.912.313.1Transverse deformation/mmCXY-014.39.911.813.315.516.817.919.821.6CXY-022.39.911.212.714.416.618.820.523.2CXY-034.710.913.014.617.118.519.721.824.0CXY-042.811.913.415.217.319.922.624.628.0CXY-053.47.99.410.612.413.414.315.817.4CXY-062.18.910.111.413.014.916.918.521.1Seepage pressure/kPaSYJ-0153.5078.9994.8089.7082.0577.9775.4270.3364.17SYJ-0261.2869.9582.5069.4466.8964.3459.7556.6950.14SYJ-0358.8586.89104.2898.6790.2685.7782.9677.3670.59SYJ-0473.5483.9499.0083.3380.2777.2171.7068.0360.34SYJ-0542.8063.1975.8471.7665.6462.3860.3456.2651.44SYJ-0655.1562.9674.2562.5060.2057.9153.7851.0244.93The preliminary analysis of the monitoring data (Table2) suggests that the vertical deformation of the dam grew slowly, while the transverse deformation increased sharply in stage 1 and recovered smoothly in stages 2–9. The seepage pressure of the dam exhibited an increasing trend in stages 1–3, reached the peak in stage 3, and declined thereafter.To verify the performance of the prediction model on data of different trends, the data of monitoring point 01 were analyzed in detail. The data collected in the first 8 stages (4–32 days) were taken as the original sequence to build the prediction model, and the data collected in stage 9 (32–36 days) were used as the test sequence to validate the model performance.During the construction of the GM-WC, the top priority is to determine the value ofα∗. The value of α was increased at a step length of Δα=0.01 in the interval α∈0,1 until α=1. In this way, a total of 100 fitted values were obtained. Then, α∗α∗ was selected according to the principle of minimum MRE (Figure 3).Figure 3
Weight coefficient and MRE of predicted value.According to the iterative results in Figure3, α∗ is not simply defined as 0.5 but irregularly distributed on both sides of the centerline. For the vertical deformation prediction model, α∗=0.40, and the corresponding MRE = 0.02190. For the transverse deformation prediction model, α∗=0.45, and the corresponding MRE = 0.02841. For the seepage pressure prediction model, α∗=0.77α∗=0.77, and the corresponding MRE = 0.04978.In addition, for either the GM-WC or GM-LT, the prediction accuracy could be improved by smoothing the original sequence. The predicted value (extrapolated value) in stage 9 was solved by the first 3 models in Table1. Tables 3–5 list the prediction results on dam deformation at monitoring point 01, using different models.Table 3
Predicted vertical deformation of monitoring point 01 using different models.
StageMeasured value/mmGM (1, 1)GM-LTGM-WCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%17.47.4000.007.4000.007.4000.00211.211.551−3.1311.520−2.8611.487−2.56311.912.186−2.4112.139−2.0112.116−1.81413.212.8572.6012.8052.9912.7793.19514.013.5653.1113.5243.4013.4793.72614.614.3111.9814.2992.0614.2162.63715.015.098−0.6615.138−0.9214.9950.04815.615.929−2.1116.044−2.8515.816−1.38916.516.806−1.8517.027−3.1916.681−1.10Table 4
Predicted transverse deformation of monitoring point 01 using different models.
StageMeasured value/mmGM (1,1)GM-LTGM-WCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%14.34.3000.004.3000.004.3000.0029.910.669−7.7710.607−7.1410.609−7.16311.811.865−0.5511.7200.6811.7920.07413.313.1940.7913.0042.2213.1061.46515.514.6735.3414.4936.4914.5666.02616.816.3172.8716.2273.4116.1903.63717.918.146−1.3718.255−1.9817.994−0.53819.820.179−1.9120.638−4.2320.000−1.01921.622.440−3.8923.455−8.5922.229−2.91Table 5
Predicted seepage pressure of monitoring point 01 using different models.
StageMeasured value/mmGM (1, 1)GM-LTGM-WCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%153.5053.5000.0053.5000.0053.5000.00278.9989.421−13.2189.569−13.3988.729−12.33394.8086.5788.6786.5648.6985.9059.38489.7083.8256.5583.6816.7183.1717.28582.0581.1591.0980.9151.3880.5241.86677.9778.578−0.7878.261−0.3777.9620.01775.4276.080−0.8775.713−0.3975.481−0.08870.3373.660−4.7473.266−4.1773.079−3.91964.1771.318−11.1470.916−10.5170.753−10.26Similarly, the dam deformation of monitoring points 02–06 in stage 9 was predicted by different models and compared with the above results (Table6).Table 6
Prediction results (extrapolated value) of monitoring points 02–06 in stage 9.
Monitoring contentMonitoring pointGM (1,1)GM-LTGM-WCRelative error/%Relative error/%Relative error/%Vertical deformation/mmSZY-01−1.85−3.19−1.10SZY-02−1.97−3.39−1.43SZY-03−1.67−3.11−0.84SZY-04−1.71−3.07−1.04SZY-05−1.58−2.97−1.27SZY-06−1.49−2.94−1.18Transverse deformation/mmCXY-01−3.89−8.59−2.91CXY-02−1.51−5.72−1.03CXY-03−2.95−7.53−2.17CXY-04−1.02−5.00−0.31CXY-05−2.99−7.90−1.83CXY-06−0.57−4.860.14Seepage pressure/kPaSYJ-01−11.14−10.51−10.26SYJ-02−10.33−10.04−9.56SYJ-03−11.13−10.50−10.25SYJ-04−10.01−9.71−9.25SYJ-05−10.92−10.31−10.04SYJ-06−10.81−10.54−10.04The following can be derived by comparing Tables3–6.(1)
Among the three models, the GM-WC achieved higher overall prediction accuracy than GM(1, 1) and GM-LT(2)
When the vertical deformation grew gently, all three models had a high prediction accuracy. Their relative errors were controlled within 0.84%∼3.39%. When the transverse deformation fluctuated significantly, the GM-WC witnessed the greatest improvement of prediction accuracy, with a minimum relative error of only 0.14%. In contrast, the GM-LT had a relatively low prediction accuracy, and the prediction effect of the model was not ideal after function recovery. When the original sequence of seepage pressure had peaks, the three models failed to make sufficiently accurate predictions, and their prediction accuracy needs to be further optimized.The predicted values of the three models are compared in Figures4–6 to better illustrate their prediction performance for the trend of monitoring data.Figure 4
Vertical deformation of monitoring point 01 predicted by different models.Figure 5
Transverse deformation of monitoring point 01 predicted by different models.Figure 6
Seepage pressure of monitoring point 01 predicted by different models.The following can be derived by comparing Figures4–6.The three models excelled in fitting the overall trend of vertical deformation and transverse deformation but did not perform so well in predicting the trend of seepage pressure. The predicted values only increased in stages 1–2, peaking in stage 2. This trend is different from the measured trend. In the subsequent stages 2–9, the predicted values declined slower than the measured values. Hence, the goodness of fit of the three models is positively correlated with the smoothness of the measured values. The fitting effects of these models need to be optimized, if there are peaks in the measured values.


### 4.2. Deformation Prediction of GDPM-DC

The above analysis shows that the GM-WC has certain advantages in overall prediction accuracy. But the model faces some defects in predicting volatile data. Therefore, the GDPM-DC was extended from the GM-WC to solve these defects.Firstly, the vertical deformation dataX0=7.4,11.2,11.9,13.2,14.0,14.6,15.0,15.6 of monitoring point 01 were divided into 4 subsequences with different number of samples:(26)X10=7.4,11.2,11.9,13.2,14.0,14.6,15.0,15.6,X20=11.2,11.9,13.2,14.0,14.6,15.0,15.6,X30=11.9,13.2,14.0,14.6,15.0,15.6,X40=13.2,14.0,14.6,15.0,15.6.Next, a GM-WC was established for each subsequence. The GM-WC was used to predict the values of the first 8 stages and obtain the extrapolated value of stage 9 of these subsequences. After that, the mean of the items corresponding to the predicted value of each subsequence can be calculated by(27)X^0=7.400,11.344,12.118,13.063,13.745,14.366,15.016,15.696.The mean extrapolated valueX^09 can be calculated by(28)X^09=16.408.The prediction sequenceX^0, which is calculated in a piecewise manner, was used to establish the GDPM-DC. Firstly, the relative error sequence ε=0,−0.013,−0.018,0.010,0.018,0.016,−0.001,−0.006 was established, after calculating the difference between the mean predicted value and the measured value. By the mean and MSE of the relative error sequence, the sequence was divided into four regions: S1=−0.0180,−0.0058, S2=−0.0058,−0.0008, S3=−0.0008,0.0073, and S4=0.0073,0.0180. The corresponding states are shown in Table 7.Table 7
Residual sequence states of vertical deformation at monitoring point 01.
StageRelative error/%State10.00S32−1.28S13−1.83S141.04S451.82S461.60S47−0.11S28−0.62S1Based on Table7, the transition probability of each state was entered into the state transition matrix P:(29)P1=1/2001/21000100001/302/3,P2=0001000010001/31/301/3,P3=0001000000011/21/200,P4=01/201/2000000011000.Then, the autocorrelation coefficientrK of the relative error sequence was calculated, and the weight wK of each transfer step was obtained after normalization. The calculation results are shown in Table 8.Table 8
Autocorrelation coefficient of each order and weight of each transfer step.
StepKrKwK10.4170.27220.3560.23230.5500.35940.2090.137Based on the relative error and state transition probability of vertical deformation in the first 8 stages, the deformation of the stage 9 was predicted. Besides, the maximumP1=0.45 of the relative error transition probability in stage 9 was derived from the results in Table 8 and the probability matrix in steps 1–4. Table 9 shows the predicted vertical deformation of monitoring point 01 in stage 9.Table 9
Predicted vertical deformation of monitoring point 01 in stage 9.
StageStepKInitial statewKProbability distribution of statesS1S2S3S481S10.2721/2001/272S20.232000063S40.3591/21/20154S40.1371000It can be inferred from Table9 that the relative error in stage 9 peaked in state S1. Then, the extrapolated value of vertical deformation in stage 9 can be calculated:(30)Y^09=16.408×1−−0.018−0.00582=16.603.Similarly, the lateral deformation and seepage pressure of monitoring point 01 were predicted using the above GDPM-DC. The results are summarized in Table10.Table 10
Predicted deformation of monitoring point 01 by the GDPM-DC.
Monitoring contentGDPM-DCStage123456789Vertical deformation/mmMeasured value/mm7.411.211.913.21414.61515.616.5Predicted value/mm7.40011.34412.11813.06313.74514.36615.01615.69616.603Relative error/%0.00−1.28−1.831.041.821.60−0.11−0.62−0.63Transverse deformation/mmMeasured value/mm4.39.911.813.315.516.817.919.821.6Predicted value/mm4.30010.25511.91213.38114.97016.44218.06019.83921.699Relative error/%0.00−3.58−0.95−0.613.422.13−0.89−0.20−0.46Seepage pressure/kPaMeasured value/mm53.578.9994.889.782.0577.9775.4270.3364.17Predicted value/mm53.50084.20692.06887.87482.74578.71674.89271.26366.382Relative error/%0.00−6.602.882.04−0.85−0.960.70−1.33−3.45Figures7–9 compare the predicted values in Table 9 with the results of GM (1,1) and GM-LT.Figure 7
Vertical deformation of monitoring point 01 predicted by different models.Figure 8
Transverse deformation of monitoring point 01 predicted by different models.Figure 9
Seepage pressure of monitoring point 01 predicted by different models.The following can be derived by comparing Figures7–9.The extrapolated value of stage 9 calculated by the GDPM-DC was very close to the measured value. It is more accurate than that calculated by GM (1, 1) and GM-LT. The three models had a smaller difference in predicting the trend of vertical and transverse deformations. However, the GDPM-DC boasted an obvious advantage in predicting the trend of seepage pressure. The prediction results increased in stages 1–3 and declined in stages 3–8. This trend is consistent with the measured values.To better reflect the optimization degree of the GDPM-DC, the prediction results of monitoring points 01–06 are summarized in Table11. The difference of relative error was adopted to measure the difference between the mean value of relative error of GM (1, 1), GM-LT, and GM-WC and the relative error of the GDPM-DC.Table 11
Comparison and verification of stage 9 prediction results (extrapolated value).
Monitoring contentMonitoring pointGM (1, 1)GM-LTGM-WCGDPM-DCDifference ofrelative error/%Relative error/%Relative error/%Relative error/%Relative error/%Vertical deformation/mmSZY-01−1.85−3.19−1.1−0.631.42SZY-02−1.97−3.39−1.43−0.981.28SZY-03−1.67−3.11−0.84−0.291.58SZY-04−1.71−3.07−1.04−0.581.36SZY-05−1.58−2.97−1.27−1.000.94SZY-06−1.49−2.94−1.18−0.651.22Transverse deformation/mmCXY-01−3.89−8.59−2.91−0.464.67CXY-02−1.51−5.72−1.03−0.732.02CXY-03−2.95−7.53−2.17−0.593.63CXY-04−1.02−5.00−0.31−0.351.76CXY-05−2.99−7.90−1.83−0.353.89CXY-06−0.57−4.860.14−0.071.79Seepage pressure/kPaSYJ-01−11.14−10.51−10.26−3.457.19SYJ-02−10.33−10.04−9.56−3.256.73SYJ-03−11.13−10.50−10.25−3.257.38SYJ-04−10.01−9.71−9.25−2.966.70SYJ-05−10.92−10.31−10.04−2.797.63SYJ-06−10.81−10.54−10.04−3.716.75By relative error, the different prediction models could be ranked as GDPM-DC < GM-WC < GM (1, 1) and GM-LT. The GDPM-DC is clearly superior to the other three models. The superiority is particularly notable in predicting seepage pressure: our model reduced the maximum relative error from 11.14% to 3.45%, the improvement of accuracy is the most obvious, and the improvement effect is the best.


### 4.3. Model Validation

In order to fully verify the superiority and applicability of the GDPM-DC, based on the monitoring data of transverse deformation in the first 8 stages of the dam in reference [22] and the monitoring data of accumulative vertical deformation in the first 8 stages of the railway in reference [23], the GDPM-DC is used to predict the deformation in stages 9 and 10.Firstly, the transverse deformation dataX0=6.36,6.94,7.31,7.80,8.03,8.17,8.46,8.64 in the first 8 stages of the dam were divided into 4 subsequences with different number of samples, and the mean extrapolated value X^09 can be calculated by(31)X^09=8.90.According to the division of residual sequence states of transverse deformation, the transition probability of each state was entered into the state transition matrixP:(32)P1=1/2001/2000110001/31/301/3,P2=00011000100001/201/2,P3=01/201/2000000011/2001/2,P4=01/201/2000000011000.Then, through the autocorrelation calculation and the probability matrix in steps 1–4, the maximumP1=0.45 of the relative error transition probability in stage 9 can be calculated (see Table 12 for the specific process).Table 12
Predicted transverse deformation of the dam in stage 9.
StageStepKInitial statewKProbability distribution of statesS1S2S3S481S10.1531/2001/272S40.49001/201/263S20.242000054S40.1151000It can be inferred from Table12 that the relative error in stage 9 peaked in state S4. Then, the extrapolated value of transverse deformation in stage 9 can be calculated:(33)Y^09=8.90×1−0.0029+0.0092=8.80.Similarly, the extrapolated value of transverse deformation in stage 10 can be calculated:(34)Y^010=9.10.In addition, by predicting the accumulative vertical deformation of the railway in stages 9–10, the extrapolated values can be calculated asY^09=3.04 and Y^010=3.14. The above calculation results are summarized in Tables 13 and 14.Table 13
Predicted transverse deformation of the dam using different models.
StageMeasured value/mmReference [22] modelGM-WCGDPM-DCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%16.366.55−2.936.360.006.360.0026.946.821.707.11−2.457.02−1.2837.317.112.817.36−0.707.37−1.8347.807.405.157.622.297.731.0458.037.713.987.891.737.961.8268.178.031.768.170.008.191.6078.468.361.178.460.028.43−0.1188.648.71−0.828.76−1.368.68−0.6298.709.07−4.229.07−4.228.80−1.15108.919.44−5.999.39−5.369.10−2.13Mean relative error/%—5.11—4.79—1.64Table 14
Predicted accumulative vertical deformation of the railway using different models.
StageMeasured value/mmReference [23] modelGM-WCGDPM-DCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%13.043.077−1.223.0400.003.0400.0023.213.0823.993.0166.053.113−1.2832.943.085−4.933.034−3.212.958−1.8342.853.088−8.353.053−7.122.9481.0453.173.0902.523.0723.113.0531.8263.033.092−2.053.090−1.993.0861.6073.283.0935.703.1095.203.118−0.1183.093.094−0.133.128−1.243.152−0.6292.953.095−4.923.148−6.693.040−3.05103.213.0953.583.1671.353.1402.18Mean relative error/%—4.25—4.02—2.62According to the prediction results of the transverse deformation of the dam in stages 9–10, the prediction accuracy difference between reference [22] model and the GM-WC is small, and the mean relative error is about 5%. The GDPM-DC has the highest prediction accuracy in which the relative errors of stages 9 and 10 are −1.15% and −2.13%, respectively, and the mean relative error is 1.64%. Compared with reference [22] model, its value decreases by 3.47%. Moreover, it can be seen from the prediction results of accumulative vertical deformation of the railway in stages 9–10 that the prediction accuracy of the GDPM-DC is still the highest. Compared with reference [23] model, the mean relative error is reduced by 1.63%, and the comparison results of relative errors are plotted as shown in Figures 10 and 11.Figure 10
Comparison of relative error of transverse deformation.Figure 11
Comparison of relative error of accumulative vertical deformation.Through comparison with reference [22, 23] model, it is found that the prediction value calculated by the GDPM-DC is more consistent with the monitoring data, and the prediction accuracy of the model is greatly improved, which verifies that the model has obvious superiority in predicting dam deformation. Moreover, the GDPM-DC can still have high prediction accuracy when predicting the data of accumulative vertical deformation in the railway with great fluctuation, so that the relative error of each stage can maintain a more gentle change trend, which indicates that the GDPM-DC also has a certain application value for processing more complex data of accumulative vertical deformation.


## 4.1. Deformation Prediction of GM-WC

During construction and operation, the tailings dam undergoes vertical and transverse displacement and variation of seepage pressure. The collective result of these processes is the deformation of the tailings dam. Taking a tailings dam in Xinjiang as an example, this paper carries out automatic monitoring of the dam deformation. The monitoring data were processed for prediction, which better reflect the dam stability under seepage. Figure2 shows the monitoring plan for the tailings dam.Figure 2
Monitoring plan for the tailings dam.Six monitoring points were mainly arranged at the abutment and middle of the tailings dam. The long-term monitoring was realized using instruments like level gauge, inclinometer, and osmometer. The monitoring data were compiled into the original sequence of the prediction model (Table2).Table 2
Monitoring data at points 01–06.
Monitoring contentMonitoring pointMonitoring period/day4812162024283236Monitoring stage123456789Vertical deformation/mmSZY-017.411.211.913.214.014.615.015.616.5SZY-026.59.810.411.612.312.813.213.714.5SZY-036.59.39.911.011.712.212.513.013.8SZY-047.811.812.513.914.815.415.816.417.4SZY-056.810.310.912.212.913.413.914.415.3SZY-065.98.89.410.411.111.511.912.313.1Transverse deformation/mmCXY-014.39.911.813.315.516.817.919.821.6CXY-022.39.911.212.714.416.618.820.523.2CXY-034.710.913.014.617.118.519.721.824.0CXY-042.811.913.415.217.319.922.624.628.0CXY-053.47.99.410.612.413.414.315.817.4CXY-062.18.910.111.413.014.916.918.521.1Seepage pressure/kPaSYJ-0153.5078.9994.8089.7082.0577.9775.4270.3364.17SYJ-0261.2869.9582.5069.4466.8964.3459.7556.6950.14SYJ-0358.8586.89104.2898.6790.2685.7782.9677.3670.59SYJ-0473.5483.9499.0083.3380.2777.2171.7068.0360.34SYJ-0542.8063.1975.8471.7665.6462.3860.3456.2651.44SYJ-0655.1562.9674.2562.5060.2057.9153.7851.0244.93The preliminary analysis of the monitoring data (Table2) suggests that the vertical deformation of the dam grew slowly, while the transverse deformation increased sharply in stage 1 and recovered smoothly in stages 2–9. The seepage pressure of the dam exhibited an increasing trend in stages 1–3, reached the peak in stage 3, and declined thereafter.To verify the performance of the prediction model on data of different trends, the data of monitoring point 01 were analyzed in detail. The data collected in the first 8 stages (4–32 days) were taken as the original sequence to build the prediction model, and the data collected in stage 9 (32–36 days) were used as the test sequence to validate the model performance.During the construction of the GM-WC, the top priority is to determine the value ofα∗. The value of α was increased at a step length of Δα=0.01 in the interval α∈0,1 until α=1. In this way, a total of 100 fitted values were obtained. Then, α∗α∗ was selected according to the principle of minimum MRE (Figure 3).Figure 3
Weight coefficient and MRE of predicted value.According to the iterative results in Figure3, α∗ is not simply defined as 0.5 but irregularly distributed on both sides of the centerline. For the vertical deformation prediction model, α∗=0.40, and the corresponding MRE = 0.02190. For the transverse deformation prediction model, α∗=0.45, and the corresponding MRE = 0.02841. For the seepage pressure prediction model, α∗=0.77α∗=0.77, and the corresponding MRE = 0.04978.In addition, for either the GM-WC or GM-LT, the prediction accuracy could be improved by smoothing the original sequence. The predicted value (extrapolated value) in stage 9 was solved by the first 3 models in Table1. Tables 3–5 list the prediction results on dam deformation at monitoring point 01, using different models.Table 3
Predicted vertical deformation of monitoring point 01 using different models.
StageMeasured value/mmGM (1, 1)GM-LTGM-WCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%17.47.4000.007.4000.007.4000.00211.211.551−3.1311.520−2.8611.487−2.56311.912.186−2.4112.139−2.0112.116−1.81413.212.8572.6012.8052.9912.7793.19514.013.5653.1113.5243.4013.4793.72614.614.3111.9814.2992.0614.2162.63715.015.098−0.6615.138−0.9214.9950.04815.615.929−2.1116.044−2.8515.816−1.38916.516.806−1.8517.027−3.1916.681−1.10Table 4
Predicted transverse deformation of monitoring point 01 using different models.
StageMeasured value/mmGM (1,1)GM-LTGM-WCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%14.34.3000.004.3000.004.3000.0029.910.669−7.7710.607−7.1410.609−7.16311.811.865−0.5511.7200.6811.7920.07413.313.1940.7913.0042.2213.1061.46515.514.6735.3414.4936.4914.5666.02616.816.3172.8716.2273.4116.1903.63717.918.146−1.3718.255−1.9817.994−0.53819.820.179−1.9120.638−4.2320.000−1.01921.622.440−3.8923.455−8.5922.229−2.91Table 5
Predicted seepage pressure of monitoring point 01 using different models.
StageMeasured value/mmGM (1, 1)GM-LTGM-WCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%153.5053.5000.0053.5000.0053.5000.00278.9989.421−13.2189.569−13.3988.729−12.33394.8086.5788.6786.5648.6985.9059.38489.7083.8256.5583.6816.7183.1717.28582.0581.1591.0980.9151.3880.5241.86677.9778.578−0.7878.261−0.3777.9620.01775.4276.080−0.8775.713−0.3975.481−0.08870.3373.660−4.7473.266−4.1773.079−3.91964.1771.318−11.1470.916−10.5170.753−10.26Similarly, the dam deformation of monitoring points 02–06 in stage 9 was predicted by different models and compared with the above results (Table6).Table 6
Prediction results (extrapolated value) of monitoring points 02–06 in stage 9.
Monitoring contentMonitoring pointGM (1,1)GM-LTGM-WCRelative error/%Relative error/%Relative error/%Vertical deformation/mmSZY-01−1.85−3.19−1.10SZY-02−1.97−3.39−1.43SZY-03−1.67−3.11−0.84SZY-04−1.71−3.07−1.04SZY-05−1.58−2.97−1.27SZY-06−1.49−2.94−1.18Transverse deformation/mmCXY-01−3.89−8.59−2.91CXY-02−1.51−5.72−1.03CXY-03−2.95−7.53−2.17CXY-04−1.02−5.00−0.31CXY-05−2.99−7.90−1.83CXY-06−0.57−4.860.14Seepage pressure/kPaSYJ-01−11.14−10.51−10.26SYJ-02−10.33−10.04−9.56SYJ-03−11.13−10.50−10.25SYJ-04−10.01−9.71−9.25SYJ-05−10.92−10.31−10.04SYJ-06−10.81−10.54−10.04The following can be derived by comparing Tables3–6.(1)
Among the three models, the GM-WC achieved higher overall prediction accuracy than GM(1, 1) and GM-LT(2)
When the vertical deformation grew gently, all three models had a high prediction accuracy. Their relative errors were controlled within 0.84%∼3.39%. When the transverse deformation fluctuated significantly, the GM-WC witnessed the greatest improvement of prediction accuracy, with a minimum relative error of only 0.14%. In contrast, the GM-LT had a relatively low prediction accuracy, and the prediction effect of the model was not ideal after function recovery. When the original sequence of seepage pressure had peaks, the three models failed to make sufficiently accurate predictions, and their prediction accuracy needs to be further optimized.The predicted values of the three models are compared in Figures4–6 to better illustrate their prediction performance for the trend of monitoring data.Figure 4
Vertical deformation of monitoring point 01 predicted by different models.Figure 5
Transverse deformation of monitoring point 01 predicted by different models.Figure 6
Seepage pressure of monitoring point 01 predicted by different models.The following can be derived by comparing Figures4–6.The three models excelled in fitting the overall trend of vertical deformation and transverse deformation but did not perform so well in predicting the trend of seepage pressure. The predicted values only increased in stages 1–2, peaking in stage 2. This trend is different from the measured trend. In the subsequent stages 2–9, the predicted values declined slower than the measured values. Hence, the goodness of fit of the three models is positively correlated with the smoothness of the measured values. The fitting effects of these models need to be optimized, if there are peaks in the measured values.


## 4.2. Deformation Prediction of GDPM-DC

The above analysis shows that the GM-WC has certain advantages in overall prediction accuracy. But the model faces some defects in predicting volatile data. Therefore, the GDPM-DC was extended from the GM-WC to solve these defects.Firstly, the vertical deformation dataX0=7.4,11.2,11.9,13.2,14.0,14.6,15.0,15.6 of monitoring point 01 were divided into 4 subsequences with different number of samples:(26)X10=7.4,11.2,11.9,13.2,14.0,14.6,15.0,15.6,X20=11.2,11.9,13.2,14.0,14.6,15.0,15.6,X30=11.9,13.2,14.0,14.6,15.0,15.6,X40=13.2,14.0,14.6,15.0,15.6.Next, a GM-WC was established for each subsequence. The GM-WC was used to predict the values of the first 8 stages and obtain the extrapolated value of stage 9 of these subsequences. After that, the mean of the items corresponding to the predicted value of each subsequence can be calculated by(27)X^0=7.400,11.344,12.118,13.063,13.745,14.366,15.016,15.696.The mean extrapolated valueX^09 can be calculated by(28)X^09=16.408.The prediction sequenceX^0, which is calculated in a piecewise manner, was used to establish the GDPM-DC. Firstly, the relative error sequence ε=0,−0.013,−0.018,0.010,0.018,0.016,−0.001,−0.006 was established, after calculating the difference between the mean predicted value and the measured value. By the mean and MSE of the relative error sequence, the sequence was divided into four regions: S1=−0.0180,−0.0058, S2=−0.0058,−0.0008, S3=−0.0008,0.0073, and S4=0.0073,0.0180. The corresponding states are shown in Table 7.Table 7
Residual sequence states of vertical deformation at monitoring point 01.
StageRelative error/%State10.00S32−1.28S13−1.83S141.04S451.82S461.60S47−0.11S28−0.62S1Based on Table7, the transition probability of each state was entered into the state transition matrix P:(29)P1=1/2001/21000100001/302/3,P2=0001000010001/31/301/3,P3=0001000000011/21/200,P4=01/201/2000000011000.Then, the autocorrelation coefficientrK of the relative error sequence was calculated, and the weight wK of each transfer step was obtained after normalization. The calculation results are shown in Table 8.Table 8
Autocorrelation coefficient of each order and weight of each transfer step.
StepKrKwK10.4170.27220.3560.23230.5500.35940.2090.137Based on the relative error and state transition probability of vertical deformation in the first 8 stages, the deformation of the stage 9 was predicted. Besides, the maximumP1=0.45 of the relative error transition probability in stage 9 was derived from the results in Table 8 and the probability matrix in steps 1–4. Table 9 shows the predicted vertical deformation of monitoring point 01 in stage 9.Table 9
Predicted vertical deformation of monitoring point 01 in stage 9.
StageStepKInitial statewKProbability distribution of statesS1S2S3S481S10.2721/2001/272S20.232000063S40.3591/21/20154S40.1371000It can be inferred from Table9 that the relative error in stage 9 peaked in state S1. Then, the extrapolated value of vertical deformation in stage 9 can be calculated:(30)Y^09=16.408×1−−0.018−0.00582=16.603.Similarly, the lateral deformation and seepage pressure of monitoring point 01 were predicted using the above GDPM-DC. The results are summarized in Table10.Table 10
Predicted deformation of monitoring point 01 by the GDPM-DC.
Monitoring contentGDPM-DCStage123456789Vertical deformation/mmMeasured value/mm7.411.211.913.21414.61515.616.5Predicted value/mm7.40011.34412.11813.06313.74514.36615.01615.69616.603Relative error/%0.00−1.28−1.831.041.821.60−0.11−0.62−0.63Transverse deformation/mmMeasured value/mm4.39.911.813.315.516.817.919.821.6Predicted value/mm4.30010.25511.91213.38114.97016.44218.06019.83921.699Relative error/%0.00−3.58−0.95−0.613.422.13−0.89−0.20−0.46Seepage pressure/kPaMeasured value/mm53.578.9994.889.782.0577.9775.4270.3364.17Predicted value/mm53.50084.20692.06887.87482.74578.71674.89271.26366.382Relative error/%0.00−6.602.882.04−0.85−0.960.70−1.33−3.45Figures7–9 compare the predicted values in Table 9 with the results of GM (1,1) and GM-LT.Figure 7
Vertical deformation of monitoring point 01 predicted by different models.Figure 8
Transverse deformation of monitoring point 01 predicted by different models.Figure 9
Seepage pressure of monitoring point 01 predicted by different models.The following can be derived by comparing Figures7–9.The extrapolated value of stage 9 calculated by the GDPM-DC was very close to the measured value. It is more accurate than that calculated by GM (1, 1) and GM-LT. The three models had a smaller difference in predicting the trend of vertical and transverse deformations. However, the GDPM-DC boasted an obvious advantage in predicting the trend of seepage pressure. The prediction results increased in stages 1–3 and declined in stages 3–8. This trend is consistent with the measured values.To better reflect the optimization degree of the GDPM-DC, the prediction results of monitoring points 01–06 are summarized in Table11. The difference of relative error was adopted to measure the difference between the mean value of relative error of GM (1, 1), GM-LT, and GM-WC and the relative error of the GDPM-DC.Table 11
Comparison and verification of stage 9 prediction results (extrapolated value).
Monitoring contentMonitoring pointGM (1, 1)GM-LTGM-WCGDPM-DCDifference ofrelative error/%Relative error/%Relative error/%Relative error/%Relative error/%Vertical deformation/mmSZY-01−1.85−3.19−1.1−0.631.42SZY-02−1.97−3.39−1.43−0.981.28SZY-03−1.67−3.11−0.84−0.291.58SZY-04−1.71−3.07−1.04−0.581.36SZY-05−1.58−2.97−1.27−1.000.94SZY-06−1.49−2.94−1.18−0.651.22Transverse deformation/mmCXY-01−3.89−8.59−2.91−0.464.67CXY-02−1.51−5.72−1.03−0.732.02CXY-03−2.95−7.53−2.17−0.593.63CXY-04−1.02−5.00−0.31−0.351.76CXY-05−2.99−7.90−1.83−0.353.89CXY-06−0.57−4.860.14−0.071.79Seepage pressure/kPaSYJ-01−11.14−10.51−10.26−3.457.19SYJ-02−10.33−10.04−9.56−3.256.73SYJ-03−11.13−10.50−10.25−3.257.38SYJ-04−10.01−9.71−9.25−2.966.70SYJ-05−10.92−10.31−10.04−2.797.63SYJ-06−10.81−10.54−10.04−3.716.75By relative error, the different prediction models could be ranked as GDPM-DC < GM-WC < GM (1, 1) and GM-LT. The GDPM-DC is clearly superior to the other three models. The superiority is particularly notable in predicting seepage pressure: our model reduced the maximum relative error from 11.14% to 3.45%, the improvement of accuracy is the most obvious, and the improvement effect is the best.


## 4.3. Model Validation

In order to fully verify the superiority and applicability of the GDPM-DC, based on the monitoring data of transverse deformation in the first 8 stages of the dam in reference [22] and the monitoring data of accumulative vertical deformation in the first 8 stages of the railway in reference [23], the GDPM-DC is used to predict the deformation in stages 9 and 10.Firstly, the transverse deformation dataX0=6.36,6.94,7.31,7.80,8.03,8.17,8.46,8.64 in the first 8 stages of the dam were divided into 4 subsequences with different number of samples, and the mean extrapolated value X^09 can be calculated by(31)X^09=8.90.According to the division of residual sequence states of transverse deformation, the transition probability of each state was entered into the state transition matrixP:(32)P1=1/2001/2000110001/31/301/3,P2=00011000100001/201/2,P3=01/201/2000000011/2001/2,P4=01/201/2000000011000.Then, through the autocorrelation calculation and the probability matrix in steps 1–4, the maximumP1=0.45 of the relative error transition probability in stage 9 can be calculated (see Table 12 for the specific process).Table 12
Predicted transverse deformation of the dam in stage 9.
StageStepKInitial statewKProbability distribution of statesS1S2S3S481S10.1531/2001/272S40.49001/201/263S20.242000054S40.1151000It can be inferred from Table12 that the relative error in stage 9 peaked in state S4. Then, the extrapolated value of transverse deformation in stage 9 can be calculated:(33)Y^09=8.90×1−0.0029+0.0092=8.80.Similarly, the extrapolated value of transverse deformation in stage 10 can be calculated:(34)Y^010=9.10.In addition, by predicting the accumulative vertical deformation of the railway in stages 9–10, the extrapolated values can be calculated asY^09=3.04 and Y^010=3.14. The above calculation results are summarized in Tables 13 and 14.Table 13
Predicted transverse deformation of the dam using different models.
StageMeasured value/mmReference [22] modelGM-WCGDPM-DCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%16.366.55−2.936.360.006.360.0026.946.821.707.11−2.457.02−1.2837.317.112.817.36−0.707.37−1.8347.807.405.157.622.297.731.0458.037.713.987.891.737.961.8268.178.031.768.170.008.191.6078.468.361.178.460.028.43−0.1188.648.71−0.828.76−1.368.68−0.6298.709.07−4.229.07−4.228.80−1.15108.919.44−5.999.39−5.369.10−2.13Mean relative error/%—5.11—4.79—1.64Table 14
Predicted accumulative vertical deformation of the railway using different models.
StageMeasured value/mmReference [23] modelGM-WCGDPM-DCPredicted value/mmRelative error/%Predicted value/mmRelative error/%Predicted value/mmRelative error/%13.043.077−1.223.0400.003.0400.0023.213.0823.993.0166.053.113−1.2832.943.085−4.933.034−3.212.958−1.8342.853.088−8.353.053−7.122.9481.0453.173.0902.523.0723.113.0531.8263.033.092−2.053.090−1.993.0861.6073.283.0935.703.1095.203.118−0.1183.093.094−0.133.128−1.243.152−0.6292.953.095−4.923.148−6.693.040−3.05103.213.0953.583.1671.353.1402.18Mean relative error/%—4.25—4.02—2.62According to the prediction results of the transverse deformation of the dam in stages 9–10, the prediction accuracy difference between reference [22] model and the GM-WC is small, and the mean relative error is about 5%. The GDPM-DC has the highest prediction accuracy in which the relative errors of stages 9 and 10 are −1.15% and −2.13%, respectively, and the mean relative error is 1.64%. Compared with reference [22] model, its value decreases by 3.47%. Moreover, it can be seen from the prediction results of accumulative vertical deformation of the railway in stages 9–10 that the prediction accuracy of the GDPM-DC is still the highest. Compared with reference [23] model, the mean relative error is reduced by 1.63%, and the comparison results of relative errors are plotted as shown in Figures 10 and 11.Figure 10
Comparison of relative error of transverse deformation.Figure 11
Comparison of relative error of accumulative vertical deformation.Through comparison with reference [22, 23] model, it is found that the prediction value calculated by the GDPM-DC is more consistent with the monitoring data, and the prediction accuracy of the model is greatly improved, which verifies that the model has obvious superiority in predicting dam deformation. Moreover, the GDPM-DC can still have high prediction accuracy when predicting the data of accumulative vertical deformation in the railway with great fluctuation, so that the relative error of each stage can maintain a more gentle change trend, which indicates that the GDPM-DC also has a certain application value for processing more complex data of accumulative vertical deformation.


## 5. Conclusions

Under the joint influence of multiple factors, the deformation of the tailings dam is highly uncertain and nonlinear. The dam deformation cannot be predicted accurately by a single grey system model. To overcome this problem, we combined the weight coefficient with the autocorrelation coefficient to establish the GDPM-DC. Through repeated comparisons and verifications, the following conclusions were drawn:(1)
Compared with the other three models, the GDPM-DC can effectively mine the internal laws of the original data and make accurate predictions. The model controlled the relative error below 1.00%, during the prediction of the vertical and transverse deformations of the tailings dam, and controlled the relative error in 2.79%∼3.71%, during the prediction of the seepage pressure. The relative error was 6.70%∼7.63% lower than the mean value of relative error of the three contrastive models.(2)
The nonaftereffect property is another advantage of the GDPM-DC. This property mitigates the impact from the large fluctuations of data and facilitates the prediction of dam deformation. Our model accurately forecasted the increase in stages 1–3 and the decline in stages 3–9. The prediction accuracy was far better than that of any single prediction model.(3)
Through the calculation and analysis, it can be shown that the GDPM-DC can quantitatively predict the dam deformation of the operating tailings dam, which fully verifies the superiority of the model and lays the evaluation basis for seepage stability and safety management of the tailings dam. Moreover, the GDPM-DC also has high prediction accuracy when dealing with more complex accumulative settlement data, which shows that the model is also suitable for deformation prediction in other engineering fields and widens its application range, and fully verifies the applicability of the model.(4)
The GDPM-DC needs less sample data, the model is easy to program and operate, and the short-term prediction accuracy is high. However, in the subsequent model optimization, it needs to perfect the theory of optimal parameter selection to improve the prediction efficiency and reduce the requirements and restrictions of the model on data samples to expand the application range. Moreover, it is necessary to further explore the combined model under the influence of multiple factors, so that it can not only effectively predict the deformation in the short-term but also realize the long-term prediction of the deformation development trend.



---

*Source: 6103860-2022-08-16.xml*