# A Novel Model Using Virtual State Variables and Bayesian Discriminant Analysis to Classify Surrounding Rock Stability

Authors: Jinglai Sun; Darui Ren; Yu Song; Mingyuan Yu; Zhaofei Chu; Baoguo Liu; Shaogang Li; Xinyang Guo
Journal: Shock and Vibration
 (2021)

---

## Abstract

To accurately classify the stability of surrounding rock masses, a novel method (VSV-BDA) based on virtual state variables (VSVs) and Bayesian discriminant analysis (BDA) is proposed. The factors influencing stability are mapped by an artificial neural network (ANN) capable of recognizing the model of rock mass classification, and the obtained output vector is treated as VSVs, which are verified as obeying a multinormal distribution with equal covariance matrixes by normal distribution testing and constructed statistics. The prediction variance ratio test method is introduced to determine the optimal dimension of the VSVs. The VSV-BDA model is constructed through the use of VSVs and the optimal dimension on the basis of the training samples, which are divided from the collected samples into three situations having different numbers. ANN and BDA models are also constructed based on the same training samples. The predictions by the three models for the testing samples are compared; the results show that the proposed VSV-BDA model has high prediction accuracy and can be applied in practical engineering.


---

## Body

## 1. Introduction

Rock mass classification is generally considered a usable and practical approach to evaluating the stability of a rock mass in underground engineering [1, 2]. Reasonable classification can reflect the mechanical characteristics of the rock mass and provide a reliable basis for the design of underground engineering excavations and supporting systems [3–5]. The classification approach has gradually evolved from the beginning with a single index to the current generation with multiple quantitative and/or qualitative indexes.Rock quality design (RQD) was an efficient method proposed for rock mass quality assessment [6, 7]. The rock structure rating (RSR) system was proposed for tunnel support design [8]. RSR was further developed into rock mass rating (RMR), a portion-rating system [9, 10]. A rock tunneling quality index called the Q-system has been widely used in underground engineering and is closely related to RMR [11]. In recent times, with the development of numerical analysis and computer science, researchers have extensively studied rock mass classification and stability evaluation approaches: fuzzy analytic hierarchy process [12], artificial neural networks (ANNs) [13, 14], distance discriminant analysis [15, 16], set pair analysis [17], and so on.The stability of surrounding rocks is determined by many factors including geological conditions, exploitation factors, rock properties, and hydrogeological conditions [18–24]. These factors have mutual effects and complicated nonlinear relationships. ANNs can represent nonlinear relationships, but a sufficient number of available training samples is needed [25, 26]. Bayesian discriminant analysis (BDA) is a statistical method to classify samples that can use prior probability to increase the accuracy, but the relationships are too simple [27, 28].A novel method combining virtual state variables (VSVs) and Bayesian discriminant analysis (BDA), denoted as the VSV-BDA model, is proposed. The characteristics include the following: (1) the VSVs and prior probability is based to increase the accuracy of prediction; (2) quantile-quantile (QQ) plots and constructed statistics are used to test the normal distribution and equality of covariance instead of assuming a normal distribution and an equal covariance.The remainder of the paper includes three parts. In Section2, the VSV-BDA approach is proposed with virtual state variables determination, Bayesian discriminant theory, and evaluation of the VSV-BDA. In Section 3, the proposed approach is used in rock mass classification and testing the accuracy. In Section 4, existing problems are analyzed and further study is analyzed.


## 2. VSV-BDA Approach

The factors influencing rock mass stability have nonlinear relationships [11], and ANNs can recognize these relationships [29–31]. The output of training samples predicted by an ANN model, denoted as Y−, is compared with the actual output by the residual variance ratio (RVR) method to determine the construction [32, 33]. Then, the first k−1 layers (except for the output layer) and the regulating unit are used to construct a recognition network [34]. The output vector of the recognition network through nonlinear transformation of the input vector is called the VSVs by ANN, with the i-th variable denoted as vi. The vector with l components is called the virtual state vector, denoted as V, which does not have physical meaning but contains particular information characteristics to classify the stability, as shown in Figure 1. In addition, these variables instead of the influencing factors are used to construct the BDA model to classify the stability.Figure 1
Forecasting method.The processes of VSV-BDA are summarized as follows: (1) collect data and group them into training samples and testing samples; (2) construct a multilayer ANN to determine the VSVs (influencing factors); (3) calculate the Mahalanobis distance of VSVs and use the QQ plot to normality test; the data should be transformed to satisfy the multinormal distribution by a Box-Cox transformation; (4) create discriminant functions based on the VSV of the training samples and cross-validation to estimate the accuracy; (5) use the constructed discriminant function to calculate the scores of the testing samples belonging to the collectivities; (6) classify the testing samples associated with the highest posterior probability and calculate the ratio of misclassification. The process of BDA is illustrated in Figure2.Figure 2
Flow chart of BDA.


