The support vector regression (SVR) is data-based prediction model improved by the support vector classification. The basic idea of SVR is to find an optimal classification surface to minimize the error of all training samples from the optimal classification surface (Li et al., 2018). Suppose that there is a training sample set {( x i ; y i )}( i =1,2, … ; n ). x i denotes the input variable vector of the i th training sample, y i ( y i ∈ R ) represents the corresponding output value, and n is the number of sample points.

The goal of SVR is to find a function relationship f x between the input vector x i and the output vector y i under the premise that the relationship between input and output variables is unknown. The function can be expressed as a linear relationship as follows: f x = w ϕ x + b (1) where w represents the weight coefficient matrix, b represents the value of offsets, and ϕ x is a non-linear mapping function which is used to transform the complex non-linear problem in to a simple linear problem.

Based on the principle of structural risk minimization, the function f x should loosely fit training data and avoid over-fitting problem by minimizing the norm of w (i.e., w ). To cope with infeasible constrains, two slack variables ξ i and ξ i ∗ are introduced. Then the seeking process of minimizing w can be converted to solving the convex optimization problem: min 1 2 w 2 + C ∑ i = 1 m ξ i + ξ i ∗ s . t . y i − w ϕ x − b ≤ ε + ξ i − y i + w ϕ x + b ≤ ε + ξ i ∗ , i = 1,2 , ⋯ m ξ i ≥ 0 , ξ i ∗ ≥ 0 (2) where ε is the insensitive loss function which represents the error requirements for the regression function, and C C > 0 is the penalty factor. A larger C indicates that a larger penalty will be exerted on the samples when the training error is bigger than ε .

The convex optimization problem can be converted to solving the extremum of Lagrangian function L through the Lagrange multiplier method. L w , b , ξ i , ξ i ∗ , α i , α i ∗ , η i , η i ∗ = 1 2 w 2 + C ∑ i = 1 n ξ i + ξ i ∗ − ∑ i = 1 n α i ξ i + ε − y i + w ϕ x i + b − ∑ i = 1 n α i ∗ ξ i ∗ + ε + y i − w ϕ x i − b − ∑ i = 1 n η i ξ i + η i ∗ ξ i ∗ (3) where α i , α i ∗ , η i , and η i ∗ are the Lagrangian multipliers, which satisfy the positivity constraints.

According to the Karush-Kuhn-Tucher (KKT) condition which describes the necessary and sufficient conditions to meet the optimal solution, the derivatives of L about the original variable must be 0 to obtain optimal results (Smola and Scholkopf, 2004). ∂ L ∂ b = ∑ i = 1 n α i − α i ∗ = 0,0 ≤ α i , α i ∗ ≤ C (4) ∂ L ∂ w = w − ∑ i = 1 n α i − α i ∗ φ x i = 0 ⇒ w = ∑ i = 1 n α i − α i ∗ φ x i (5) ∂ L ∂ ξ I ∗ = C − α i ∗ − η i ∗ = 0 ⇒ C = α i ∗ + η i ∗ (6)

According to Eq. 5, the regression function of the SVR model can be transformed as follows. f x = ∑ i = 1 n α i − α i ∗ ϕ x i · ϕ x + b = ∑ i = 1 n α i − α i ∗ K x , x i + b (7) where K x , x i represents the kernel function.

When the SVR is used to solve the non-linear regression problem in practice, the non-linear problem is mapped to a high-dimensional space and the linear function is constructed in this space by selecting an appropriate kernel function. The selection of kernel function has a significant influence on the performance of the SVR because different kernel function is suitable for different data types (Huang et al., 2012). A radial basis kernel function is more favored than other kernel functions due to its facilitating implementation and strong mapping performance (Rasmussen, 2003). The expression of a radial basis kernel function is shown in Eq. 8. K x , x i = exp − x − x i 2 2 σ 2 (8) where σ denotes the parameter related to the width of kernel in statistics.

