RNN is a special kind of network structure used to process sequence data. Different from the traditional neural network model, the RNN network has a high-dimensional hidden nonlinear internal structure, which can memorize the information of the previous moment and apply it to the current output vector. Theoretically, the RNN network can process sequence data of any length. But in practice, the information that the RNN network can store is limited, and there may be varying degrees of information forgetting during long-distance data transmission, and the input of the hidden layer has a significant impact on the output of the network. The influence decays as the network loop continues to recurse, making it difficult to deal with long-term dependencies such as vanishing gradients.

To solve this problem, related scholars introduced a gating mechanism upon the RNN network to form a LSTM network. Compared with the RNN network, LSTM adds three control gates, namely input gate, output gate, and forget gate. The unit structure of LSTM is shown in Figure 1, which is mainly composed of one memory unit and three gating units. Through the gated unit structure, the memory unit can effectively filter the historical information. The output vector of the hidden layer at the previous moment and the information at the current moment are combined into the input vector at the current moment, and then be stored, forgotten and output adaptively through the LSTM unit.

LSTM network can effectively deal with long-term dependence problems such as gradient disappearance in variable-length sequences, but TBM tunneling parameters are time-dependent, and the data change law is not only related to past information, but also considers the impact of future information on the model prediction performance. For changing the limitation that the traditional LSTM network only draws conclusions from past information, some scholars have improved the propagation method of the LSTM network and proposed a Bi-LSTM network.

Bi-LSTM network is composed of forward LSTM and backward LSTM. Compared with traditional LSTM network, Bi-LSTM network can better adapt to sequence data with strong round-trip correlation and effectively capture the dynamic characteristics of TBM tunneling parameters. Figure 2 shows the Bi-LSTM network structure, which adds a backpropagation layer upon LSTM. The forward training sequence and the backward training sequence are two independent LSTM structures. The two structures are symmetrical, and the information transmission directions are opposite. The forward training sequence calculates the current moment information forward, and the backward training sequence calculates the same sequence backwards. Bi-LSTM predicts all input vectors based on the timing of the input features, the hidden layer integrates the past and future information and outputs to the output layer, and finally integrates the forward and reverse information to output the prediction result. h t = f ( w 1 x t + w 2 h t − 1 ) (1) h t ′ = f ( w 3 x t + w 5 h t + 1 ′ ) (2) y t = g ( w 4 h t + w 6 h t ′ ) (3) where, h t is the current state of the forward hidden layer; h t ' is the current state of the reverse hidden layer; x t is the input value of the input layer at the current moment; y t is the output value of the output layer at the current moment; w 1 、 w 2 、 w 3 、 w 4 、 w 5 、 w 6 are the corresponding weight matrices of information propagating along the time series respectively.

PSO is a population-based computational method for intelligent bionic evolution, first proposed by American professors Eberhart and Kennedy in 1995. The concept of PSO originates from the research on the predation behavior of birds and realizes the intelligence of problem solving through cooperation and information sharing among individuals in the group. Compared with optimization algorithms such as genetic algorithm and ant colony algorithm, PSO has a relatively simple structure, fewer adjustment parameters, and has the advantage of fast convergence.

PSO uses the velocity-position iterative search method to determine the global optimal solution when solving the optimization problem. First, some particles are randomly initialized in the solution space. Each particle can be regarded as a search individual in the search space. The current velocity and position of the particle correspond to a candidate set of the optimization problem and are dynamically adjusted in accordance with the fitness function. The point where the fitness function is optimal is recorded as the current individual extreme value, also known as the local optimal solution. The local optimal solution is shared with other particles in the particle swarm, and the optimal individual extreme value is obtained by comparative analysis as the current global optimal solution of the whole particle swarm. The particle speed and position adjust the update state in the light of the two optimal solutions in the optimization process, and finally find the optimal parameters. According to the principle of the PSO above, the state update equation of the particle at time t+1 is as follows: v ( t + 1 ) = ω v t + c 1 r 1 ( P t − x t ) + c 2 r 2 ( G t − x t ) (4) x ( t + 1 ) = x t + v ( t + 1 ) (5) where, v is the particle swarm velocity; ω is the inertia weight (non-negative), the larger the value of ω , the stronger the global optimization ability and the weaker the local optimization; x is the particle swarm position; r 1 , r 2 are a random function in the range [0, 1]; c 1 , c 2 are the learning factor; P t is the local optimal solution; G t is the global optimal solution.

