# A Novel Intelligent Method for Predicting the Penetration Rate of the Tunnel Boring Machine in Rocks

Authors: Yan Zhang; Mingdong Wei; Guoshao Su; Yao Li; Jianbin Zeng; Xueqin Deng
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

In the construction of rock tunnels, the penetration rate of the tunnel boring machine (TBM) is influenced by many factors (e.g., geomechanical parameters), some of which are highly uncertain. It is difficult to establish a precise model for predicting the penetration rate on the basis of the influencing factors. Thus, this work proposed a useful method, based on the relevance vector machine (RVM) and particle swarm optimization (PSO), for the prediction of the TBM penetration rate. In this method, the RVM played a vital role in establishing a nonlinear mapping relationship between the penetration rate and its influencing factors through training-related samples. Then, the penetration rate could be predicted using some collected data of the influencing factors. As for the PSO, it helped to find the optimum value of a key parameter (called the basis function width) that was needed in the RVM model. Subsequently, the validity of the proposed RVM-PSO method was checked with the data monitored from a rock tunnel. The results showed that the RVM-PSO method could estimate the penetration rate of the TBM, and it proved superior to the back-propagation artificial neural network, the least-squares support vector machine, and the conventional RVM methods, in terms of the prediction performance. Moreover, the proposed RVM-PSO method could be applied to identify the difference in the importance of the various factors affecting the TBM penetration rate prediction for a tunnel.


---

## Body

## 1. Introduction

As a consequence of sustained interest in the development of underground space, tunnels play a more and more important role in human production activities and life. For example, there are many tunnels under construction in hydraulic engineering (e.g., transportation tunnels, drainage tunnels, and grouting tunnels [1–3]) and in traffic engineering (road tunnels). Thus, considerable scientific efforts have been made to develop advanced methods, techniques, and theories for improving the safety and efficiency of tunnel constructions [4–11]. In the recent decades, the tunnel boring machine (TBM) has been widely applied to tunnel constructions, in general, and long tunnels in particular, due to its many superior performances over the traditional new Austrian tunnelling method (NATM) and drilling-blasting methods, such as less pollution to the surrounding environment and higher safety, efficiency, and quality. The penetration rate of the TBM affects not only the construction schedule and the amount of capital cost but also the stability of surrounding rock mass during tunnelling. Hence, it is significant to accurately predict the penetration rate for improving construction efficiencies and reducing construction risks.Researchers have never ceased to develop empirical and theoretical models for evaluating the penetration rate/performance of the TBM since it was invented. Some typical models are the QTBM model [12], the Field Penetration Index (FPI) model [13], the Specific Rock Mass Boreability Index (SRMBI) model [14], the Colorado School of Mines (CSM) model [15], and the Norwegian University of Science and Technology (NTNU) model [16]. Among them, the CSM and NTNU are the most famous models and are widely used. The CSM model is established on the basis of laboratory test results and by analyzing the rock-breaking mechanism of the TBM. This model has the advantage that it can take into account the influence of the strength and brittleness of rock on the penetration rate. However, the original CSM model has not considered the effect of joint cracks in rock masses, which can significantly affect the TBM performance. For this reason, much research attention has been paid to replacing the compressive strength of rock with that of the rock mass, in order to consider the effect of joint cracks. Even still, the modified CSM mode, however, cannot meet the accuracy requirement for many complicated rock structures. In addition, there are only a few engineering factors considered in the derivation process of the CSM model, restricting the ability of the model to solve complex engineering problems. The NTNU model not only considers the properties of the TBM but also involves the interaction between hobs and rock mass. Moreover, all empirical parameters involved in this model are obtained from practical tunnel engineering by using back analysis; thus, the NTNU model has received a certain degree of engineering application. Nevertheless, the experiments to determine the rock drillability and the friction coefficient between rock mass and hobs are not commonly conducted in the rock mechanics community; the acquisition of these parameters limits the application of the model somewhat.In the recent years, with the rapid development of artificial intelligence (AI) technology, many intelligent models, such as the artificial neural network (ANN) [17, 18], support vector machine (SVM) [19, 20], group method of data handling (GMDH) [21], and adaptive neurofuzzy inference system (ANFIS) [22], have been applied to predicting the penetration rate. For example, there are a number of intelligent prediction models based on the dataset from the Queens Water Tunnel #3 (New York City, USA) (see Table 1), which is also used in this study to assess the proposed method. These models have a typical advantage that they can establish the complicated nonlinear relationship between the penetration rate and various influencing factors with no need to consider too many indeterminate mechanical mechanisms, which simplifies the problem of estimating the penetration rate. However, they also have some deficiencies. For example, it is difficult for the ANN to determine an appropriate network topology and reasonable hyperparameter values. As for the SVM, its parameters are always obtained by manual trials, which often contain a certain degree of subjectivity and blindness. To circumvent these issues, many researchers prefer using the optimization methods/algorithms, i.e., particle swarm optimization (PSO), cuckoo search (CS), and genetic programming (GP), to search the optimal hyperparameters for the prediction models [33]. Then, some hybrid models have been developed by combining the optimization methods with the intelligent models, including the PSO-ANN [34], PSO-SVM [35], and GP-ANN [36]. Although these hybrid models can solve the problem of parameter determination, they still have limitations, such as relatively low accuracy of ANN-based models and comparatively high computational cost/time of SVM-based models. Until now, a well-recognized model that can overcome all the problems mentioned above has not emerged. For the penetration rate prediction, it is, thus, necessary to explore a high-efficiency and accurate-enough model that can not only solve the parameter determination problem but also obtain better prediction results.Table 1
Review of several TBM performance prediction models using the dataset from the Queens Water Tunnel#3.
ReferenceModelMain factorsPredicted valueEvaluation parameterYagiz [13]SAUCS, BTS, BI, DPW,α, FNFPIR2, rGholamnejad and Tayarani [18]ANNDPW, UCS, RQDPRMSE,R2Mahdevari et al. [19]SVMUCS, BTS, BI, DPW,α, SE, TF, CP, CTPRMSE,R2Yagiz et al. [23]ANNDPW, UCS, BI,αPRRMSE, VAF,rMikaeil et al. [24]FISDPW, UCS, BTS,α, PSIPRDPDOWen et al. [25]BP-ANNUCS, BTS, SWP, PSI,αPRMREYagiz and Karahan [26]PSOUCS, BTS, BI, DPW,αPRMSEGhasemi et al. [27]FISUCS, BI, DPW,αPRVAF, RMSEYagiz and Karahan [28]DE, HS-BFGS, GWODPW, UCS, BI,αPRMSEMikaeil et al. [29]EMFEUCS, PSI, DPW,αPRR2Naghadehi et al. [30]GEPUCS, BTS, BI, DPW,αPR, FPIR2, RMSE, MAPENaghadehi et al. [31]ICA, CART, GEPUCS, BTS, DPW, BI,αPRR2, RMSE, MAPENaghadehi et al. [32]GEP, ICAUCS, BTS, DPW, BI,αPRR2, RMSE, MAPENomenclature in this table: brittleness index (BI); Brazilian tensile strength (BTS); classification and regression tree (CART); cutterhead power (CP); cutterhead torque (CT); differential evolution (DE); Dubois–Prade decision operator (DPDO); distance between planes of weakness (DPW); extended multifactorial fuzzy evaluation (EMFE); field penetration index (FPI); gene expression programming (GEP); grey wolf optimizer (GWO); hybrid harmony search (HS-BFGS); imperialist competitive algorithm (ICA); individual cutter force (FN); mean absolute percentage error (MAPE); mean relative error (MRE); mean squared error (MSE); penetration rate (PR); peak slope index (PSI); coefficient of correlation (R2); root mean square error (RMSE); rock quality designation (RQD); coefficient of the cross correlation (r); statistical analysis (SA); specific energy (SE); the spacing between the weakness planes (SWP); thrust force (TF); uniaxial compressive strength of intact rock (UCS); variance accounted for (VAF); angle between the longitudinal tunnel axis and the plane of weakness (α).To address this issue, a useful method is developed in this study for predicting the TBM penetration rate by integrating the PSO optimization algorithm into the relevance vector machine (RVM) method. In the proposed method, the RVM contributes to establishing a nonlinear mapping relationship between the penetration rate and its various influencing factors, and PSO helps to optimize the parameter of the RVM. There is only one parameter (i.e., the basis function width) needed to determine during the use of the RVM, for which fewer parameters are required than the ANN and SVM. Compared with the ANN, the RVM has a higher prediction precision, and moreover, it utilizes fewer basis functions and requires lower computational cost than the SVM. Besides, the RVM has many additional advantages, such as a small number of samples required and non-“Mercer” kernels [37]. Moreover, PSO is an excellent stochastic search optimization method with the characteristics of simple computation and rapid convergence capability, and it is widely applied to various complicated optimization problems [38]. Thus, the RVM-PSO method makes full use of the merits of both PSO and the RVM and possesses a better prediction performance for the TBM penetration rate.