Using the Lagrangian multiplier method, duality principle, and the kernel function, the problem is transformed into a quadratic programming optimization one. max − ∑ i , j = 1 n α i − α i ∗ α j − α j ∗ e x p − x i − x j 2 2 σ 2 − ε ∑ i = 1 n α i + α i ∗ + ∑ i = 1 n y i α i − α i ∗ s u b j e c t   t o ∑ i = 1 n α i − α i ∗ = 0   a n d   α i , α i ∗ ∈ 0 , C (9)

Obtaining the Lagrange multiplier α i and α i ∗ from the above quadratic optimization problem, the regression function of support vector machine can be expressed as Eq. 10. f x = ∑ i = 1 n α i − α i ∗ exp − x − x i 2 2 σ 2 + b (10) where α i and α i ∗ are the Lagrange multipliers, and σ denotes the parameter related to the width of kernel in statistics.

There are two essential parameters (the penalty parameter C and kernel parameter σ ) in a SVR. The penalty parameter C controls the trade-off between the complexity of the function and the frequency in which errors are allowed. The parameter σ affects the mapping transformation of the input data to the feature space and controls the complexity of the model. Thus, it is important to select suitable parameters in the SVR.

As mentioned above, the penalty parameter C and the kernel parameter σ are essential in a SVR. Many swarm intelligence algorithms mentioned above were presented to optimize those two parameters. However, these algorithms are more likely to fall into local optimum solutions. To cope with the problem, a novel algorithm called Levy flight-based Grey Wolf Optimizer (LGWO) is introduced to select the suitable parameters for the SVR.

GWO is a swarm intelligence meta-heuristic algorithm given by Mirjalili et al. (2014). The inspiration of the GWO algorithm is based on the social hierarchy and hunting strategy of grey wolves in nature (Searemi et al., 2014). In each group of grey wolves, there is a very strict social dominant hierarchy shown in Figure 1.

To simulate the social hierarchy of grey wolves, four categories of wolves are defined--alpha ( α ), beta ( β ), delta ( δ ), and omega ( ω ). In the iterative calculation process, the first three best solutions are considered as alpha, beta, and delta, respectively. The rest of the candidate solutions are called omega. The wolves need to encircle the first three optimal solutions (alpha, beta, and delta) to find better solution for the problem (Mirjalili et al., 2014), which is modelled as: D → = C → · X → P t − X → t (11) X → t + 1 = X → P t − A → · D → (12) where t denotes the current iteration, X → p denotes the position vector of the prey, X → represents the position vector of a grey wolf, and A → , C → are the random vectors.

The random vectors A → and C → are formulated as: C → = 2 r → 2 , A → = 2 a → g r → 1 − a → (13) where a is gradually linearly decreased from 2 to 0, and r → 1 , Combination forecast model for concrete dam displacement considering residual correction r → 2 are random vectors in 0,1 .

Different mathematical operators are defined in Eqs 11–13, which can be summarized as follows: D → g B → = d 1 b 1 , d 2 b 2 , L , d n b n E → = e 1 , e 2 , L , e n where D → , B → , and E → are N-dimensional vectors. D → = d 1 , d 2 , ⋯ , d n , B → = b 1 , b 2 , ⋯ , b n , and E → = e 1 , e 2 , ⋯ , e n .

Figure 2 shows how a grey wolf updates its position X , Y according to the position of the prey X ∗ , Y ∗ . In the process of encircling prey, the grey wolf can reach different places around the best agent by adjusting the parameter values of A and C in Eqs 11, 12.

To simulate the hunting behavior of the grey wolves, the alpha, beta, and delta wolves in the GWO algorithm are three best solutions obtained so far. The omega wolves are obliged to update their positions according to the positions of the above best wolves. The hunting process can be mathematically described (Song et al., 2015): D → α = C → 1 g X → α − X → , D → β = C → 2 g X → β − X → , D → δ = C → 3 g X → δ − X → (14) X → 1 = X → α − A → 1 g D → α , X → 2 = X → β − A → 2 g D → β , X → 3 = X → δ − A → 3 g D → δ (15) X → t + 1 = X → 1 + X → 2 + X → 3 3 (16) where C 1 , C 2 , C 3 , A 1 , A 2 , and A 3 represent the random vectors, X α , X β , and X δ represent the positions of the alpha, beta, and delta wolves respectively, t represents the number of iteration, and X r is the position of current solution.