The basic idea of constructing the PSO-Bi-LSTM model is to use the PSO to optimize the bias and weights of the Bi-LSTM neural network, thereby improving its convergence speed and prediction accuracy. Compared with some deterministic algorithms, the PSO, as an uncertain algorithm of bionic optimization, has self-organization, robustness, and essential parallelism. Each particle cooperates with each other to better adapt to the environment, ensuring the effectiveness and practicability of the algorithm under different conditions and environments. Figure 3 is the specific algorithm flow of PSO.

The specific steps of PSO to optimize Bi-LSTM are as follows:

1) The method of constructing binary state discriminant function is used to preprocess the original data of TBM tunneling, and the data and outliers of the stop section are eliminated to obtain valid data samples.

The Songhua River water conveyance project in Jilin Province has a total length of 263.01 km and a water diversion tunnel of 133.98 km. It is a large-scale cross-regional water diversion project with the largest scale, the longest water transmission line, and the most difficult construction in Jilin Province. The total length of the TBM3 bid section is 22,955 m, and the construction projects are mainly water diversion tunnels, shaft excavation, support, and drainage work. Among them, the tunnel excavation section with a total length of 20,198 m and a maximum buried depth of 260 m is mainly constructed by TBM and supplemented by the drilling and blasting method, and its engineering geological conditions are shown in Figure 4. The lithology of the stratum involved in this section is mainly granite and limestone, and the surrounding rock grades II to V are distributed. The groundwater is divided into five grades in the light of the degree of humidity. The distribution ratio of the granite area to the limestone area is 1:1.58; the proportions of grades II to V of the surrounding rocks in the granite area are 7.03, 65.23, 25.14 and 2.60% respectively; the proportions of grades II to V of the surrounding rocks in the limestone area are 5.19, 66.70, 22.83 and 5.28% respectively.

The data of the TBM3 bit section is taken as the research object. The data collection frequency is 1 Hz, and an average of 86,400 data are collected every day. Each data includes 199-dimensional TBM tunneling parameters such as the propulsion speed, total thrust, cutter head torque and cutter head power. As shown in Figure 5, the TBM tunneling process is composed of several tunneling cycles. Generally, the cylinder stroke is used as one tunneling cycle, with about 10–20 tunneling cycles per day. The tunneling data of a complete tunneling cycle is selected to predict the tunneling parameters of the stable segment of the next tunneling cycle, thereby judging the tunneling state of TBM.