## 2. Methodology

The two main parts (i.e., RVM and PSO) of the RVM-PSO method developed in this study will be introduced herein.


### 2.1. The RVM Model

The RVM model is firstly proposed by Tipping [39] on the basis of the Sparse Bayesian Learning theory, which is a probabilistic model to remedy the defects of the SVM. The principle of the RVM is as follows.It is well known that many problems involving machine learning would meet the supervised learning, in which a set of learning examples (number:N) is given as xn,tnn=1N, where xn denotes an input vector and tn represents a scalar output or target. Then, the relationship between xn and tn is described as(1)tn=yxn,ω+ξn,where yxn,ω is expressed as equation (2), and it has a linear relationship with basis functions ϕxn. ω is the weight matrix. ξn is the Gaussian noise obeying the normal distribution with a mean of 0 and a variance of σ2. Hence, the expression ptnx=Ntnyxn,σ2 belongs to a Gaussian distribution determined by tn, yxn and σ2.(2)yxn,ω=∑n=1Nωnϕxn=ωTφ,where ϕxn is often identified with the Gaussian radial basis function parameterized by training vectors, as shown in the following equation:(3)ϕxn=exp−x−xn22r2,where r is the width of radial basis function, which is an important parameter affecting the accuracy of prediction results. With the assumption of independence to tn, the complete probability expression for the sample dataset can be written as(4)ptω,σ2=2πσ2−N/2exp−12σ2t−Φω2,where t=t1t2,…,tNT and ω=ω0ω1,…,ωNT. Φ=ϕx1ϕx2,…,ϕxNT is a N×N+1 design matrix. The maximum-likelihood estimations of ω and σ2 can lead to overfitting for training because many parameters of learning samples are used. Thus, some additional constraints are generally exerted to solve this problem. The Bayesian perspective is usually adopted for the RVM to constrain the parameters by defining a simpler function of ω, which is shown as follows:(5)pωα=∏i=0NNωi0,αi−1,where α is a vector with N + 1 hyperparameters. Then, there are individual distributed parameters for every weight. Thus, the complexity of the prior distribution can be greatly reduced. In order to match the final function, a noise variance σ2 is added for α. These quantities belong to the scale parameters, for which the priors are the Gamma distribution:(6)pα=∏i=0NGammaαia,b,pβ=Gammaβc,d,where β=σ2 and(7)Gammaαa,b=Γa−1baαa−1e−bα.In equation (7), Γa=∫0∞ta−1e−tdt is the Gamma function. To ensure that these parameters have no prior information, it is better to set them to small values. If they are fixed to a = b = c = d = 0, more uniform hyperparameters can be obtained, and therefore, most of the related literature adopt this parameter setting strategy. To help to understand the general idea of this introduction, a visual implementation process is shown in Figure 1.Figure 1
A schematic of the implementation process of the RVM.On the basis of the prior, the probabilistic distribution of the posterior is given as equation (8) by following Bayes’ rule.(8)pω,α,σ2t=Ptω,α,σ2Pω,α,σ2pt.By assuming that there is a new prediction samplex∗,∗, the predictive distribution can be written as(9)pt∗t=∫pt∗ω,α,σ2pω,α,σ2tdwdαdσ2,where pω,α,σ2t cannot be directly computed, so it is decomposed as follows:(10)pω,α,σ2t=pωt,α,σ2pα,σ2t.Based on the abovementioned introduction, the posterior distribution ofω is given by(11)pωt,α,σ2=ptω,σ2pωαptα,σ2=Nμ,∑,where the posterior mean μ and covariance ∑ are given in equations (12) and (13), respectively:(12)μ=σ−2∑ΦTt.(13)∑=σ−2ΦTΦ+A−1,where A is the diagonal matrix A=diagα0,α1,…,αN.When the hyperparameters are updated to learn, many ofωi (i = 1, 2, …, N) tend to be zero, so that only a few basis functions (relevance vectors) can work on the results. Thus, the goal of sparseness is achieved.By training the learning samplesxn,tnn=1N, the prediction could proceed on the basis of the posterior distribution over the weights. After maximizing the hyperparameters, the optimal parameters αMP and σMP2 can be found. The predictive distribution can be computed by equation (9). For a new prediction sample x∗,∗, the predictive probability is given by(14)pt∗t,αMP,σMP2=∫pt∗ω,σMP2pωt,αMP,σMP2dω.Both of the two integrands on the right-hand side of equation (14) are Gaussian distributions, so ∗ also belongs to a Gaussian distribution, which is computed by(15)pt∗t,αMP,σMP2=Nt∗y∗,σ∗2.Then, the meany∗ and the covariance σ∗2 are, respectively, given as(16)y∗=μTϕx∗,σ∗2=σMP2+ϕx∗T∑ϕx∗.


### 2.2. The PSO Algorithm