As showed in Figure 3, a grey wolf can update its position according to the positions of alpha, beta, and delta wolves in a 2D search space. The final possible position of the prey is distributed in a circle determined by the positions of alpha, beta, and delta in the search space. The hunting process can be summarized that the position of the prey is estimated by the best three wolves, whereas the other wolves update their positions randomly around the prey.

When the grey wolves start to attack the prey, the encirclement of wolves became smaller and smaller. The GWO algorithm are more likely to fall into local optimum solutions under the small search encirclement (Mirjalili et al., 2015). Considering the small encirclement, the Levy flight is introduced to increase the ability of global and local search simultaneously.

The Levy flight is a category of random search process (Amirsadri et al., 2017). In this method, the short-range exploratory hopping and occasional long-distance walking are combined to result in a more effective search. The hopping behavior ensures that the search agents can search the small areas carefully, whereas the long-distance walking behavior ensures that the search agents can enter into another areas and search a wider range. The jump size in the Levy flight follows the Levy probability distribution function (Yang and Deb, 2009). Considering the difficulty of calculating the search path, a simple mathematical definition of the Levy distribution is described: s = u / v 1 β (17) where s is the random step size obeying the Levy distribution, and u , v are the random numbers produced by normal distribution. u : N 0 , σ 2 , v : N 0,1 (18) with σ = f 1 + β g sin π β 2 β g f 0.5 + β 2 g 2 β / 2 − 0.5 1 / β (19) where f x is the standard gamma function, and the range of β is from 0 to 3.

In this study, a hybrid optimization algorithm which combines the GWO algorithm with the Levy flight is presented. In the proposed algorithm, all the wolves except the three leading wolves update their positions through the Levy flight. Therefore, the following equations can be used to update the position. X → t + 1 L G W O = X → t + 1 + S (20) where S is the step size determined by Eqs 17–19, X → L G W O t + 1 is the updated position of the wolf after the Levy flight, and X → t + 1 is the updated position of the wolf without the Levy flight calculated by Eqs 14–16.

As mentioned above, the penalty parameter C and the kernel function parameter σ in a SVR have a significant influence on the prediction performance. The LGWO algorithm is introduced to obtain the best series of parameters for the SVR. Therefore, the position of each wolf in the LGWO algorithm represents a parameter pair C , σ and the root mean squared error between the measured and predicted values is served as the fitness of each wolf.

Figure 4 shows the overall process of training the SVR using the LGWO algorithm at each iteration. At the beginning, the positions of wolves are obtained from the last iteration and served as the parameters of SVR. After giving the position of each wolf in SVR successively and training SVR, the predicted values of the testing sample are generated in SVR. Then the root mean squared error (RMSE) as the fitness of each wolf is given to the LGWO algorithm from SVR. The positions of the wolves are updated in the LGWO algorithm according to the fitness given by SVR. Finally, the best parameter pair C , σ is obtained from the LGWO algorithm and the predicted value with highest prediction accuracy is generated in SVR after reaching the maximum iteration.

In this paper, a novel hybrid model composed of the Levy flight-based grey wolf optimizer (LGWO) and support vector regression (SVR) was proposed and applied to make prediction for dam deformation. The structure of the proposed hybrid method for dam deformation prediction is shown in Figure 5. To avoid calculation error caused by numerical differences, the data is normalized and all the samples are divided into the training and testing samples. Main steps are as followed.

Step 1:

Normalize the parameters of the LGWO-SVR model.

Step 2:

Initialize the population of the wolf pack.

Step 3:

Enter the position of each wolf to the SVR model and obtain the fitness of each wolf.

Step 4:

Select the alpha, beta, and delta wolves in the wolf pack.

Step 5:

Update the positions of the omega wolves.

Step 6:

If the maximum iteration number reaches, the iteration is terminated and the position of the alpha wolf is outputted. Otherwise, repeat steps 3–6;

Step 7:

Train the SVR model according to the outputted parameters and obtain the prediction values.

To evaluate the performance of the proposed model, three widely used quantitative evaluation indicators are introduced. The specific expression of these indicators is shown as follows:

Squared correlation coefficient ( R 2 ) R 2 = n ∑ i = 1 n y ^ i y i − ∑ i = 1 n y ^ i ∑ i = 1 n y i 2 n ∑ i = 1 n y ^ i 2 − ∑ i = 1 n y i 2 2 n ∑ i = 1 n y ^ i 2 − ∑ i = 1 n y ^ i 2 (21)

Mean absolute error (MAE) M A E = 1 n ∑ i = 1 n y i − y ^ i (22)

Root mean squared error (RMSE) R M S E = 1 n ∑ i = 1 n y i − y ^ i 2 (23) where n is the number of testing samples, y ^ i is the predicted value of the i th testing sample, and y i is the measured value of the i th testing sample.

The project used in this paper is situated on the Luo River (a tributary of the Huaihe River) in Anhui province, China. It is a multiple-arch dam consisting of 20 sections and 21 arches showed in Figure 6. The total height of the dam is 75.9 m, and the total length of the dam is 510 m. To understand the real-time working status of the dam during operation, the dam is installed with pendulum monitoring system. The aim is to monitor and assess the horizontal displacements of the dam. A total of 21 pendulum systems are installed in the arches. The pendulum monitoring system consists of 20 pendulum lines (PL) and three inverted pendulum lines (IP). The distribution of the pendulums is showed in Figure 7. The raw data was recorded by manual and automated equipment every day.

At the same time, some environmental data are also monitored, such as reservoir level, air temperature, water temperature. There are 57 thermometers embedded in the dam body, which are used to measure air, water, and concrete temperatures. In this study, the dam section 13 is selected for testing the model. The thermometer distribution of the dam section 13 is showed in Figure 8.

In this study, the environmental monitoring data and displacements of the dam section 13 are used. Time series with 1700 data points from January 2009 to September 2013 are selected. The time series are divided into training and testing samples. The training time series are from January 2009 to June 2013. The testing time series are from June 2013 to September 2013. The amounts of training and testing samples are 1,000 and 260. The time series of air temperature, reservoir level, and displacements are presented in Figure 9.

Dam deformation is mainly influenced by hydraulic effect, thermal effect, and time-varying effect (Wei et al., 2019). For the hydraulic effect, it is usually considered as a reversible effect and can be represented by the polynomials( H 4 ; H 3 ; H 2 ; H ) of the reservoir level. To consider the thermal effect, the measured temperatures including air, water, dam foundation, and concrete temperatures( T a i r ; T w a t e r ; T f o u n d a t i o n ,; T c o n c t r e t e ) are used. The time-varying effect represents the irreversible deformation of dam over the time effect. The effect could be represented by the combination of θ , ln 1 + θ , and θ / θ + 1 where θ = t / 100 and t is the cumulative days from current data to initial monitoring date. Therefore, a total of 11variables are used as the input variables to construct the model and the displacements are the outputs.

As mentioned above, there are two parameters (the penalty parameter C and kernel function parameter σ ) in the SVR. The range of values of C and σ is [0.01,100] in SVR. In the LGWO algorithm, the number of the grey wolves is 20 and the maximum iteration is 100.

To describe the prediction performance of the LGWO-SVR model, four other algorithms are combined with SVR to predict the displacements: GS, PSO, GWO, and GA. Table 1 shows the initial parameters of these algorithms. To test the stability for each algorithm, 20 independent runs are carried out.