Figure 5 shows that the TBM tunneling parameter-time variation in a certain day. The adjacent tunneling cycles are divided by the shutdown state. The shutdown data (that is, the data with the TBM tunneling parameter value of 0) belongs to useless data for the machine learning algorithm. For eliminating the influence of the shutdown state data, the binary state discriminant function is constructed to preprocess the original working data and extract the effective tunneling data of each tunneling cycle from the record file. The formula for judging whether the TBM tunneling parameter value is the data of the shutdown section is as follows: I = f ( N ) f ( T ) f ( F ) f ( V ) (6) f ( x ) = { 1 , x ≠ 0 0 , x = 0 (7) where, I is the binary state discriminant function; N is the propulsion speed; F is total thrust; T is cutter head torque; V is cutter head power. When the binary state discriminant function I = 0, it means that this segment of data is the shutdown segment data; when I = 1, it means that this segment of data is the TBM normal running state data.

Due to the influence of equipment factors, construction environment and construction experience of operators, there are abnormal working conditions in the collected data samples, which affect the prediction accuracy of the model. Therefore, the 3 σ criterion is used to identify and remove outliers. The data processing flow is shown in Figure 6.

The original data is read row by row, and the input data is marked and judged based on total thrust. The first total thrust (non-zero value) data is marked as P1; then continue to read the next data, and the second total thrust (zero value) data is marked as P2. Subsequently, it is necessary to determine whether the data between P1 and P2 is within the range of 500 s to 5,000 s. If not, discard the data; if it is, output the data of P1 and P2 to the corresponding file. According to the above method, all data are sequentially extracted, and the data preprocessing process ends. Based on data preprocessing, a total of 9,900 groups of valid tunneling cycles were extracted, 9,350 groups of tunneling cycles were randomly selected as the model training set, and the remaining 550 groups of tunneling cycles were used as the model test set.

In accordance with the data preprocessing results, the working state of the TBM at each moment is determined, and valid data samples are extracted, wherein each group of continuous tunneling sequences corresponds to a complete tunneling cycle, and the data between adjacent tunneling cycles corresponds to a set of stop data. As shown in Figure 7, the TBM tunneling process can generally be divided into an ascending segment and a stable segment.

In the existing TBM tunneling parameter prediction research, several scholars such as Hou et al. (Hou et al., 2020) and Zhou et al. (Zhou et al., 2020) usually analyze and predict the parameters of the stable segment based on the data of the ascending segment, which helps the driver to judge the tunneling state and optimize the adjustment of parameters to a certain extent. However, the method of predicting TBM tunneling parameters by using the ascending segment data requires a lot of time to identify and extract data. The accuracy of the ascending segment identification directly affects the model prediction performance. Therefore, how to avoid the extraction of invalid features of the ascending segment is still an urgent problem to be solved. Secondly, the duration of the ascending segment is short, and the amount of data contained is limited, which cannot completely and accurately display the law of rock-machine interaction. In addition, the stable segment occupies most time of the total tunneling section. For the geological conditions that may be different from the ascending segment, when only the rock-machine interaction law of the ascending section is used to reflect the working state of the TBM stable segment, there will be large error. Therefore, this paper uses the data of the complete tunneling period to train the prediction model, which can effectively reduce the complexity and difficulty of the huge data processing generated during the TBM tunneling process, greatly improve the data preprocessing efficiency, and more accurately display the rock-machine interaction. It can provide a higher-quality intelligent decision-making method for TBM tunneling construction (Figure 8).

The construction conditions of TBM tunneling are often complex, and the data generated during the working process of the equipment contains 199- dimensional parameters. However, when predicting TBM tunneling parameters, the selection of characteristic variables is not the better. On the one hand, too many parameters will cause the model dimension to be too high and the calculation will be slow; on the other hand, some parameters are not highly correlated with the prediction parameters, and too many selections will reduce the model prediction accuracy. Therefore, it is necessary to perform feature selection on the 199-dimensional tunneling parameters in the TBM operating data, with the aim of assisting the optimization and adjustment of the TBM tunneling parameters.

There are many factors affecting the prediction of TBM tunneling parameters, but most of the parameters that characterize the law of rock-machine interaction are greatly influenced by human factors and have poor correlation with rock-machine interaction, so they cannot be directly used for model prediction and analysis. Therefore, this paper only selects the propulsion speed, total thrust, cutter head torque and cutter head power of the TBM tunneling stable segment as the output characteristics of the model.

The selection of input features is of great significance to the parameter optimization of the model and the efficiency of intelligent decision-making. Based on the complexity of the model, the accuracy of prediction and the comprehensive consideration of previous studies, this paper uses Pearson correlation analysis to identify the 199-dimensional TBM tunneling parameters, eliminates irrelevant parameter variables, and selects key model parameters as model inputs. The formula for calculating the Pearson correlation coefficient is as follows:

1) Calculate the average of the data:

Based on the Pearson correlation analysis results, the PSO-Bi-LSTM model constructed in this paper selects the 21 key parameters with the highest correlation with the predicted parameters in the complete tunneling cycle data as the input features of the model, to realize the prediction and analysis of the four key parameter indicators in the TBM stable segment. The 21-dimensional input feature identification results corresponding to the four prediction parameters of the model divided according to the size of the correlation coefficient of each parameter are shown in Figure 9.