Particle swarm optimization is an evolutionary computation technology based on the swarm intelligence method, which was proposed by Kennedy and Eberhart [40] and first intended for simulating social behavior. PSO is similar to the genetic algorithm (GA); both of them are optimization tools based on population information. However, unlike GA, PSO does not need crossover and mutation, and it finds the best solution by tracking the optimum particles. Compared with GA, PSO has the merits of a simpler implementation process and a better intelligent background, and it is not only used in scientific research but also applied to engineering practices in the industry [41].To solve the optimization problem, PSO is initialized with a generation of random particles first (the number of particles is denoted as NP) and, then, searches for the optimum by updating generations. In PSO, the potential solutions, i.e., particles, fly through the problem space by following the current optimum particles. Each particle keeps track of its coordinates (in the problem space), which are associated with the best solution that the current-generation particles have achieved so far. This value is calledpbest. Another “best” value tracked by the particle swarm optimizer is the best value obtained so far by all particles; this value is a global best called gbest. After finding the two best values, the particles update their velocities and positions according to the following equation:(17)vid=wvid+c1r1pid−xid+c2r2pgd−xidxid=xid+vid,where vid is the velocity vector. xid is the position vector. pid represents the best position in the current stage of particle I, and the corresponding fitness is pbest. pgd denotes the global best position in the process of particle iteration on the dth dimension, and its fitness is gbest. The parameters r1 and r2 are two random values uniformly distributed in [0, 1] [42]. c1 and c2 are acceleration constants, which can guide the particle on how to follow the pbest and the gbest. Small values of c1 and c2 may cause the particle to vibrate out of the target region; conversely, the large values can make the particle fly rapidly towards the target region and even surpass the target region. c1 and c2 usually have the same value between 1.8 and 2.0. The parameter w represents the inertia weight, which controls the influence of the previous velocity on the new velocity. When w is large, the algorithm has a strong ability of global search, while a small w value can make the algorithm have a strong ability of local search. The inertia weight w can be determined by the following equation:(18)w=wmax−wmax−wmintmaxt,where t is the current iteration step, tmax is the maximum iteration step, wmax is the maximum inertia weight, and wmin is the minimum inertia weight. In general, wmax = 0.9 and wmin = 0.4 [43].A schematic diagram is given in Figure2 to show how PSO solves the global optimization problem. In Figure 2, a two-dimensional problem is assumed. The red point represents pgd. The velocity of particle i is updated by three parts, denoted as v1,v2, and v3. v1=c1r1pid−xid is the velocity of particle i caused by pid; v2=wvid is the velocity resulted from the initial velocity; and v3=c2r2pgd−xid is the velocity induced by pgd. Then, the particle i is updated to a new position. By following the same way, the location of particle i continues to be updated and gradually approaches the global optimal solution.Figure 2
A schematic of PSO for solving the global optimization problem.


### 2.3. The RVM-PSO Method

Since the RVM contributes to establishing a nonlinear mapping relationship between the output/target and input vectors and PSO helps significantly to optimize the parameter (i.e., the basis function width) of the RVM, the RVM-PSO method is developed by integrating PSO into the RVM. The implementation steps of the RVM-PSO method are as follows.First, one should analyze and process the collected data to get the input and output datasets for the RVM-PSO model. For example, the output parameter is the penetration rate in this study, and the input parameters are the influencing factors. Given that the input parameters are usually not in the same order of magnitude and this fact will affect the learning result, the data should be normalized. The principles for normalization are expressed as(19)ai=xis,bi=yi−y¯,where ai and bi are the normalized values, xi denotes an i-dimensional input vector, yi represents the output of the ith sample, y¯ is the mean of all yi values, and s is determined using the following equation:(20)s=1n−1∑i=1nxi−x¯21/2,x=1n∑i=1nxi,y=1n∑i=1nyi.Subsequently, the normalized values of input and output data are employed to establish training sets, which can be used for training and fitting to find the nonlinear relation between the input and output data.In the next step, PSO program should be started to generate the particles, which are treated as the widths of the basis function to be involved in the RVM model. The RVM model should be optimized automatically by updating the particles in PSO to reduce the error between the predicted results and the expected results until the prediction corresponding to a certain hyperparameter meets the accuracy requirement. A flow chart of the detailed implementation process of the RVM-PSO method is shown in Figure3.Figure 3
A flow chart of the RVM-PSO method used for predicting a target variable.


## 2.1. The RVM Model

The RVM model is firstly proposed by Tipping [39] on the basis of the Sparse Bayesian Learning theory, which is a probabilistic model to remedy the defects of the SVM. The principle of the RVM is as follows.It is well known that many problems involving machine learning would meet the supervised learning, in which a set of learning examples (number:N) is given as xn,tnn=1N, where xn denotes an input vector and tn represents a scalar output or target. Then, the relationship between xn and tn is described as(1)tn=yxn,ω+ξn,where yxn,ω is expressed as equation (2), and it has a linear relationship with basis functions ϕxn. ω is the weight matrix. ξn is the Gaussian noise obeying the normal distribution with a mean of 0 and a variance of σ2. Hence, the expression ptnx=Ntnyxn,σ2 belongs to a Gaussian distribution determined by tn, yxn and σ2.(2)yxn,ω=∑n=1Nωnϕxn=ωTφ,where ϕxn is often identified with the Gaussian radial basis function parameterized by training vectors, as shown in the following equation:(3)ϕxn=exp−x−xn22r2,where r is the width of radial basis function, which is an important parameter affecting the accuracy of prediction results. With the assumption of independence to tn, the complete probability expression for the sample dataset can be written as(4)ptω,σ2=2πσ2−N/2exp−12σ2t−Φω2,where t=t1t2,…,tNT and ω=ω0ω1,…,ωNT. Φ=ϕx1ϕx2,…,ϕxNT is a N×N+1 design matrix. The maximum-likelihood estimations of ω and σ2 can lead to overfitting for training because many parameters of learning samples are used. Thus, some additional constraints are generally exerted to solve this problem. The Bayesian perspective is usually adopted for the RVM to constrain the parameters by defining a simpler function of ω, which is shown as follows:(5)pωα=∏i=0NNωi0,αi−1,where α is a vector with N + 1 hyperparameters. Then, there are individual distributed parameters for every weight. Thus, the complexity of the prior distribution can be greatly reduced. In order to match the final function, a noise variance σ2 is added for α. These quantities belong to the scale parameters, for which the priors are the Gamma distribution:(6)pα=∏i=0NGammaαia,b,pβ=Gammaβc,d,where β=σ2 and(7)Gammaαa,b=Γa−1baαa−1e−bα.In equation (7), Γa=∫0∞ta−1e−tdt is the Gamma function. To ensure that these parameters have no prior information, it is better to set them to small values. If they are fixed to a = b = c = d = 0, more uniform hyperparameters can be obtained, and therefore, most of the related literature adopt this parameter setting strategy. To help to understand the general idea of this introduction, a visual implementation process is shown in Figure 1.Figure 1
A schematic of the implementation process of the RVM.On the basis of the prior, the probabilistic distribution of the posterior is given as equation (8) by following Bayes’ rule.(8)pω,α,σ2t=Ptω,α,σ2Pω,α,σ2pt.By assuming that there is a new prediction samplex∗,∗, the predictive distribution can be written as(9)pt∗t=∫pt∗ω,α,σ2pω,α,σ2tdwdαdσ2,where pω,α,σ2t cannot be directly computed, so it is decomposed as follows:(10)pω,α,σ2t=pωt,α,σ2pα,σ2t.Based on the abovementioned introduction, the posterior distribution ofω is given by(11)pωt,α,σ2=ptω,σ2pωαptα,σ2=Nμ,∑,where the posterior mean μ and covariance ∑ are given in equations (12) and (13), respectively:(12)μ=σ−2∑ΦTt.(13)∑=σ−2ΦTΦ+A−1,where A is the diagonal matrix A=diagα0,α1,…,αN.When the hyperparameters are updated to learn, many ofωi (i = 1, 2, …, N) tend to be zero, so that only a few basis functions (relevance vectors) can work on the results. Thus, the goal of sparseness is achieved.By training the learning samplesxn,tnn=1N, the prediction could proceed on the basis of the posterior distribution over the weights. After maximizing the hyperparameters, the optimal parameters αMP and σMP2 can be found. The predictive distribution can be computed by equation (9). For a new prediction sample x∗,∗, the predictive probability is given by(14)pt∗t,αMP,σMP2=∫pt∗ω,σMP2pωt,αMP,σMP2dω.Both of the two integrands on the right-hand side of equation (14) are Gaussian distributions, so ∗ also belongs to a Gaussian distribution, which is computed by(15)pt∗t,αMP,σMP2=Nt∗y∗,σ∗2.Then, the meany∗ and the covariance σ∗2 are, respectively, given as(16)y∗=μTϕx∗,σ∗2=σMP2+ϕx∗T∑ϕx∗.