### 2.1. Virtual State Variable Determination

Interconnected processing elements, called neurons or cells, are used to construct the ANN, which offers a computational mechanism to acquire, compute, and represent a mapping from one information space layer to another, and to obtain a dataset to represent the relationships. An ANN identifies relationships by focusing on the parallel processing of many simple units, which can describe a complex function when combined. Essentially, the ANN is a gathering of simple processing units that exchange information that can be modified and filtered by the processing units’ connections [35].A multilayer ANN, which consists of an input layer, an output layer, and one or more hidden layers, is used to explain how to determine the VSVs. Factors influencing rock stability are taken as the input vector, represented byx=x1,x2,…,xl. To simplify, a typical representation with one hidden layer whose function is denoted by fu is illustrated in Figure 2. The output is the resulting rock mass classification, denoted as Y− (Figure 3).Figure 3
Pattern recognition network.The input dataset of the hidden layer isuj=∑i=1lUijxi, and the output dataset is fjuj with a weight index wj. When the hidden layer function is a wavelet function or a linear function, the output Y− can be denoted as equation (1).(1)Y−=∑i=1lwjfjuj.The common transfer function in this paper is the sigmoid function [36], so equation (1) is replaced by equation (2).(2)Y−=wS∑fjuj+θ,where w=1 and θ is a threshold value. These parameters can be calculated by the ANN training process [29, 37]. As aforementioned, the output of the hidden layer fjuj is used as the VSV (vj), and the virtual state vector is represented by V=v1,v2,v3,…,vlT.The optimal dimension is determined by RVR, which is illustrated by equation (3). The initial state of V dimension is Λ, and the calculated RVR is denoted as SeΛ. When the dimension is decreased or increased, called a state change and denoted as Λ↔, the RVR is represented by SeΛ↔.(3)Se=1N−1∑k=1NYANN−k−Yk2,where YANN−k is the prediction by the ANN model, Yk is the actual stability level, and N is the number of training samples.WhenSe is supposed to follow a normal distribution, to test the availability, a testing variable F is illustrated by(4)F=SeΛ↔SeΛ.For the given significance levelα, the processes of RVR are shown as follows:(1)
WhenF>FαN−1,N−1, the RVR has a significant increase, and the prediction performance is worse, so the new state is invalid(2)
WhenF<F1−αN−1,N−1, the RVR has a significant decrease, and the prediction performance is better, so the new state is valid(3)
In the third situation, the RVR is changing, but it cannot decide for the better or the worse, so the state remains unchanged


### 2.2. Bayesian Discriminant Theory (BDA)

Bayesian discriminant analysis is a probability analysis method with various types of distribution density functions that should be obtained at the beginning. The prior distribution is used to describe the awareness level of training samples before extracting the testing samples; then, the posterior distribution is obtained by modifying the prior distribution from the testing samples.Suppose there arek collectivities with p member indexes (considering p indexes):G1,G2,…,Gkk≥2, and the covariance matrix Σi>0i=1,2,…,k. The prior probability of Gi is denoted as pii=1,2,…,k and allocated by the proportion of Gi to all collectivities with ∑i=1kpi=1, as shown in(5)pi=nin,i=1,2,…,k,where ni is the number of the training samples belonging to Gi and n is the number of training samples.According to Bayesian theory, the posterior probability of sample (VSV from ANN)V=v1,v2,v3,…,vlT belonging to collectivity Gi is(6)PGi|V=pifiV∑i=1kpifiV,where fiX is the distribution function of Gi.Then, the optimal belonging can be obtained as(7)V∈Gj:PjV=maxPiV1≤i≤k.


#### 2.2.1. Normal Distribution Testing

The square difference of the Mahalanobis distance, denoted asdi2V, is used to represent the distance between sample V and collectivities Gi.(8)di2V=V−μiT∑i−1V−μi,where μi is the expect vector of Gi.The square differences of the Mahalanobis distances are sorted from the smallest to the largest.(9)d12≤d22≤⋯≤dn2.A QQ plot, which is represented by the pointsdi2,χi2, where χi2 is the chi-square distribution, is used to test whether the collectivities obey the multinormal distribution. When the points are all near the line passing through the origin with slope equal to one, the collectivities are regarded as obeying the multinormal distribution; otherwise, the data should be transformed to satisfy the multinormal distribution by a Box-Cox transformation. [33].


#### 2.2.2. Discriminant Criterion