To solve the problems of low prediction efficiency and slow convergence speed caused by the dimensional difference between the tunneling parameters, this paper normalizes the input data of the model, and the formula is as follows: x ′ = x − x min x max − x min (12) where, x ' is the normalized input parameter; x is the actual input value; x min is the minimum value of the input parameter; x max is the maximum value of the input parameter.

Figure 10 shows the TBM tunneling parameter prediction model based on the PSO-Bi-LSTM. Before using the PSO to optimize, the Bi-LSTM model structure needs to be determined first. Compared with the unidirectional LSTM, a back-propagation layer is added. Each layer of this structure contains two parallel layers, which control the forward and backward outputs respectively, and can utilize the contextual information in the meanwhile.

After data preprocessing, the original data is divided into training set and test set in a ratio of 17:1. The 21 key parameters with the highest correlation with a certain tunneling parameter to be predicted in the complete tunneling cycle extracted from the training set are input into the Bi-LSTM. The Bi-LSTM is optimized by PSO, and the global optimal solution obtained by iterative optimization is input into the Bi-LSTM.

To verify the generalization ability of the model, this paper adopts the goodness of fit R 2 and the mean absolute percentage error (MAPE) as the performance indicators for evaluating the intelligent prediction model. M A P E = 1 n ∑ i = 1 n | y i ∗ − y i | y i (13) R 2 = 1 − ∑ i − 1 n ( y i − y i ∗ ) 2 ∑ i = 1 n y i 2 (14) where, y i ∗ is the predicted value of the model; y i is the true value; i is the sample number; n is the total number of samples in the test set.

The preprocessed data set is divided into training set and test set proportionally, in which the data of the training set is continuously used to iteratively optimize the model and update the model parameters; the test set is used to test the performance of the model and verify the generalization ability of the model. To improve the accuracy of the training model, this paper uses the Adam optimizer to update and optimize the network, and the loss function is the MSE function. The number of hidden layer nodes and layers will be discussed in Model Parameter Analysis. A fully connected layer is added after the hidden layer, and the learned distributed feature representation is mapped to the sample label space, so that the output feature dimension is represented as one dimension, which speeds up the model calculation efficiency. Considering the prevention of the model from overfitting, each layer of unit is regularized, and the dropout probability is set to 0.3. The maximum number of iterations is set to 50. When the number of iterations reaches the requirement or the average loss value no longer decreases with the iterations, the training is stopped.

Before using PSO to optimize Bi-LSTM, the network structure of Bi-LSTM needs to be determined. Parameters such as the number of hidden layer nodes, the number of hidden layers, and the training learning rate α of the Bi-LSTM network need to be adjusted many times, so it is difficult to effectively judge the optimal state of the model. To improve the prediction accuracy, this paper takes the propulsion speed as an example to discuss the hyperparameters and determines the optimal parameters of the model through two indicators, namely the MAPE and the goodness of fit R 2 .

1) Discussion on number of hidden layer nodes in Bi-LSTM network

In the Bi-LSTM neural network, the selection of the number of hidden layer nodes is extremely important, which not only has a great impact on the performance of the neural network model, but also is the direct cause of overfitting. However, the theory of determining the optimal number of hidden layer nodes in neural networks is not yet mature. Research shows that the number of hidden layer nodes is not only related to the number of input/output layer nodes, but also depends on the complexity of the problem to be solved and the characteristics of the sample. If the number of nodes is too large, it is prone to over-fitting, which reduces the generalization ability of the model; if the number of nodes is too small, the model cannot fully learn the characteristics of time series data and reduces the performance of the model. To solve the above problems, scholars have proposed a variety of discrimination methods to determine the number of hidden layer units based on many experimental studies (Liu et al., 2012; Han et al., 2020; Zheng et al., 2020). In this study, the empirical formula is used to determine the value range of the number of hidden layer units: N h i d = ( N i n + N o u t ) 1 2 + a (15) where, N i n is the number of nodes in the input layer; N h i d is the number of hidden layer nodes; N o u t is the number of output layer nodes; a is a constant, 1 ≤ a ≤ 10 .