## 2.2. The PSO Algorithm

Particle swarm optimization is an evolutionary computation technology based on the swarm intelligence method, which was proposed by Kennedy and Eberhart [40] and first intended for simulating social behavior. PSO is similar to the genetic algorithm (GA); both of them are optimization tools based on population information. However, unlike GA, PSO does not need crossover and mutation, and it finds the best solution by tracking the optimum particles. Compared with GA, PSO has the merits of a simpler implementation process and a better intelligent background, and it is not only used in scientific research but also applied to engineering practices in the industry [41].To solve the optimization problem, PSO is initialized with a generation of random particles first (the number of particles is denoted as NP) and, then, searches for the optimum by updating generations. In PSO, the potential solutions, i.e., particles, fly through the problem space by following the current optimum particles. Each particle keeps track of its coordinates (in the problem space), which are associated with the best solution that the current-generation particles have achieved so far. This value is calledpbest. Another “best” value tracked by the particle swarm optimizer is the best value obtained so far by all particles; this value is a global best called gbest. After finding the two best values, the particles update their velocities and positions according to the following equation:(17)vid=wvid+c1r1pid−xid+c2r2pgd−xidxid=xid+vid,where vid is the velocity vector. xid is the position vector. pid represents the best position in the current stage of particle I, and the corresponding fitness is pbest. pgd denotes the global best position in the process of particle iteration on the dth dimension, and its fitness is gbest. The parameters r1 and r2 are two random values uniformly distributed in [0, 1] [42]. c1 and c2 are acceleration constants, which can guide the particle on how to follow the pbest and the gbest. Small values of c1 and c2 may cause the particle to vibrate out of the target region; conversely, the large values can make the particle fly rapidly towards the target region and even surpass the target region. c1 and c2 usually have the same value between 1.8 and 2.0. The parameter w represents the inertia weight, which controls the influence of the previous velocity on the new velocity. When w is large, the algorithm has a strong ability of global search, while a small w value can make the algorithm have a strong ability of local search. The inertia weight w can be determined by the following equation:(18)w=wmax−wmax−wmintmaxt,where t is the current iteration step, tmax is the maximum iteration step, wmax is the maximum inertia weight, and wmin is the minimum inertia weight. In general, wmax = 0.9 and wmin = 0.4 [43].A schematic diagram is given in Figure2 to show how PSO solves the global optimization problem. In Figure 2, a two-dimensional problem is assumed. The red point represents pgd. The velocity of particle i is updated by three parts, denoted as v1,v2, and v3. v1=c1r1pid−xid is the velocity of particle i caused by pid; v2=wvid is the velocity resulted from the initial velocity; and v3=c2r2pgd−xid is the velocity induced by pgd. Then, the particle i is updated to a new position. By following the same way, the location of particle i continues to be updated and gradually approaches the global optimal solution.Figure 2
A schematic of PSO for solving the global optimization problem.


## 2.3. The RVM-PSO Method

Since the RVM contributes to establishing a nonlinear mapping relationship between the output/target and input vectors and PSO helps significantly to optimize the parameter (i.e., the basis function width) of the RVM, the RVM-PSO method is developed by integrating PSO into the RVM. The implementation steps of the RVM-PSO method are as follows.First, one should analyze and process the collected data to get the input and output datasets for the RVM-PSO model. For example, the output parameter is the penetration rate in this study, and the input parameters are the influencing factors. Given that the input parameters are usually not in the same order of magnitude and this fact will affect the learning result, the data should be normalized. The principles for normalization are expressed as(19)ai=xis,bi=yi−y¯,where ai and bi are the normalized values, xi denotes an i-dimensional input vector, yi represents the output of the ith sample, y¯ is the mean of all yi values, and s is determined using the following equation:(20)s=1n−1∑i=1nxi−x¯21/2,x=1n∑i=1nxi,y=1n∑i=1nyi.Subsequently, the normalized values of input and output data are employed to establish training sets, which can be used for training and fitting to find the nonlinear relation between the input and output data.In the next step, PSO program should be started to generate the particles, which are treated as the widths of the basis function to be involved in the RVM model. The RVM model should be optimized automatically by updating the particles in PSO to reduce the error between the predicted results and the expected results until the prediction corresponding to a certain hyperparameter meets the accuracy requirement. A flow chart of the detailed implementation process of the RVM-PSO method is shown in Figure3.Figure 3
A flow chart of the RVM-PSO method used for predicting a target variable.


## 3. Application of the Prediction Methods to an Engineering Practice

To check the validity of the developed RVM-PSO method on the TBM penetration rate prediction, the method will be applied to estimate the penetration rate data collected from the Queens Water Tunnel#3 (New York City, USA) [25, 44], and then, the predicted values will be compared with the actual values. Moreover, the prediction performance of the RVM-PSO method will be compared with that of the back-propagation artificial neural network (BP-ANN) [25], the least-squares support vector machines (LS-SVM), and the conventional RVM methods, which have already been widely applied to nonparametric and nonlinear classification and forecasting problems. For assessing our proposed RVM-PSO method by comparing the results with those from the BP-ANN method, this study adopted the same data as those used in the BP-ANN method in reference [25], that is, the data in Table 2.Table 2
Dataset of the TBM penetration rate (PR) collected from the Queens Water Tunnel#3.
Sample no.α (°)UCS (MPa)BTS (MPa)SPW (m)PSI (kN/mm)PR (m/h)L152141.49.41.6261.74L228153.49.12231.94L344147.610.82302.09L4111319.90.4321.78L528140.79.11.6382.1L633173.19.80.8311.84L77117610.32301.6L881150.710.81.6301.98L9831258.31.6272.21L1069129.28.20.8291.84L1173134.590.2262.01L1255192.711.40.4422.18L1366182.410.20.8392L14741598.70.2362.28L1551156.910.70.8302.16T14018992562.2∗T23419110.40.4542.95∗T370137.28.81.6422.2∗T432136.29.20.4402.03∗T510128.69.91.6321.75∗T6111319.90.4321.78∗T750138.6100.2302.39∗T818133.410.51.6301.6∗T952154.510.12331.58∗T1033173.19.80.8311.84∗


### 3.1. Implementation of the RVM-PSO Method