Through the testing by the QQ plot and data transformation, all collectivities would obey the normal distribution; i.e.,Gi∼Npμi,Σi,i=1,2,…,kk≥2.(1)
WhenΣ1=Σ2=⋯=Σk=Σ, the distribution function is(10)fiV=2π−0.5pΣ−0.5exp−0.5V−μiTΣ−1V−μi.The discriminant function (WiV) of collectivity Gi is linear, and the best divisions (Rj) are obtained by equations (10), (11), and (14).(11)WiV=aiTV+bi,whereai=Σ−1μi, bi=−0.5μiTΣ−1μi+lnpi.(2)
WhenΣ1≠Σ2≠⋯≠Σk≠Σ, the distribution function is(12)fiX=2π−0.5pΣi−0.5exp−0.5V−μiTΣi−1V−μi.The discriminant function (WiX) of collectivity Gi is quadratic, and the best divisions are obtained by equations (12)–(14).(13)di2V=V−μiTΣi−1V−μi+lnΣi−2lnpi,(14)Rj=V:dj2V=mindj2V,1≤i≤k=V:PjGj|V=maxPiGi|V,1≤i≤k,j=1,2,…,k.


#### 2.2.3. Estimation of Parameters

In fact,μi and Σi are unknown, and an unbiased estimation can be obtained from the training samples. Supposing that training samples Vi=V1i,V2i,…,VniiT belong to Gi, μi is defined as V¯i.(15)V¯i=1ni∑j=1niVji.Σi is defined as Si.(16)Si=1ni−1∑j=1niVi−V¯iVi−V¯iT,i=1,2,…,k.The statistic used to examine the equality of the covariance matrix is defined as(17)ξ=1−DM∼χ2f,(18)M=n−klnS−∑i=1kni−1lnSi,(19)S=∑i=1kni−1Sin−k,(20)f=pp+1k−12,(21)D=2p2+3p−16p+1k−1∑i=1k1ni−1−1n−k,niare not completely equal,2p2+3p−1k+16p+1k−1,niare completely equal.For the givenα, the probability p=Pξ>χα2f is calculated. If p>α, Si are completely equal; otherwise, they are not completely equal.


### 2.3. Evaluation of the VSV-BDA

The resubstitution method and cross-validation method can also be used to estimate the reliability of the constructed discriminant criterion [33, 38–40]. For high accuracy, cross-validation is used, the principle of which is to choose one sample as a testing sample and use the rest of the collectivities as training samples to construct the discriminant criteria. The constructed criteria are used to classify the testing sample. The processes are shown as follows:(1)
Choose one sample fromG1 as the testing sample and construct the discriminant criterion with the other n−1 samples.(2)
Classify the testing sample with the criterion constructed in process (1).(3)
Repeat processes (1) and (2), and define the misclassification number asN1 after all the samples of G1 are tested.(4)
Repeat the processes (1), (2), and (3) for collectivityGi1≤i≤k and define the misclassification number as Ni after all the samples of Gi are tested.(5)
The ratio of misclassificationη can be calculated by(22)η=N1+N2+⋯+Nkn.


## 2.1. Virtual State Variable Determination

Interconnected processing elements, called neurons or cells, are used to construct the ANN, which offers a computational mechanism to acquire, compute, and represent a mapping from one information space layer to another, and to obtain a dataset to represent the relationships. An ANN identifies relationships by focusing on the parallel processing of many simple units, which can describe a complex function when combined. Essentially, the ANN is a gathering of simple processing units that exchange information that can be modified and filtered by the processing units’ connections [35].A multilayer ANN, which consists of an input layer, an output layer, and one or more hidden layers, is used to explain how to determine the VSVs. Factors influencing rock stability are taken as the input vector, represented byx=x1,x2,…,xl. To simplify, a typical representation with one hidden layer whose function is denoted by fu is illustrated in Figure 2. The output is the resulting rock mass classification, denoted as Y− (Figure 3).Figure 3
Pattern recognition network.The input dataset of the hidden layer isuj=∑i=1lUijxi, and the output dataset is fjuj with a weight index wj. When the hidden layer function is a wavelet function or a linear function, the output Y− can be denoted as equation (1).(1)Y−=∑i=1lwjfjuj.The common transfer function in this paper is the sigmoid function [36], so equation (1) is replaced by equation (2).(2)Y−=wS∑fjuj+θ,where w=1 and θ is a threshold value. These parameters can be calculated by the ANN training process [29, 37]. As aforementioned, the output of the hidden layer fjuj is used as the VSV (vj), and the virtual state vector is represented by V=v1,v2,v3,…,vlT.The optimal dimension is determined by RVR, which is illustrated by equation (3). The initial state of V dimension is Λ, and the calculated RVR is denoted as SeΛ. When the dimension is decreased or increased, called a state change and denoted as Λ↔, the RVR is represented by SeΛ↔.(3)Se=1N−1∑k=1NYANN−k−Yk2,where YANN−k is the prediction by the ANN model, Yk is the actual stability level, and N is the number of training samples.WhenSe is supposed to follow a normal distribution, to test the availability, a testing variable F is illustrated by(4)F=SeΛ↔SeΛ.For the given significance levelα, the processes of RVR are shown as follows:(1)
WhenF>FαN−1,N−1, the RVR has a significant increase, and the prediction performance is worse, so the new state is invalid(2)
WhenF<F1−αN−1,N−1, the RVR has a significant decrease, and the prediction performance is better, so the new state is valid(3)
In the third situation, the RVR is changing, but it cannot decide for the better or the worse, so the state remains unchanged