To prevent the model from overfitting and reduce the complexity of the model, the maximum number of iterations of the model in this section is 50, and the number of particle swarm populations is 30; other parameters such as the number of hidden layers is set to two according to experience, the learning rate is set to 0.001, and the value range of the number of the hidden layer is 3–13. Figure 11 shows the prediction of the propulsion speed with different numbers of hidden layer units based on lots of training samples. As the number of hidden layer units increases, the mean absolute percent error first decreases and then increases, and the goodness of fit first increases and then decreases. Too many or too few hidden layer nodes will reduce the prediction performance of the model. To sum up, when the number of nodes is 5, the prediction effect is optimal, so the number of hidden layer nodes in this model is set to 5.

2) Discussion on hidden layers of Bi-LSTM network

For the sake of exploring the optimal network structure, the number of layers of the model is further analyzed and studied. Theoretically, the prediction performance of the model is proportional to the number of hidden layers, but the time cost of model training increases gradually with the increase of the number of layers, the convergence effect and efficiency drop sharply, and the model often falls into the dilemma of local optimal solution, which leads to training long-term dependency problems such as information loss and gradient disappearance in the process.

Taking the propulsion speed as an example, the number of hidden layer nodes is set to 5, the maximum number of iterations is 50, the number of particle swarm populations is 30, the range of layers is set to 1∼8 based on experience, and the learning rate is 0.001.

As shown in Figure 12, when the number of hidden layers of Bi-LSTM is 2, the model accuracy is higher than that of single layer. However, the increase of the number of layers did not improve the prediction performance of the model. On the contrary, it fell into the problem of unable to find the global optimal solution, which greatly increased the risk of overfitting. Therefore, the number of hidden layers in this paper is selected as 2.

3) Discussion on learning rate α for Bi-LSTM network training

After the Bi-LSTM neural network structure is determined, the training learning rate of the model needs to be discussed. As an important hyperparameter in deep learning, the learning rate determines whether the objective function can converge to the local minimum and the convergence efficiency. If the learning rate is set too small, the convergence process is very slow, and even the algorithm does not converge; if the learning rate is set too large, the algorithm is prone to oscillation or divergence, resulting in gradient explosion (Shen et al., 2020). Therefore, the maximum number of iterations of the established model is 50, the number of particle swarm populations is 30, the number of hidden layer nodes is set to 5, the number of layers is set to 2, and the learning rate is set to 0.0001, 0.0005, 0.001, 0.01, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5 according to experience, respectively.

Taking the propulsion speed as an example, Figure 13 shows the prediction results of the model under different learning rates. When the learning rate is 0.001, the prediction result of the model is the most stable, the MAPE is relatively low compared to other cases, and the goodness of fit is relatively high, so this paper selects the model learning rate of 0.001.

In summary, the final structural parameters of the model are determined through experimental research: the number of hidden layers is 2, the number of nodes is 5, and the learning rate is 0.001.

In the actual tunneling process, the interaction between the TBM and the surrounding rock produces a complex rock-machine relationship. How to realize the real-time perception of the current rock mass information and ensure that the current TBM tunneling state matches the complex and changeable rock mass conditions is still a hot issue in TBM construction. This section studies the prediction performance of PSO-Bi-LSTM model under different surrounding rock grades based on the construction database of TBM3 bid section of Songhua River water conveyance project in Jilin Province. The 21-dimensional parameters with the highest correlation with the predicted parameters in the complete tunneling cycle are selected as the input features of the model to realize the prediction and analysis of the four key parameters (total thrust, the propulsion speed, cutter head torque and cutter head power) in the stable segment of the TBM tunneling; and the goodness of fit R 2 and the MAPE are used as evaluation indicators to judge the validity of the model prediction.

The data samples in this paper are collected from a certain construction section of the TBM3 bid section. Based on the data preprocessing, a total of 9,900 valid tunneling cycles are extracted, and 9,350 tunneling cycles are randomly selected as the model training set, and the remaining data are used as the model test set. It is used to verify the applicability of the trained model to different surrounding rock levels. Among them, the surrounding rock grades are mainly distributed in grades II to V, among which grade III surrounding rock accounts for the largest proportion in this tunneling section, and there are relatively few grades II and V surrounding rock.