Note that the factors affecting the TBM penetration rate are too many to consider all of them in a prediction model. For this reason, it is, thus, common to consider only the relatively significant influencing factors in the applications of the current intelligent prediction methods. In addition, since the TBM machines used during the whole tunnelling process of the Queens Water Tunnel #3 are identical, some mechanical and geometrical parameters of TMB machines, such as the hob radius, the thickness of the hob end, and the rated thrust of a single hob, are not the main factors to induce the variation of penetration rate data recorded in that tunnelling practice. Hence, these parameters are not considered in this case study. Five parameters are tentatively speculated as to the main factors affecting the TBM penetration rate of the Queens Water Tunnel#3, and they are involved in the penetration rate prediction model established here. The five main factors are the angle between the longitudinal tunnel axis and the plane of weakness (α), the uniaxial compressive strength of intact rock (UCS), the Brazilian tensile strength (BTS), the spacing between the weakness planes (SWP) [45, 46], and the peak slope index (PSI, i.e., the ratio of the maximum force to the corresponding displacement in the impact test [47–54]).Table2 lists the penetration rate data corresponding to different values of the main influencing factors. Fifteen groups of the data are randomly selected as learning samples (denoted as L1 to L15), and the remaining ten groups are assigned as test samples (represented as T1 to T10). In the RVM-PSO model established here, the penetration rate and the main influencing factors are the output and input parameters, respectively. PSO is employed for the optimization of the Gaussian kernel width. The termination condition is set to f < 0.03 (where f is the mean relative error), and the maximum number of iterations is specified as 1000. Parameter settings in PSO are d = 1, NP = 10, c1 = c2 = 2, and Vmax = 1. After establishing the nonlinear relationship between the penetration rate and the influencing factors by training the learning samples, the prediction can be made by inputting some values of the key influencing factors that are measured from the tunnel. In this study, PSO finds the best value of the Gaussian kernel width through the search of three-generation particles. The detailed search process of the particles is shown in Figure 4, which indicates that the particles of a newer generation are closer to the termination condition, and the mean relative error corresponding to one of the third-generation particles is only 0.028, achieving the required accuracy.Figure 4
The detailed search process of the particles in the RVM-PSO model established for predicting the TBM penetration rate in the Queens Water Tunnel#3.


### 3.2. Brief Introduction of the BP-ANN and the LS-SVM Methods

Since the prediction made by the RVM-PSO method will be compared with those by the BP-ANN and the LS-SVM methods, the latter two methods are briefly introduced here.


#### 3.2.1. The BP-ANN Method

The artificial neural network (ANN) is an information system proposed by abstracting, simplifying, and simulating the structure, function, and basic characteristics of the human brain neural network. The back-propagation artificial neural network (BP-ANN), originally developed by McClelland and Rumelhart [55], is a multilayer forward neural network based on the error back-propagation algorithm. As a typical learning algorithm of the ANN, the main structure of the BP-ANN is composed of an input layer, one or more hidden layers, and an output layer. Each layer is composed of a number of neurons (nodes), and the output value of each node is determined by input values, action functions, and thresholds. The BP learning algorithm is used to train multilayer networks and implement a gradient search to minimize the squared error between realized and desired outputs.The learning process of the network includes two processes: forward information propagation and backward error propagation. In the forward propagation process, the input information is transmitted from the input layer to the output layer through the hidden layer. After the action functions are operated, the output value is compared with the expected value. If there is an error, the error propagates back and returns along the original connection path. By modifying the weight of neurons layer by layer, the error can be reduced, and the loop proceeds until the output result meets the specified accuracy requirement.The outputOj of neuron j on the hidden or output layer of the BP neural network is determined by the following equation:(21)Oj=fj∑ωjxi+θ,where fj is the excitation function corresponding to neuron j, θ is the threshold of neuron j, xi denotes the inputs to neuron j, and ωj represents the connection weight corresponding to the input and neuron j.


#### 3.2.2. The LS-SVM Method

In the least-squares support vector machine method proposed by SuyKens et al. [56], the least-squares linear system is used as the loss function, and the solution process is transformed into solving a set of equations, which dramatically speeds up the solution procedure. Thus, the LS-SVM has achieved good results on pattern recognition and nonlinear estimation. In the use of the LS-SVM, the samples of training data can be expressed as (X1, Y1), (X2, Y2), …, (Xl, Yl), where yi is the target value and xi is the input vector. The function estimation problem can be described as solving the following question:(22)minω,b,eJω,e=12ωTω+γ12∑i=1lei2.The constraint condition is as follows:(23)yi=Φxi×w+b+ξi,i=1,2,…,l,where Φ(): Rn⟶Rnh is the kernel function, ω (∈Rnh) is the weight vector, ξi (∈R) is the error variable, b is the deviation value, and γ an adjustable parameter.The kernel function can map the samples in the original space to a vector in a high-dimensional feature space to solve the linear indivisibility problem, and thus, the optimization problem can be solved using the Lagrangian method:(24)Lω,b,e,α=Jω,e−∑i=1lαiωTΦxk+b+ei−yi,where αi (i = 1, 2, …, l) is the Lagrange multiplier. The optimization conditions are as follows:(25)∂L∂ω=0⟶ω=∑i=1lαiΦxi,∂L∂b=0⟶∑i=1lαi=0,∂L∂ei=0⟶αi=γei,∂L∂αi=0⟶ωTΦxi+b+ei−yi=0.Subsequently, using these optimization conditions, one can get(26)01vT−1vΩ+1γIbα=0y,where x = [x1, x2, …, xl], y = [y1, y2, …, yl], and α = [α1, α2, …, αl]. The kernel function K (xk,xl) (=Φ(xk)TΦ(xl)) is a symmetric function satisfying Mercer’s condition.The regression estimate of the least-squares support vector machine is as follows:(27)yx=∑i=1lαiKx,xi+b,where the parameters α and b are determined from equation (26) and the kernel function K (x, xi) is an arbitrary symmetric function satisfying Mercer’s condition, such as K (x, y) = exp(−(x − y)2/(2σ2)).


## 3.1. Implementation of the RVM-PSO Method

Note that the factors affecting the TBM penetration rate are too many to consider all of them in a prediction model. For this reason, it is, thus, common to consider only the relatively significant influencing factors in the applications of the current intelligent prediction methods. In addition, since the TBM machines used during the whole tunnelling process of the Queens Water Tunnel #3 are identical, some mechanical and geometrical parameters of TMB machines, such as the hob radius, the thickness of the hob end, and the rated thrust of a single hob, are not the main factors to induce the variation of penetration rate data recorded in that tunnelling practice. Hence, these parameters are not considered in this case study. Five parameters are tentatively speculated as to the main factors affecting the TBM penetration rate of the Queens Water Tunnel#3, and they are involved in the penetration rate prediction model established here. The five main factors are the angle between the longitudinal tunnel axis and the plane of weakness (α), the uniaxial compressive strength of intact rock (UCS), the Brazilian tensile strength (BTS), the spacing between the weakness planes (SWP) [45, 46], and the peak slope index (PSI, i.e., the ratio of the maximum force to the corresponding displacement in the impact test [47–54]).Table2 lists the penetration rate data corresponding to different values of the main influencing factors. Fifteen groups of the data are randomly selected as learning samples (denoted as L1 to L15), and the remaining ten groups are assigned as test samples (represented as T1 to T10). In the RVM-PSO model established here, the penetration rate and the main influencing factors are the output and input parameters, respectively. PSO is employed for the optimization of the Gaussian kernel width. The termination condition is set to f < 0.03 (where f is the mean relative error), and the maximum number of iterations is specified as 1000. Parameter settings in PSO are d = 1, NP = 10, c1 = c2 = 2, and Vmax = 1. After establishing the nonlinear relationship between the penetration rate and the influencing factors by training the learning samples, the prediction can be made by inputting some values of the key influencing factors that are measured from the tunnel. In this study, PSO finds the best value of the Gaussian kernel width through the search of three-generation particles. The detailed search process of the particles is shown in Figure 4, which indicates that the particles of a newer generation are closer to the termination condition, and the mean relative error corresponding to one of the third-generation particles is only 0.028, achieving the required accuracy.Figure 4
The detailed search process of the particles in the RVM-PSO model established for predicting the TBM penetration rate in the Queens Water Tunnel#3.