## 2.2. Bayesian Discriminant Theory (BDA)

Bayesian discriminant analysis is a probability analysis method with various types of distribution density functions that should be obtained at the beginning. The prior distribution is used to describe the awareness level of training samples before extracting the testing samples; then, the posterior distribution is obtained by modifying the prior distribution from the testing samples.Suppose there arek collectivities with p member indexes (considering p indexes):G1,G2,…,Gkk≥2, and the covariance matrix Σi>0i=1,2,…,k. The prior probability of Gi is denoted as pii=1,2,…,k and allocated by the proportion of Gi to all collectivities with ∑i=1kpi=1, as shown in(5)pi=nin,i=1,2,…,k,where ni is the number of the training samples belonging to Gi and n is the number of training samples.According to Bayesian theory, the posterior probability of sample (VSV from ANN)V=v1,v2,v3,…,vlT belonging to collectivity Gi is(6)PGi|V=pifiV∑i=1kpifiV,where fiX is the distribution function of Gi.Then, the optimal belonging can be obtained as(7)V∈Gj:PjV=maxPiV1≤i≤k.


### 2.2.1. Normal Distribution Testing

The square difference of the Mahalanobis distance, denoted asdi2V, is used to represent the distance between sample V and collectivities Gi.(8)di2V=V−μiT∑i−1V−μi,where μi is the expect vector of Gi.The square differences of the Mahalanobis distances are sorted from the smallest to the largest.(9)d12≤d22≤⋯≤dn2.A QQ plot, which is represented by the pointsdi2,χi2, where χi2 is the chi-square distribution, is used to test whether the collectivities obey the multinormal distribution. When the points are all near the line passing through the origin with slope equal to one, the collectivities are regarded as obeying the multinormal distribution; otherwise, the data should be transformed to satisfy the multinormal distribution by a Box-Cox transformation. [33].


### 2.2.2. Discriminant Criterion

Through the testing by the QQ plot and data transformation, all collectivities would obey the normal distribution; i.e.,Gi∼Npμi,Σi,i=1,2,…,kk≥2.(1)
WhenΣ1=Σ2=⋯=Σk=Σ, the distribution function is(10)fiV=2π−0.5pΣ−0.5exp−0.5V−μiTΣ−1V−μi.The discriminant function (WiV) of collectivity Gi is linear, and the best divisions (Rj) are obtained by equations (10), (11), and (14).(11)WiV=aiTV+bi,whereai=Σ−1μi, bi=−0.5μiTΣ−1μi+lnpi.(2)
WhenΣ1≠Σ2≠⋯≠Σk≠Σ, the distribution function is(12)fiX=2π−0.5pΣi−0.5exp−0.5V−μiTΣi−1V−μi.The discriminant function (WiX) of collectivity Gi is quadratic, and the best divisions are obtained by equations (12)–(14).(13)di2V=V−μiTΣi−1V−μi+lnΣi−2lnpi,(14)Rj=V:dj2V=mindj2V,1≤i≤k=V:PjGj|V=maxPiGi|V,1≤i≤k,j=1,2,…,k.


### 2.2.3. Estimation of Parameters

In fact,μi and Σi are unknown, and an unbiased estimation can be obtained from the training samples. Supposing that training samples Vi=V1i,V2i,…,VniiT belong to Gi, μi is defined as V¯i.(15)V¯i=1ni∑j=1niVji.Σi is defined as Si.(16)Si=1ni−1∑j=1niVi−V¯iVi−V¯iT,i=1,2,…,k.The statistic used to examine the equality of the covariance matrix is defined as(17)ξ=1−DM∼χ2f,(18)M=n−klnS−∑i=1kni−1lnSi,(19)S=∑i=1kni−1Sin−k,(20)f=pp+1k−12,(21)D=2p2+3p−16p+1k−1∑i=1k1ni−1−1n−k,niare not completely equal,2p2+3p−1k+16p+1k−1,niare completely equal.For the givenα, the probability p=Pξ>χα2f is calculated. If p>α, Si are completely equal; otherwise, they are not completely equal.