Figure 14 shows the prediction results of the PSO-Bi-LSTM model under different surrounding rock levels. It can be seen from the analysis of the figure that the predicted result and the measured value of each tunneling parameter have a high degree of fitting, and there is a large deviation only in individual samples. In addition, the prediction agreement between the propulsion speed and total thrust is the highest under the conditions of different surrounding rock grades; the prediction agreement between cutter head power and cutter head torque is relatively low.

Table 1 shows the statistical analysis results of the evaluation indicators of the four key tunneling parameters under different surrounding rock grades. From the analysis of the goodness of fit and error of the prediction results, the PSO-Bi-LSTM model has high prediction accuracy. The higher the goodness of fit R 2 of the TBM tunneling parameters, the lower the MAPE, and the better the prediction effect of the model. Comparing the goodness of fit, the average values of the goodness of fit of the cutter head power, cutter head torque, the propulsion speed and total thrust are 0.912, 0.933, 0.975 and 0.972, respectively, and the goodness of fit of all reached above 0.91. This shows that the PSO-Bi-LSTM model has excellent prediction effect for each parameter. The prediction accuracy of total thrust in surrounding rocks at all levels is relatively high, and the fluctuation of the goodness of fit of its samples is much lower than that of other parameters, and the goodness of fit of other parameters fluctuates greatly under different surrounding rock grades.

In consideration of the influence of different surrounding rock grades on the prediction effect of the model, it can be found that the average values of R 2 of the tunneling parameters of grade II surrounding rock, grade III surrounding rock, and grade IV surrounding rock are 0.971, 0.962, 0.943, and 0.916, respectively. This shows that the PSO-Bi-LSTM model has a good fitting effect and can sensitively reflect the interaction information between the surrounding rock at all levels and the TBM. By comparing the goodness of fit under different surrounding rock grades, it is found that the higher the surrounding rock grade, the better the stability. Similarly, among the four types of surrounding rock grades, the MAPE of each tunneling parameter of grade V surrounding rock is also relatively high, reaching 7.685%. The average MAPE of the other three types of surrounding rocks is not much different, and they are 2.850, 3.139 and 5.344% respectively according to the grades from high to low. This shows that the stability of low-grade surrounding rock is relatively poor, and the rock-mechanism relationship is greatly affected by human factors, so the prediction accuracy is lower than that of high-grade surrounding rock. In addition, the fitting fluctuations of the parameters of the low-grade surrounding rock are relatively large, and it is difficult for the model to learn its variation law. Meanwhile, compared with high-grade surrounding rock, the problem of the number of samples of low-grade surrounding rock may cause larger accidental errors, which is also one of the reasons for the decline of model accuracy.

The same set of data was selected to quantitatively analyze the prediction performance of the Bi-LSTM model under different surrounding rock levels. Since total thrust is less affected by external factors, and the propulsion speed and cutter head torque are greatly affected by external rock mass information factors and have high sensitivity to the rock-machine relationship, this section takes cutter head torque as the research object and predict and analyze different surrounding rock grades in turn, and the prediction result is shown in Figure 15.

Figure 15 shows that the fluctuation of the prediction result curve of the Bi-LSTM model for cutter head torque under all levels of surrounding rock is significantly larger than that of the PSO-Bi-LSTM model. Combining with Table 2, among the grade II, III, IV and V surrounding rocks, the MAPE of the Bi-LSTM model is 1.938, 3.765, 1.890, and 3.408%, respectively, which is much larger than that of the PSO-Bi-LSTM model. On the one hand, it shows that the prediction accuracy of the traditional Bi-LSTM model is relatively poor; on the other hand, the Bi-LSTM model optimized by PSO can almost accurately predict the change law of tunneling parameters under different surrounding rock levels. This shows that the PSO-Bi-LSTM model has higher prediction accuracy and faster convergence speed for time sample series and can more effectively predict TBM tunneling parameters.

Given all of that, compared with the Bi-LSTM, the prediction model based on PSO-Bi-LSTM has better prediction performance and generalization ability. It can assist the safe and efficient tunneling of TBM under different stratum conditions and provide certain guidance for the optimization and adjustment of parameters.