## 3.2. Brief Introduction of the BP-ANN and the LS-SVM Methods

Since the prediction made by the RVM-PSO method will be compared with those by the BP-ANN and the LS-SVM methods, the latter two methods are briefly introduced here.


### 3.2.1. The BP-ANN Method

The artificial neural network (ANN) is an information system proposed by abstracting, simplifying, and simulating the structure, function, and basic characteristics of the human brain neural network. The back-propagation artificial neural network (BP-ANN), originally developed by McClelland and Rumelhart [55], is a multilayer forward neural network based on the error back-propagation algorithm. As a typical learning algorithm of the ANN, the main structure of the BP-ANN is composed of an input layer, one or more hidden layers, and an output layer. Each layer is composed of a number of neurons (nodes), and the output value of each node is determined by input values, action functions, and thresholds. The BP learning algorithm is used to train multilayer networks and implement a gradient search to minimize the squared error between realized and desired outputs.The learning process of the network includes two processes: forward information propagation and backward error propagation. In the forward propagation process, the input information is transmitted from the input layer to the output layer through the hidden layer. After the action functions are operated, the output value is compared with the expected value. If there is an error, the error propagates back and returns along the original connection path. By modifying the weight of neurons layer by layer, the error can be reduced, and the loop proceeds until the output result meets the specified accuracy requirement.The outputOj of neuron j on the hidden or output layer of the BP neural network is determined by the following equation:(21)Oj=fj∑ωjxi+θ,where fj is the excitation function corresponding to neuron j, θ is the threshold of neuron j, xi denotes the inputs to neuron j, and ωj represents the connection weight corresponding to the input and neuron j.


### 3.2.2. The LS-SVM Method

In the least-squares support vector machine method proposed by SuyKens et al. [56], the least-squares linear system is used as the loss function, and the solution process is transformed into solving a set of equations, which dramatically speeds up the solution procedure. Thus, the LS-SVM has achieved good results on pattern recognition and nonlinear estimation. In the use of the LS-SVM, the samples of training data can be expressed as (X1, Y1), (X2, Y2), …, (Xl, Yl), where yi is the target value and xi is the input vector. The function estimation problem can be described as solving the following question:(22)minω,b,eJω,e=12ωTω+γ12∑i=1lei2.The constraint condition is as follows:(23)yi=Φxi×w+b+ξi,i=1,2,…,l,where Φ(): Rn⟶Rnh is the kernel function, ω (∈Rnh) is the weight vector, ξi (∈R) is the error variable, b is the deviation value, and γ an adjustable parameter.The kernel function can map the samples in the original space to a vector in a high-dimensional feature space to solve the linear indivisibility problem, and thus, the optimization problem can be solved using the Lagrangian method:(24)Lω,b,e,α=Jω,e−∑i=1lαiωTΦxk+b+ei−yi,where αi (i = 1, 2, …, l) is the Lagrange multiplier. The optimization conditions are as follows:(25)∂L∂ω=0⟶ω=∑i=1lαiΦxi,∂L∂b=0⟶∑i=1lαi=0,∂L∂ei=0⟶αi=γei,∂L∂αi=0⟶ωTΦxi+b+ei−yi=0.Subsequently, using these optimization conditions, one can get(26)01vT−1vΩ+1γIbα=0y,where x = [x1, x2, …, xl], y = [y1, y2, …, yl], and α = [α1, α2, …, αl]. The kernel function K (xk,xl) (=Φ(xk)TΦ(xl)) is a symmetric function satisfying Mercer’s condition.The regression estimate of the least-squares support vector machine is as follows:(27)yx=∑i=1lαiKx,xi+b,where the parameters α and b are determined from equation (26) and the kernel function K (x, xi) is an arbitrary symmetric function satisfying Mercer’s condition, such as K (x, y) = exp(−(x − y)2/(2σ2)).


## 3.2.1. The BP-ANN Method

The artificial neural network (ANN) is an information system proposed by abstracting, simplifying, and simulating the structure, function, and basic characteristics of the human brain neural network. The back-propagation artificial neural network (BP-ANN), originally developed by McClelland and Rumelhart [55], is a multilayer forward neural network based on the error back-propagation algorithm. As a typical learning algorithm of the ANN, the main structure of the BP-ANN is composed of an input layer, one or more hidden layers, and an output layer. Each layer is composed of a number of neurons (nodes), and the output value of each node is determined by input values, action functions, and thresholds. The BP learning algorithm is used to train multilayer networks and implement a gradient search to minimize the squared error between realized and desired outputs.The learning process of the network includes two processes: forward information propagation and backward error propagation. In the forward propagation process, the input information is transmitted from the input layer to the output layer through the hidden layer. After the action functions are operated, the output value is compared with the expected value. If there is an error, the error propagates back and returns along the original connection path. By modifying the weight of neurons layer by layer, the error can be reduced, and the loop proceeds until the output result meets the specified accuracy requirement.The outputOj of neuron j on the hidden or output layer of the BP neural network is determined by the following equation:(21)Oj=fj∑ωjxi+θ,where fj is the excitation function corresponding to neuron j, θ is the threshold of neuron j, xi denotes the inputs to neuron j, and ωj represents the connection weight corresponding to the input and neuron j.


## 3.2.2. The LS-SVM Method

In the least-squares support vector machine method proposed by SuyKens et al. [56], the least-squares linear system is used as the loss function, and the solution process is transformed into solving a set of equations, which dramatically speeds up the solution procedure. Thus, the LS-SVM has achieved good results on pattern recognition and nonlinear estimation. In the use of the LS-SVM, the samples of training data can be expressed as (X1, Y1), (X2, Y2), …, (Xl, Yl), where yi is the target value and xi is the input vector. The function estimation problem can be described as solving the following question:(22)minω,b,eJω,e=12ωTω+γ12∑i=1lei2.The constraint condition is as follows:(23)yi=Φxi×w+b+ξi,i=1,2,…,l,where Φ(): Rn⟶Rnh is the kernel function, ω (∈Rnh) is the weight vector, ξi (∈R) is the error variable, b is the deviation value, and γ an adjustable parameter.The kernel function can map the samples in the original space to a vector in a high-dimensional feature space to solve the linear indivisibility problem, and thus, the optimization problem can be solved using the Lagrangian method:(24)Lω,b,e,α=Jω,e−∑i=1lαiωTΦxk+b+ei−yi,where αi (i = 1, 2, …, l) is the Lagrange multiplier. The optimization conditions are as follows:(25)∂L∂ω=0⟶ω=∑i=1lαiΦxi,∂L∂b=0⟶∑i=1lαi=0,∂L∂ei=0⟶αi=γei,∂L∂αi=0⟶ωTΦxi+b+ei−yi=0.Subsequently, using these optimization conditions, one can get(26)01vT−1vΩ+1γIbα=0y,where x = [x1, x2, …, xl], y = [y1, y2, …, yl], and α = [α1, α2, …, αl]. The kernel function K (xk,xl) (=Φ(xk)TΦ(xl)) is a symmetric function satisfying Mercer’s condition.The regression estimate of the least-squares support vector machine is as follows:(27)yx=∑i=1lαiKx,xi+b,where the parameters α and b are determined from equation (26) and the kernel function K (x, xi) is an arbitrary symmetric function satisfying Mercer’s condition, such as K (x, y) = exp(−(x − y)2/(2σ2)).