## 2.2.1. Normal Distribution Testing

The square difference of the Mahalanobis distance, denoted asdi2V, is used to represent the distance between sample V and collectivities Gi.(8)di2V=V−μiT∑i−1V−μi,where μi is the expect vector of Gi.The square differences of the Mahalanobis distances are sorted from the smallest to the largest.(9)d12≤d22≤⋯≤dn2.A QQ plot, which is represented by the pointsdi2,χi2, where χi2 is the chi-square distribution, is used to test whether the collectivities obey the multinormal distribution. When the points are all near the line passing through the origin with slope equal to one, the collectivities are regarded as obeying the multinormal distribution; otherwise, the data should be transformed to satisfy the multinormal distribution by a Box-Cox transformation. [33].


## 2.2.2. Discriminant Criterion

Through the testing by the QQ plot and data transformation, all collectivities would obey the normal distribution; i.e.,Gi∼Npμi,Σi,i=1,2,…,kk≥2.(1)
WhenΣ1=Σ2=⋯=Σk=Σ, the distribution function is(10)fiV=2π−0.5pΣ−0.5exp−0.5V−μiTΣ−1V−μi.The discriminant function (WiV) of collectivity Gi is linear, and the best divisions (Rj) are obtained by equations (10), (11), and (14).(11)WiV=aiTV+bi,whereai=Σ−1μi, bi=−0.5μiTΣ−1μi+lnpi.(2)
WhenΣ1≠Σ2≠⋯≠Σk≠Σ, the distribution function is(12)fiX=2π−0.5pΣi−0.5exp−0.5V−μiTΣi−1V−μi.The discriminant function (WiX) of collectivity Gi is quadratic, and the best divisions are obtained by equations (12)–(14).(13)di2V=V−μiTΣi−1V−μi+lnΣi−2lnpi,(14)Rj=V:dj2V=mindj2V,1≤i≤k=V:PjGj|V=maxPiGi|V,1≤i≤k,j=1,2,…,k.


## 2.2.3. Estimation of Parameters

In fact,μi and Σi are unknown, and an unbiased estimation can be obtained from the training samples. Supposing that training samples Vi=V1i,V2i,…,VniiT belong to Gi, μi is defined as V¯i.(15)V¯i=1ni∑j=1niVji.Σi is defined as Si.(16)Si=1ni−1∑j=1niVi−V¯iVi−V¯iT,i=1,2,…,k.The statistic used to examine the equality of the covariance matrix is defined as(17)ξ=1−DM∼χ2f,(18)M=n−klnS−∑i=1kni−1lnSi,(19)S=∑i=1kni−1Sin−k,(20)f=pp+1k−12,(21)D=2p2+3p−16p+1k−1∑i=1k1ni−1−1n−k,niare not completely equal,2p2+3p−1k+16p+1k−1,niare completely equal.For the givenα, the probability p=Pξ>χα2f is calculated. If p>α, Si are completely equal; otherwise, they are not completely equal.


## 2.3. Evaluation of the VSV-BDA

The resubstitution method and cross-validation method can also be used to estimate the reliability of the constructed discriminant criterion [33, 38–40]. For high accuracy, cross-validation is used, the principle of which is to choose one sample as a testing sample and use the rest of the collectivities as training samples to construct the discriminant criteria. The constructed criteria are used to classify the testing sample. The processes are shown as follows:(1)
Choose one sample fromG1 as the testing sample and construct the discriminant criterion with the other n−1 samples.(2)
Classify the testing sample with the criterion constructed in process (1).(3)
Repeat processes (1) and (2), and define the misclassification number asN1 after all the samples of G1 are tested.(4)
Repeat the processes (1), (2), and (3) for collectivityGi1≤i≤k and define the misclassification number as Ni after all the samples of Gi are tested.(5)
The ratio of misclassificationη can be calculated by(22)η=N1+N2+⋯+Nkn.


## 3. VSV-BDA Model for Rock Mass Classification

### 3.1. Influencing Factors of Rock Stability and Sample Collection