## 4. Comparison of the Prediction Results Obtained from the Methods

After learning the training samples, the prediction sample set was predicted. The prediction results for the ten prediction samples are summarized in Table3, in which the prediction results based on other methods were also listed for comparison. The maximum relative error of the BP-ANN method is up to 17.24%, and that of the LS-SVM and the RVM methods are as high as −14.7% and −14.27%, respectively. However, the RVM-PSO method has the minimum value of the maximum relative error, which is only 12.7%. To compare the prediction performances of the four methods (i.e., RVM-PSO, BP-ANN, LS-SVM, and RVM) in an intuitive manner, the predicted penetration rates and the actual values for all the test samples are drawn in Figure 5. It can be seen that the BP-ANN method provides terrible estimates for the penetration rates in samples T1, T2, T4, and T6. As for the LS-SVM method, the predictions on samples T7 and T9 significantly deviate from the measured values. When the conventional RVM method is used, there are also relatively significant differences between the predicted and actual penetration rates for samples T7 and T9. As with the RVM-PSO method, the prediction results seem close to the actual penetration rates, except for sample T9. These observations indicate that the capabilities of those four methods to predict the penetration rate are not the same, and moreover, if the test samples are sorted by the accuracies of their corresponding prediction results, the sequences are not consistent for all of those methods. Most importantly, it appears that the RVM-PSO method overall yields the best prediction. However, it is still difficult to sufficiently distinguish the capacities of those four methods to estimate the TBM penetration rate (for example, Figure 5 shows that the predictions by the LS-SVM and the conventional RVM methods are quite close for a given test sample). Hence, more reasonable and in-depth analysis is required to assess the four methods further. For this reason, some performance indices, including the mean absolute error (MAE), maximum absolute error (MXAE), mean relative error (MRE), root mean square error (RMSE), and Theil’s inequality coefficient (TIC) [57], are calculated here for evaluating and comparing the prediction performances of those methods. The formulas for determining these parameters are as follows:(28)MAE=1n∑i=1nyi−y^i,MXAE=maxyi−y^i,MRE=1n∑i=1nyi−y^iyi×100%,RMSE=1n∑i=1nyi−y^i2,TIC=1/n∑i=1nyi−y^i2∑i=1nyi2+∑i=1ny^i2,where n denotes the total number of samples; yi and y^i represent the measured value and the predicted value, respectively. The terms MAE, MXAE, MRE, and RMSE can reflect the prediction accuracy of the methods. The index TIC represents the level of agreement between the prediction and the real penetration rate.Table 3
Comparison of the predicted results with different methods.
Sample no.In situ measurement (m/h)BP-ANN [25]LS-SVMRVMRVM-PSOThe predicted results (m/h)The relative errors (%)The predicted results (m/h)The relative errors (%)The predicted results (m/h)The relative errors (%)The predicted results (m/h)The relative errors (%)T12.2∗2.4611.822.220.982.18−0.82.18−0.73T22.95∗3.136.12.78−5.722.84−3.732.93−0.68T32.2∗2.19−0.452.1−4.662.06−6.222.13−3.05T42.03∗2.3817.242.061.532.040.552.13.48T51.75∗1.74−0.571.781.481.71−2.461.771.37T61.78∗1.61−9.551.896.371.917.251.780.16T72.39∗2.27−5.022.04−14.72.05−14.272.36−1.31T81.6∗1.716.871.8012.651.7610.181.653.13T91.58∗1.57−0.631.7410.351.7410.251.7812.7T101.84∗1.71−7.071.955.841.998.221.861.08Figure 5
Comparison between the penetration rates measured in situ and those predicted using the BP-ANN, LS-SVM, RVM-PSO, and conventional RVM methods.The calculated values of these five parameters are delineated in Figure6. As can be seen from Figure 6(a), there are some differences between the mean absolute errors of the predictions using the BP-ANN, LS-SVM, RVM, and RVM-PSO methods. The mean absolute error of the RVM-PSO method is only about 0.05 m/h, significantly lower than those of the other three methods. In addition, the BP-ANN method yields the highest mean absolute error, and the MAE value corresponding to the RVM method is lower than that generated by the LS-SVM method. Figure 6(b) shows that the RVM-PSO and LS-SVM methods give the smallest and largest values of the maximum absolute error, respectively, and the BP-ANN method produces a slightly lower MXAE than the LS-SVM method. Figure 6(c) depicts that the mean relative error (MRE) of the prediction results of the BP-ANN method is 6.53%, and that of both the LS-SVM and the RVM methods is 6.4%, while the MRE value of the RVM-PSO method is only 2.8%, even less than half of those of the former three methods. As for the remained two indices, the RVM-PSO method still generates much lower RMSE and TIC values than the other three methods, for which both the RMSE and TIC values decline from the BP-ANN to LS-SVM to RVM. Overall, it can be found from Figure 6 that the distinctions between the prediction performances of the BP-ANN, LS-SVM, and RVM methods are not very significant. More specifically, the RVM method is slightly better than the LS-SVM method, which can modify the prediction a little compared with the BP-ANN method. Most importantly, Figure 6 indicates that no matter from which index the viewpoint is, the prediction using the RVM-PSO method is always more satisfactory than that using the other three methods. Therefore, the proposed RVM-PSO method proves valid and superior for the TBM penetration rate prediction.Figure 6
Comparison between the prediction performances of the four methods in terms of (a) MAE, (b) MXAE, (c) MRE, (d) RMSE, and (e) TIC.
(a)(b)(c)(d)(e)


## 5. Comparison of the Effects of the Main Influencing Factors