Considering the typical factors, the convenience of factors can be obtained and compared in practical engineering. Five factors are selected as influencing factors, including the RQD, uniaxial axial compressive strength of the rock (Rw), rock mass integrity index (Kv), coefficient of structural surface strength (Kf), and groundwater discharge (ω) [41]. The rock stability is divided into five levels presented by collectivities G1∼G5: G1 is level I of rock stability; G2 is level II; G3 is level III; G4 is level IV; and G5 is level V. The discriminant criterion belonging to Gi is Wii=1,2,…,5.The surrounding rock classifications in the second-stage project of the Guangzhou pumped storage power station are collected to construct the VSV-BDA model and validation, as shown in Table1 [41].Table 1
The collected samples.
No.RQD (%)x1Rw (MPa) x2Kvx3Kfx4ω (L/min.10 m) x5Actual output143.826.30.280.3550IV237.522.50.250.375IV331.318.80.230.25100IV426360.220.355IV55040.20.50.510III652250.20.55III771900.350.318III875950.70.50II977.5900.570.4510II1050700.50.255III1150.9340.320.3521III1231.5200.230.2546IV1350.6260.260.3520III1475.5900.450.528II1585.5940.650.550II1680950.50.450II173570.50.350.310III1850900.50.255III1985930.60.50II2078.5920.550.56II2130.2700.40.210III2285920.70.510II2351350.320.3515III2487950.50.450II2582960.750.350II2652.530.50.380.1525IV2778800.750.50II2885920.70.510II2950450.120.35III3031.5200.230.2550IV3176.5900.550.510II3285950.650.50II3350.9340.320.3520III3475900.550.57II3575900.550.510II3650450.150.355III3775700.720.5312II


### 3.2. VSV-BDA Construction

To test the generalization of the VSV-BDA model and ensure the reliability of the results, the situation with too many training samples and too few testing samples should be avoided. The collected samples in Table1 are divided into three different situations to construct VSV-BDA models with ratios of training samples to testing samples equal to 20 : 17, 25 : 12, and 30 : 7.Taking the first ratio as an example, no. 1-20 are assigned as the training samples, and no. 21-37 are assigned as the testing samples. The method inSection 2 is applied. An ANN with 4 layers is constructed to determine the VSVs. The significance level is 0.05 (α=0.05). The input layer and two hidden layers of ANN compose the first 3 layers of the recognition network. The five influencing factors listed in Table 1 are taken as the input vector. The transfer function implemented in the two hidden layers is the sigmoid function. When the number of training samples is smaller than the dimensions of the VSVs, the covariance estimation may be ill-posed [42]. Due to the limitation of collected samples, the number of VSVs can be decreased by increasing the number of hidden layers [32]. The optimal dimension of VSVs is determined by RVR, which has a structure of 5×6×5×1. Then, the five VSVs obtained from the second hidden layer are used to construct the VSV-BDA model by the BDA theory.An ANN [41] and BDA are also used to create the model with the training samples to obtain the classification functions or variables.


### 3.3. Validation and Comparison

The testing samples are predicted by the three models to gain the output and compare the accuracy. The outputs of these 3 situations are shown in Figures4–6. The predictions of the training samples through cross-validation of the three models are the same as the actual output. From Figures 5 and 6, the outputs of the VSV-BDA model in situations 2 and 3 are completely identical to the actual output. Obviously, the VSV-BDA model has higher prediction accuracy and estimation accuracy than the other approaches when the samples are limited. All models have the same accuracy when the samples are sufficient. We conclude that the VSV-BDA model can be applied to surrounding rock mass classification with high prediction accuracy.Figure 4
Training effects and the prediction results for situation 1.Figure 5
Training effects and the prediction results for situation 2.Figure 6
Training effects and the prediction results for situation 3.In situation 1, the inaccurate predictions of testing samples occur for (1) samples 21, 26, and 31 by ANN, (2) samples 26 and 31 by BDA, and (3) sample 26 by VSV-BDA. The average relative error ratios of testing samples are 5.88% for ANN, 3.92% for BDA, and 1.96% for VSV-BDA, as shown in Figure4.In situation 2, the inaccurate prediction of samples occurs for sample 26 by ANN and by BDA. The average relative error ratios of testing samples are 3.03% for ANN and 3.03% for BDA, as shown in Figure5.In situation 3, the predictions by all three methods are correct. The three models have the same accuracy, which indicates that they can yield the correct result, as shown in Figure6.Compared with the ANN, the VSV-BDA model can use prior probability to increase the accuracy. Compared with BDA, the proposed method can represent complex relationships between the factors. From the aforementioned points, the VSV-BDA model has higher accuracy than the ANN and BDA models, especially when the training samples are insufficient. The more training samples are available, the more accurate the prediction. For the three methods, as the number of training samples increases, the accuracy increases.


### 3.4. Application of VSV-BDA Discriminant Criteria