As mentioned above, it is basically impractical to obtain the data of all the factors affecting the TBM penetration rate, and thus, it is quite interesting and useful to identify the difference in the importance of the various influencing factors for the TBM penetration rate prediction of a tunnel because this work may help researchers/engineers to ascertain or pay attention to the parameters that are significant for the penetration rate of a tunnel. In the abovementioned application example of the RVM-PSO method, five parameters are considered as the main influencing factors of the TBM penetration rate. Their importance would be evaluated and compared in this section. To this end, five new prediction models without considering each of the five influencing factors, respectively (i.e., the prediction models correspond to different combinations of just four of the five influencing factors), are established on the basis of the RVM-PSO method to predict the penetration rate. The reason for selecting the proposed RVM-PSO method is that it is already manifested to be capable of the TBM penetration rate prediction. Note that similar to the settings used in the previous section, the first fifteen groups of data listed in Table2 are deemed as learning samples for the five new prediction models, and the remaining ten groups work as test samples. Moreover, the structures and parameters (e.g., the Gaussian kernel width and the number of iterations) of each model are consistent with those used in the previous model, in which all the five influencing factors are taken into account. Then, the five four-factor models can be applied to estimate the penetration rates in the test samples, and the importance of the five influencing factors on the penetration rate can be determined by comparing the prediction accuracies of the four-factor models with that of the previous five-factor model. By following the implementation steps described in the previous sections, the four-factor models established in this section can also make predictions of the penetration rates in the test samples. The predicted penetration rates are tabulated in Table 4.Table 4
The predictions yielded by the five four-factor RVM-PSO models for the penetration rates (PR) in the test samples.
Sample no.The PR (m/h) predicted using the RVM-PSO method with one missing main influencing factorαUCSBTSSWPPSIT12.092.212.102.182.12T22.892.902.912.922.88T32.052.132.072.102.06T42.072.142.011.961.98T51.801.711.771.791.80T62.021.771.851.791.87T72.042.272.151.972.30T81.691.721.741.851.80T91.901.721.842.001.83T102.021.871.951.981.95Since the difference in the importance of the influencing factors can be reflected by the different prediction performances of the four-factor models, it is useful to compare the values of the statistic indices yielded by the models. To facilitate the comparison, the following five parameters are defined:(29)T1j=MAEjMAE0,T2j=MXAEjMXAE0,T3j=MREjMRE0,T4j=RMSEjRMSE0,T5j=TICjTIC0,where T1j, T2j, …, T5j represent the ratio of the five indices (MAE, MXAE, MRE, RMSE, and TIC) yielded by the new models to those by the previous five-factor model, respectively. Values of j (=1, 2, …, 5) correspond to the five influencing factors (α, UCS, BTS, SWP, and PSI) that are not considered in the new models in turn. For example, T31 represents the ratio of the MRE value obtained from the new model without considering α to that from the previous model. If T31 > T32, this means that the mean relative error produced by the model that does not consider α, which is larger than that obtained from the one without considering UCS.To make theTij (i, j = 1, 2, …, 5) values yielded by the five new models easy to compare, they are delineated in Figure 7, and the detailed Tij results are listed in Table 5. As can be seen from Figure 7 and Table 5, the values of Tij are basically larger than 1 (except for T22). This means that the prediction results of the four-factor models are basically inferior to those of the previous five-factor model, indicating that each of the five main influencing factors indeed affects the prediction results and the five main influencing factors selected in the previous sections are reasonable. Moreover, for a given value of i, the values of Tij are usually not the same when j varies, suggesting that the importance of the five influencing factors is different for the TBM penetration rate prediction of the Queens Water Tunnel #3. It is observed that Ti1 and Ti4 are always larger than Ti2, Ti3, and Ti5, and the difference between Ti1 and Ti4 is not very conspicuous for any value of i; therefore, the angle between the longitudinal tunnel axis and the plane of weakness (i.e., α) and the spacing between the weakness planes (SWP) have significant effects.Figure 7
TheTij values induced by the five four-factor RVM-PSO models for the cases of (a) i = 1, (b) i = 2, (c) i = 3, (d) i = 4, and (e) i = 5.
(a)(b)(c)(d)(e)Table 5
TheTij values induced by the five four-factor RVM-PSO models for the predictions of the penetration rates in the test samples.
iTijj = 1j = 2j = 3j = 4j = 512.951.342.152.922.2222.561.141.892.891.7633.141.412.2732.2741.740.721.312.111.2652.571.151.902.891.76Moreover, for a fixed value ofi, Ti3 and Ti5 are basically equal; hence, the importance of the Brazilian tensile strength (BTS) and that of the peak slope index are nearly the same for this tunnel practice. In addition, Figure 7 clearly shows that the model without considering the UCS has the closest prediction performance to the previous five-factor model, and the mean relative error yielded by the former is only 3.71%. These observations imply that the UCS has a smaller influence on the TBM penetration rate prediction of the Queens Water Tunnel #3 than the other four factors, while α and SWP are quite important for the prediction. Therefore, in the penetration rate prediction of that tunnel, data of UCS might be able to leave out for the purpose of reducing workload, but the acquisition of α and SWP should be paid more attention.In practical engineering, it is quite difficult to acquire the detailed data of all the factors affecting the TBM penetration rate, and even the data of some important factors might not be collected. Hence, it is meaningful to assess the predictive ability of the proposed RVM-PSO model under the working condition that not the data of all the main influencing factors are available. For comparison, both the RVM-PSO and the conventional RVM methods are employed to estimate the penetration rates in the test samples listed in Table2, under the assumption that the values of α or SWP are not known. In other words, both the RVM-PSO and RVM models compared here belong to the four-parameter model without considering α or SWP. The reason for selecting the RVM method (instead of the BP-ANN and LS-SVM methods) to be compared with the RVM-PSO method here is that the RVM is already proved to be more accurate than the BP-ANN and LS-SVM for TBM penetration rate prediction, as shown in the previous section. Figure 8 presents the relative errors of the penetration rates predicted by the RVM-PSO and the conventional RVM methods when α or SWP is not taken into account in the predictions. It can be found that the RVM-PSO method always generates smaller errors in the prediction than the RVM for each of the test samples, no matter for the absence of α or SWP. Thus, the RVM-PSO method is superior to the RVM method. Moreover, even when one of the two most important influencing factors identified for the TBM penetration rate prediction of the Queens Water Tunnel #3 is not considered, the mean relative error yielded by the RVM-PSO method is about 8%, less than 10%. Therefore, the RVM-PSO method is useful to estimate the TBM penetration rate in rock tunnels, and it has a relatively good tolerance for missing or lacking the data of a few influencing factors.Figure 8
The relative errors (%) of the penetration rates predicted using the RVM-PSO and the conventional RVM methods without considering (a)α and (b) SWP.
(a)(b)


## 6. Conclusions

Accurate prediction of the TBM penetration rate is of great significance to tunnel constructions, for purposes of optimizing design schemes, saving economic costs, and evaluating the feasibility of engineering constructions. Thus, a useful intelligent method, RVM-PSO, is developed to estimate the penetration rate of the TBM. To check the validity of the developed method, it was applied to the prediction of the TBM penetration rate in the Queens Water Tunnel#3. Then, the prediction performance of the RVM-PSO method was compared with that of the BP-ANN, LS-SVM, and conventional RVM methods. The main conclusions are drawn as follows:(1)
The penetration rate of the TBM is affected by various factors, and there are complex nonlinear relationships between the penetration rate and its influencing factors. The RVM-PSO method developed in this study can effectively establish a nonlinear correlation between the penetration rate and its main influencing factors and can achieve a relatively good prediction of the TBM penetration rate.(2)
The RVM-PSO method proves able to provide much better prediction than the BP-ANN, LS-SVM, and RVM. In addition, the BP-ANN and the LS-SVM methods have close prediction performances, and they are inferior to the RVM. Moreover, even if there is a lack of data on an important factor affecting the penetration rate, the prediction yielded by the RVM-PSO method is still more accurate than that using the RVM.(3)
The RVM-PSO method is useful in identifying the difference in the importance of the various factors affecting the TBM penetration rate prediction. It can help researchers/engineers to pay attention to the parameters that are significant for the penetration rate prediction and to leave out the data acquisition of the nonsignificant influencing factors in order to reduce workload.



---

*Source: 3268694-2020-09-04.xml*