Due to the limitation of space, only situation 1 is explained in this paper. The samples in Table1 are normalized. Through the constructed ANN, from the second hidden layer, five VSVs are obtained, denoted as V=v1,v2,v3,v4,v5. After calculating the Mahalanobis distance of the VSVs of training samples, the QQ plot shown in Figure 7 indicates that the points are all near the line passing through the origin with slope equal to one. The plot illustrates that the samples all obey the multinormal distribution.Figure 7
QQ plot of training samples.The prior probability is allocated by the proportion of the training samples:p1=0/20,p2=7/20, p3=8/20, p4=5/20, and p5=0/20. Due to the limitation of the samples, three collectivities and discriminant criteria have been collected and constructed.Equations (15)–(21) are calculated with the training samples in Table 1. The results illustrate, through equation (23), that the covariances of the three collectivities are equal, so the discriminant functions of rock stability are linear functions.(23)p=PT=375>χ260=79=1.Based on equations (10), (16), and (17), the distribution function of the rock stability is(24)fiV=2π−0.5pS−0.5exp−0.5V−V¯iTS−1V−V¯i.Equation (11) is used to obtain the discriminant criteria, illustrated by equations (25)–(27).(25)W2V=−11.42v1+6.05v2−26.69v3−6.3v4+5.5626v5−26.18,(26)W3V=−10.11v1+4.06v2−18.53v3+27.36v4+20.08v5−18.91,(27)W4V=−7.85v1−2.44v2−15.07v3+26.58v4+14.59v5−15.73.The training samples are used to examine the accuracy of the discriminant criteria by the cross-validation method, as shown in Figure4. The probability of misclassification is 0, which is perfect for classifying the rock mass stability.The testing samples are no. 21-37. The discriminant criteria, equations (25)–(27), are used to classify the testing samples. The output of sample no. 26 is level III but should be level IV, so the average relative error ratio is 1.96%.


## 3.1. Influencing Factors of Rock Stability and Sample Collection

Considering the typical factors, the convenience of factors can be obtained and compared in practical engineering. Five factors are selected as influencing factors, including the RQD, uniaxial axial compressive strength of the rock (Rw), rock mass integrity index (Kv), coefficient of structural surface strength (Kf), and groundwater discharge (ω) [41]. The rock stability is divided into five levels presented by collectivities G1∼G5: G1 is level I of rock stability; G2 is level II; G3 is level III; G4 is level IV; and G5 is level V. The discriminant criterion belonging to Gi is Wii=1,2,…,5.The surrounding rock classifications in the second-stage project of the Guangzhou pumped storage power station are collected to construct the VSV-BDA model and validation, as shown in Table1 [41].Table 1
The collected samples.
No.RQD (%)x1Rw (MPa) x2Kvx3Kfx4ω (L/min.10 m) x5Actual output143.826.30.280.3550IV237.522.50.250.375IV331.318.80.230.25100IV426360.220.355IV55040.20.50.510III652250.20.55III771900.350.318III875950.70.50II977.5900.570.4510II1050700.50.255III1150.9340.320.3521III1231.5200.230.2546IV1350.6260.260.3520III1475.5900.450.528II1585.5940.650.550II1680950.50.450II173570.50.350.310III1850900.50.255III1985930.60.50II2078.5920.550.56II2130.2700.40.210III2285920.70.510II2351350.320.3515III2487950.50.450II2582960.750.350II2652.530.50.380.1525IV2778800.750.50II2885920.70.510II2950450.120.35III3031.5200.230.2550IV3176.5900.550.510II3285950.650.50II3350.9340.320.3520III3475900.550.57II3575900.550.510II3650450.150.355III3775700.720.5312II


## 3.2. VSV-BDA Construction

To test the generalization of the VSV-BDA model and ensure the reliability of the results, the situation with too many training samples and too few testing samples should be avoided. The collected samples in Table1 are divided into three different situations to construct VSV-BDA models with ratios of training samples to testing samples equal to 20 : 17, 25 : 12, and 30 : 7.Taking the first ratio as an example, no. 1-20 are assigned as the training samples, and no. 21-37 are assigned as the testing samples. The method inSection 2 is applied. An ANN with 4 layers is constructed to determine the VSVs. The significance level is 0.05 (α=0.05). The input layer and two hidden layers of ANN compose the first 3 layers of the recognition network. The five influencing factors listed in Table 1 are taken as the input vector. The transfer function implemented in the two hidden layers is the sigmoid function. When the number of training samples is smaller than the dimensions of the VSVs, the covariance estimation may be ill-posed [42]. Due to the limitation of collected samples, the number of VSVs can be decreased by increasing the number of hidden layers [32]. The optimal dimension of VSVs is determined by RVR, which has a structure of 5×6×5×1. Then, the five VSVs obtained from the second hidden layer are used to construct the VSV-BDA model by the BDA theory.An ANN [41] and BDA are also used to create the model with the training samples to obtain the classification functions or variables.


## 3.3. Validation and Comparison

The testing samples are predicted by the three models to gain the output and compare the accuracy. The outputs of these 3 situations are shown in Figures4–6. The predictions of the training samples through cross-validation of the three models are the same as the actual output. From Figures 5 and 6, the outputs of the VSV-BDA model in situations 2 and 3 are completely identical to the actual output. Obviously, the VSV-BDA model has higher prediction accuracy and estimation accuracy than the other approaches when the samples are limited. All models have the same accuracy when the samples are sufficient. We conclude that the VSV-BDA model can be applied to surrounding rock mass classification with high prediction accuracy.Figure 4
Training effects and the prediction results for situation 1.Figure 5
Training effects and the prediction results for situation 2.Figure 6
Training effects and the prediction results for situation 3.In situation 1, the inaccurate predictions of testing samples occur for (1) samples 21, 26, and 31 by ANN, (2) samples 26 and 31 by BDA, and (3) sample 26 by VSV-BDA. The average relative error ratios of testing samples are 5.88% for ANN, 3.92% for BDA, and 1.96% for VSV-BDA, as shown in Figure4.In situation 2, the inaccurate prediction of samples occurs for sample 26 by ANN and by BDA. The average relative error ratios of testing samples are 3.03% for ANN and 3.03% for BDA, as shown in Figure5.In situation 3, the predictions by all three methods are correct. The three models have the same accuracy, which indicates that they can yield the correct result, as shown in Figure6.Compared with the ANN, the VSV-BDA model can use prior probability to increase the accuracy. Compared with BDA, the proposed method can represent complex relationships between the factors. From the aforementioned points, the VSV-BDA model has higher accuracy than the ANN and BDA models, especially when the training samples are insufficient. The more training samples are available, the more accurate the prediction. For the three methods, as the number of training samples increases, the accuracy increases.


## 3.4. Application of VSV-BDA Discriminant Criteria

Due to the limitation of space, only situation 1 is explained in this paper. The samples in Table1 are normalized. Through the constructed ANN, from the second hidden layer, five VSVs are obtained, denoted as V=v1,v2,v3,v4,v5. After calculating the Mahalanobis distance of the VSVs of training samples, the QQ plot shown in Figure 7 indicates that the points are all near the line passing through the origin with slope equal to one. The plot illustrates that the samples all obey the multinormal distribution.Figure 7
QQ plot of training samples.The prior probability is allocated by the proportion of the training samples:p1=0/20,p2=7/20, p3=8/20, p4=5/20, and p5=0/20. Due to the limitation of the samples, three collectivities and discriminant criteria have been collected and constructed.Equations (15)–(21) are calculated with the training samples in Table 1. The results illustrate, through equation (23), that the covariances of the three collectivities are equal, so the discriminant functions of rock stability are linear functions.(23)p=PT=375>χ260=79=1.Based on equations (10), (16), and (17), the distribution function of the rock stability is(24)fiV=2π−0.5pS−0.5exp−0.5V−V¯iTS−1V−V¯i.Equation (11) is used to obtain the discriminant criteria, illustrated by equations (25)–(27).(25)W2V=−11.42v1+6.05v2−26.69v3−6.3v4+5.5626v5−26.18,(26)W3V=−10.11v1+4.06v2−18.53v3+27.36v4+20.08v5−18.91,(27)W4V=−7.85v1−2.44v2−15.07v3+26.58v4+14.59v5−15.73.The training samples are used to examine the accuracy of the discriminant criteria by the cross-validation method, as shown in Figure4. The probability of misclassification is 0, which is perfect for classifying the rock mass stability.The testing samples are no. 21-37. The discriminant criteria, equations (25)–(27), are used to classify the testing samples. The output of sample no. 26 is level III but should be level IV, so the average relative error ratio is 1.96%.


## 4. Conclusions

(1)
The accuracy of the VSV-BDA model in rock stability classification depends on the prior probability, probability density, and complex relationships among the VSVs. To satisfy the actual conditions, the QQ plot and statistic are introduced to test the probability distribution and the equality of the covariance instead of the assumed multinormal distribution and equality used in other studies, which could yield a more reasonable distribution and criterion to classify the stability.(2)
Due to the limitation of the samples, only three collectivities yield the discriminant criterion based on the VSV-BDA. More samples should be collected to construct the other levels’ criteria.(3)
Compared with the results obtained with the ANN model and BDA, the results from the VSV-BDA model indicate that it has high prediction accuracy and can be used in practice. The VSV-BDA model provides a perfect approach to classify the surrounding rock mass stability.In the future, further study should be carried out to select the influencing factor and enhance the VSV-BDA model in practical engineering.



---

*Source: 6656882-2021-03-12.xml*