# The Rock Burst Hazard Evaluation Using Statistical Learning Approaches

Authors: Jie Chen; Jingkuan Gao; Yuanyuan Pu; Mingzhong Gao; Like Wei; Chong Wang; Bo Peng; Xusheng Zhao; Guangchao Zhang; Zhigang Zhang
Journal: Shock and Vibration
 (2021)

---

## Abstract

The great threat and destructiveness brought by a rock burst make its prediction and prevention crucial in engineering. The rock burst hazard evaluation at project locations is an effective way of preventing rock burst since currently real-time prediction is not available. Since different control factors and discrimination conditions of rock burst were accepted by conventional risk determination methods, the rock burst risk determination in the same area may produce conflicting results. In this study, Naive Bayes statistical learning models based on different model prior distributions representing highly complicated nonlinear relationship between rock burst hazard and impact factors were built to evaluate the rock burst hazards. The results suggested that the Bayes statistical learning model based on a Gaussian prior has the strongest performance over four preset prior distributions. Combining the rock mechanics parameters measured in the laboratory and the stress data collected on the project sites, the proposed model was successfully employed to evaluate the kimberlite rock burst risk of a diamond mine in Canada. The Bayes statistical learning model exhibits its robustness and generalization in rock burst hazard evaluation, which can be generalized for similar engineering cases with enough supported data.


---

## Body

## 1. Introduction

A rock burst is a kind of sudden and severe rock instability, referring to a dynamic geological disaster caused by the sudden release of elastic strain energy accumulated in the rock mass of the underground excavations [1, 2]. Rock burst ejects a large amount of rocks in a short time, which seriously endangers the safety of construction equipment and personnel. With the increase of the mining depth and the tunnel construction, the initial ground stress has seen a drastic increase, which is more likely to induce rock burst [3]. Since the first recorded rock burst that occurred in the United Kingdom in 1738, all mining countries over the world have recorded rock bursts, including China, Canada, the United States, and Australia [4]. The Brunswick lead-zinc mine in Canada, the Macassa gold mine, the Kalgoorlie gold mine in Australia, and the Idaho lead-zinc-silver mine in the United States have all experienced rock bursts that caused serious fatalities [5–7]. In 2018, 21 workers were killed and 4 were injured in a rock burst that occurred in the connecting lane at 1303 face of Long Yun Coal Mine, He Ze, Shandong Province, China. In addition to mining, rock bursts are also likely to occur in other underground excavations. For example, More than 1000 rock bursts of different levels occurred during the excavation of the No. 2 diversion tunnel in Jinping II Hydropower Station on the Ya Long River, which seriously affected the progress of the project.The hazard of the rock burst makes the research on it a hot topic in the field of rock mechanics and engineering. Researchers have been committed to predicting the exact time and location of rock burst based on the mechanism, so as to completely eliminate the threat of rock burst. However, the mechanism of rock burst is complicated, and its control factors are numerous including the mechanical properties of the rock, the field stress environment, construction support parameters, etc. Furthermore, the coupling relationships and interactions among those various influencing factors are complicated. Currently, there is no universal and reliable method for predicting rock burst in real time [8, 9]. Hence, before the construction of a project, the assessment for the rock burst risk in certain high-risk areas has become a major preventive method when rock burst cannot be accurately predicted in real time and space. Adjusting on-site excavation methods and support parameters with reference to the results of rock burst risk assessment is also an important prevention and control measure for rock burst. The evaluation of the rock burst risk is mainly to use the control factors of rock burst to comprehensively evaluate the possibility of the rock burst and the potential rock burst intensity. Conventional rock burst risk assessment methods can be divided into two categories, single index methods and comprehensive index methods. The single index methods use the selected control factors of rock burst to calculate a result value through the preset formula and then compare it with threshold value to determine the risk of a rock burst. The comprehensive index methods mainly use mathematical and statistical models to carry out a weighted mapping of the control factors of rock burst and calculate a comprehensive value to judge the risk of rock burst. The common synthetic index methods to evaluate rock burst risk include principal component analysis [10], fuzzy mathematics [11], and analytic hierarchy process [12]. Researchers have made achievements in evaluating rock burst hazards using single index methods and comprehensive index methods. However, due to some inherent defects of these two approaches, the universality is difficult to meet engineering needs, which is difficult to extend a successful risk assessment case to more engineering sites. For example, for single index methods, the calculated index is not only difficult to comprehensively reflect the mechanical behavior of the rock mass and the instability characteristics under extreme conditions but also requires presetting artificial threshold values for the risk assessment. However, in different rock burst cases, the threshold is determined by the constructor mainly based on experience. For example, both Kidybinski and Singh proposed to use strain energy storage index to judge the risk of coal rock burst. For coal mines in Silesia, Poland, Kidybinski’s judgment criteria are that the strain energy storage index less than 2.5 refers to no rock burst hazard; the strain energy storage index greater than 2.5 and less than 3.5 refers to medium rock burst hazard; the energy index greater than 3.5 and less than 5 is regarded as a strong rock burst hazard; the strain energy storage index greater than 5 refers to the violent rock burst hazard [13]. By contrast, Singh puts forward the rock burst hazard criteria for hard rock as follows: value less than 10 is regarded as no rock burst hazard; value greater than 10 and less than 15 is regarded as medium rock burst hazard; and value greater than 15 is regarded as strong rock burst hazard [14]. In terms of the comprehensive index method, it has to set a threshold value for risk judgment, as well as to assign weights to the selected control factors of rock burst. The aforementioned traditional rock burst risk assessment methods are greatly affected by the subjective judgment of researchers, which is difficult to accurately and objectively judge the risk of rock burst. Hence, the methods accepted in some cases are difficult to generalize to other rock burst cases.In view of the shortcomings of traditional risk evaluation approaches, it raises a new possibility to use the statistical learning method to judge the risk of rock burst only from case data. Statistical learning is a way of obtaining knowledge from existing data and predicting using new data. Commonly used statistical learning models include support vector machines, feedforward neural networks, logistic regression, and Bayesian methods. Rock burst risk assessment can be regarded as a multiclass supervised-learning task, which means to use the collected rock burst case data (including the values of the control factors of rock bursts and the corresponding intensity level) to train a statistical learning model and then feed new input data control factors for prediction. As an end-to-end data-driven method, statistical learning does not need any prior information about the data as well as considering the complex mapping relationship for the intermediate process in the process of model training and prediction. According to universal approximation theorem, a feedforward neural network with only one hidden layer, which embeds enough neurons, can approximate any continuous function on a compact subset ofRn with arbitrary precision [15]. Therefore, the statistical learning model is suitable for multiclassification tasks such as rock burst hazard evaluation which is not clear enough, involving many variables and complex mapping relationships.This study collected rock burst case data to train Naive Bayes classifiers under various prior distributions. By comparing the classification accuracy, the Bayesian classifier based on Gaussian distribution is selected as the rock burst risk assessment model. This model is used to evaluate the rock burst risk of kimberlite in a diamond mine in Canada. Compared with other statistical learning models, the Naive Bayes classifier does not significantly reduce the model accuracy when the training sample size has been reduced, which has strong adaptability to small sample tasks [16]. Rock burst cases are all over the world, and most of them are difficult to access case reports. Due to the lack of control factor data of some cases, the complete and high-quality data that can be collected are limited, which means that the rock burst hazard evaluation task is a small sample training task. Therefore, it is reasonable to use the Naive Bayes classifier to evaluate the rock burst hazard. The remainder of this article is organized as follows: Section 2 discusses the background theory of the Naive Bayes model; Section 3 gives the complete process of model construction, including data collection, model training, and validation; Section 4 uses the trained model to conduct rock burst hazard evaluation in a diamond mine; and finally, the conclusion and discussion are given in Section 5.


## 2. Bayesian Statistical Learning Model

The Naive Bayesian statistical classification method has excellent performance in many practical applications, such as file classification and spam filtering.Compared with most random processes, the Naive Bayes classification method can extract data features of each dimension more quickly, reducing the difficulty of high-dimensional data calculation. However, the probability calculation of Naive Bayes in practice is significantly different from the true probability, but the impact of this difference on the classification task can be ignored, and the classification effect is well represented. This method assumes that all attributes independently affect the results, given the sample attributesX = {x1, x2, …, xn}, and the Bayesian class conditional probability corresponding to the class label y is(1)Py|X=Py|x1,⋯,xn=PyPx1,⋯,xn|yPx1,⋯,xn.According to the independence hypothesis of Naive Bayes:(2)Pxi|y,x1,⋯,xi−1,xi+1,⋯,xn=Pxi|y.For attributeX, the probability of class Y is simplified as(3)Py|x1,⋯,xn=Py∏i=1nPxi|yPx1,⋯,xn=PyPX∏i=1dPxi|y,where n is the number of attributes and xi is the value of X on the i-th attribute.Generally, for all categories, the probabilityPX has the same value, so the Bayesian decision criterion based on the Bayesian optimization classifier is(4)hnbx=argmaxyPy∏i=1nPxi|y.Obviously, the judgment process of the Naive classifier of expression (4) is to estimate the prior probability Py of the class based on the dataset, and each attribute estimates the conditional probability Pxi|y.Through the basic formula of the above Bayes algorithm, the basic process is as follows:(1)
LetX=x1,x2,…,xn be the sample to be classified, where xi is a characteristic attribute of X.(2)
The set of categories to be classified isC=y1,y2,…,yn.(3)
Calculate the probability ofPy1|XPy2|X,…,Pyn|X, respectively.(4)
IfPyk|X=maxPy1|X,Pyk2|X,…,Pyn|X, then X is considered to be of type yk.According to the prior distribution that the data obey, the Naive Bayes classifier can be divided into Gaussian, multinomial, complement, Bernoulli, and other models. If the distribution of sample features is continuous, the attribute features can be assumed to obey Gaussian distribution:(5)Pxi|y=12πσy2exp−xi−μy22σy2.The parametersσy and μy can be calculated using maximum likelihood estimation (MLE).Schütze H et al. supposed that when the characteristic attributes are continuous values, which distribution obeys the Bernoulli distribution, the Bernoulli’s Naive Bayes classification algorithm can be introduced for calculating posterior distribution [4–6]. The classification decision rule is calculated based on the Bernoulli distribution probability as follows:(6)Pxi|y=Pi|yxi+1−Pi|y1−xi.Similarly, when the data prior distribution is polynomial distribution or other distribution, the corresponding Bayesian classifier can also be constructed correspondingly.The control factors of rock burst are numerous and the relationship is complicated. It is difficult to develop the mathematical laws that the control factors obey from the mechanism, which is difficult to predict the predetermined distribution of its characteristics. In this paper, Naive Bayesian models based on Gaussian, polynomial, complement, and Bernoulli distributions are established, respectively, performing statistical learning on the same dataset and finally adopting the prior distribution with the best performance as that prior distribution of the Bayesian classifier.In general, the Naive Bayes classifier can only conduct binary classification task. The risk assessment of rock burst is a multiclass classification task, which has to construct a Naive Bayes classifier suitable for multiclass classification. This article establishes multiple basic binary classifiers and finally completes the multiclass classification task through a “voting” mechanism. The basic idea is to establish a binary classifier between each pair of all categories, that is to say, there aren categories to be classified, and n (n − 1)/2 classifiers are established. Let Cij represent the classifier between category i and category j. For a single training sample X, if the classification result belongs to i, then class i gets “one vote.” Otherwise, category j gets “one vote.” Until all n (n − 1)/2 classifiers have voted for X, X will belong to the category with the most votes.


## 3. Model Construction and Training

### 3.1. Sample Data Collection

Rock burst is caused by the concentration of on-site stress exceeding the energy storage limit of the rock mass. There are two main control factors of a rock burst, the field stress indexes and the rock mass properties. In order to reflect the nature and characteristics of rock burst comprehensively, this article adopted the tangential stress of surrounding rockσθ (MPa), the ratio of tangential stress to uniaxial tensile strength σθ/σt, the ratio of uniaxial compressive strength to uniaxial tensile strength σc/σt, and elastic energy index WET as the attributes of the data sample to evaluate the risk of rock burst. According to the commonly accepted standards for ranking a rock burst grades, a rock burst can be divided into four grades [17]. That is, no rock burst, moderate rock burst, strong rock burst, and severe rock burst. These four levels are used as data labels to participate in the construction of statistical learning models by digitizing. 0 means no rock burst risk, 1 means medium rock burst risk, 2 means strong rock burst risk, and 3 means violent rock burst risk.This article collected 111 rock burst cases from different underground excavations over the world [4, 11, 12, 18–20]. Each data sample comprises four control factors of the rock burst mentioned above in this research as well as the corresponding rock burst grades, of which the dataset recorded a total of 13 samples labelled 0, 29 samples labelled 1, 55 samples labelled 2, and 15 samples labelled 3. The details of the data samples used to build the Bayes models are shown in Table 1.Table 1
Sample data.
NumberCaseσθσθ/σcσc/σtWETRiskLabel1Yu Zixi Hydropower Station Water Diversion Tunnel, China900.5315.049Strong22Yi Chuan No. 2 Tunnel, China900.4129.737.3Medium13Cavern Highway Tunnel, Iran62.60.3817.539Medium14Takiyami Power Station Headrace tunnel, Japan55.40.3224.119.3Strong25Yutoo Hydropower Station, Norway300.3423.976.6Strong26Yutoo Hydropower Station, Norway48.750.2721.695Strong27Yutoo Hydropower Station, Norway800.4426.875.5Medium18Guanyuk Road Tunnel, Japan890.3828.435Strong29Jin Ping No. 2 Diversion Tunnel, China98.60.8218.463.8Strong210Jin Ping No. 2 Diversion Tunnel, China350.2614.342.9Medium111Tasmanian River Hydropower Station, Australia157.30.5813.186.27Violent312Tasmanian River Hydropower Station, Australia148.40.4517.535.08Medium113Tasmanian River Hydropower Station, Australia132.10.3920.864.63Strong214Tasmanian River Hydropower Station, Australia127.90.2828.93.67Medium115Tasmanian River Hydropower Station, Australia107.50.236.042.29No016Tasmanian River Hydropower Station, Australia96.410.1947.931.87No017Tasmanian River Hydropower Station, Australia167.20.6613.26.83Violent318Tasmanian River Hydropower Station, Australia118.50.2233.752.89Medium119Hui Ze Lead-Zinc Mine, China34.150.634.483.17Medium120Jin Chuan No. 2 Mine, China600.4448.9764.86Medium121Jin Chuan No. 2 Mine, China600.9026.8412.15Medium122Jin Chuan No. 2 Mine, China600.5649.4986.11Medium123Jin Chuan No. 2 Mine, China600.69712.052.85Medium124Jin Chuan No. 2 Mine, China600.40216.043.5Medium125Jin Chuan No. 2 Mine, China600.43913.132.12Medium126Ma Luping lead-zinc mine, China63.80.5824.46.31Strong227Ma Luping lead-zinc mine, China2.60.136.671.39No028Ma Luping lead-zinc mine, China44.40.37245.1Medium129Ma Luping lead-zinc mine, China13.50.4511.22.03Medium130Ma Luping lead-zinc mine, China70.40.6424.46.31Strong231Ma Luping lead-zinc mine, China3.80.196.671.39No032Ma Luping lead-zinc mine, China57.60.48245.1Strong233Ma Luping lead-zinc mine, China19.50.6511.22.03Strong234Ma Luping lead-zinc mine, China81.40.7424.46.31Violent335Ma Luping lead-zinc mine, China4.60.236.671.39No036Ma Luping lead-zinc mine, China73.20.61245.1Strong237Ma Luping lead-zinc mine, China30111.22.03Violent338Birmingham Iron Mine, America15.20.2839.681.92No039Birmingham Iron Mine, America88.90.62710.73.62Violent340Birmingham Iron Mine, America59.820.69711.72.78Strong241Birmingham Iron Mine, America32.30.47910.11.1No042Birmingham Iron Mine, America30.10.3423.976.6Violent343Birmingham Iron Mine, America18.80.1127.227No044Birmingham Iron Mine, America340.2325.257.6Medium145Birmingham Iron Mine, America38.20.7213.591.6No046Birmingham Iron Mine, America11.30.1318.753.6No047Birmingham Iron Mine, America920.3524.588Medium148Birmingham Iron Mine, America62.40.2724.749Violent349Birmingham Iron Mine, America43.40.3218.965.6Violent350Birmingham Iron Mine, America110.121.434.7No051Cheng Chao Iron Mine, China18.70.237.521.5No052Cheng Chao Iron Mine, China28.60.2310.222.5Strong253Cheng Chao Iron Mine, China29.80.2311.524.6Strong254Cheng Chao Iron Mine, China33.60.2214.455.2Strong255Cheng Chao Iron Mine, China26.90.299.83.7Strong256Cheng Chao Iron Mine, China55.90.4420.38.1Violent357Cheng Chao Iron Mine, China59.90.628.261.8Medium158Cheng Chao Iron Mine, China680.6417.517.2Violent359Dong Guashan Copper Mine, China105.50.569.747.27Strong260Dong Guashan Copper Mine, China105.50.6214.055.76Strong261Dong Guashan Copper Mine, China105.50.5511.113.97Strong262Tong Yu Tunnel K21 + 680, China47.560.8116.715Medium163Tong Yu Tunnel K21 + 680, China43.620.5624.416Medium164Da Xiangling Tunnel YK55 + 119, China25.70.4345.91.7No065Da Xiangling Tunnel ZK55 + 154, China26.90.4229.92.4Medium166Da Xiangling Tunnel YK55 + 819, China40.40.5634.31.9Medium167Da Xiangling Tunnel ZK55 + 854, China39.40.628.33.4Strong268Da Xiangling Tunnel YK56 + 080, China38.20.53213.6Strong269Da Xiangling Tunnel YK56 + 109, China45.70.6621.54.1Strong270Da Xiangling Tunnel YK56 + 177, China35.80.5217.84.3Strong271Da Xiangling Tunnel YK56 + 343, China39.40.5725.63.8Strong272Da Xiangling Tunnel ZK56 + 374, China40.60.6125.63.7Strong273Da Xiangling Tunnel YK56 + 421, China390.5629.24.8Strong274Da Xiangling Tunnel YK61 + 305, China57.20.7132.25.5Violent375Da Xiangling Tunnel YK61 + 382, China55.60.4949.54.7Strong276Da Xiangling Tunnel YK61 + 400, China56.90.4645.55.2Strong277Da Xiangling Tunnel ZK61 + 440, China62.10.47555Strong278Da Xiangling Tunnel YK61 + 445, China29.70.2642.93.7Medium179Da Xiangling Tunnel YK61 + 450, China29.10.3136.13.2Medium180Da Xiangling Tunnel YK61 + 493, China27.80.3142.81.8No081Da Xiangling Tunnel YK61 + 827, China30.30.3428.33Medium182Da Xiangling Tunnel YK61 + 382, China55.60.4949.54.7Strong283Da Xiangling Tunnel ZK56 + 451, China41.60.61253.7Strong284Da Xiangling Tunnel ZK56 + 479, China40.10.5531.34.6Strong285Da Xiangling Tunnel ZK61 + 201, China58.20.6932.15.9Violent386Da Xiangling Tunnel ZK61 + 352, China56.80.550.95.2Strong287The Hegula Tunnel, Norway54.40.3421.924.9Strong288West Code Tunnel, Norway47.30.2721.345Strong289Sewage Tunnel, Norway73.50.4221.64.8Strong290La Xiwa Hydropower Station Underground Chamber ,China56.30.3223.788.9Strong291Sheng Qiao Hydropower Station Diversion Tunnel, China29.70.33256.5Strong292Er Tan Hydropower Station Branch tunnel, China88.20.4129.456.9Medium193Jin Ping Hydropower Station Diversion Tunnel, China97.20.8118.754Strong294Yu Zixi Water Hydropower Station Diversion Tunnel, China87.50.5315.148.8Strong295A gold mine roadway, South Africa1080.525.787.1Violent396Hoist underground cavern, South Africa93.90.4725.066.3Strong297Galena Gold Mine, America520.3255.2Strong298Sewage Tunnel, Norway109.40.5831.156.9Strong299Sewage Tunnel, Norway85.70.4624.035.8Medium1100Sewage Tunnel, Norway79.80.4428.574.7Strong2101Long Yangxia Hydropower Station Underground Chamber, China21.60.1231.037.6No0102Tai Pingyi Hydropower Station Underground Chamber, China59.40.3617.748.9Medium1103Li Jiaxia Hydropower Station Underground Chamber, China10.80.0923.535.5No0104Kotoyo Tunnel, Japan790.3626.834.8Strong2105Vietas Water Power Station Diversion Tunnel, Sweden810.4527.275.1Medium1106Raibl metal mine tunnel, Italy103.60.7417.725.2Violent3107Brugge Hydropower Station Underground Chamber, Switzerland35.60.2329.257.9No0108Waterfall Ditch Hydropower Station Underground Chamber, China43.40.3520.55Strong2109Qin Ling Tunnel of Xi Kang Railway, China60.70.5414.26.2Strong2110Bayu Tunnel, China74.20.3921.37.1Violent3111Kamchik Tunnel, Uzbekistan59.60.418.96.9Violent3


### 3.2. Model Training and Verification

This article uses the Scikit-Learn statistical learning platform with Python language to build and train the Naive Bayes model. In order to eliminate the dimensional difference of the sample attributes, the sample attributes are standardized according to formula (7). All attributes are scaled between [0, 1], x is the original value of the attribute, and x′ is the standardized sample attribute value.(7)f:x⟶x′=x−xminxmax−xmin.In order to compare the classification effect of the Naive Bayes models under the four prior distributions, trained models should use unseen samples that have not participated in the training for model verification. In statistical learning, the number of training samples have a great impact on the model’s performance. This article collected 111 learning samples, which are typically small. If a part of samples are retained for model validation, it will further reduce the number of samples participating in model training and harm the model performance. Hence, this paper uses a 10-fold cross-validation method to verify the accuracy of the model. That is, all training samples are randomly divided into 10 disjoint subsets. In each round of training, nine subsets are selected for model training, and the remaining one is used for model validation. The final classification accuracy is the average accuracy of 10 round validations. This strategy ensures that all samples participate in model training and also ensures that the validation samples do not participate in, that is, the independence of verification. Figure1 shows a schematic diagram of the 10-fold cross-validation.Figure 1
Schematic diagram of the 10-fold cross-validation process of the model.Figure2 shows the training process of four Bayesian models with different prior distributions.Figure 2
Model training process.
(a)(b)(c)(d)We can see from Figure2 that as the training samples enter the model, the accuracy of the Gaussian Bayes classifier is steadily increasing. Until all samples are involved in training, the model accuracy stabilizes at 0.4. Complement Naive Bayes classifier has the worst accuracy, only 0.25, which is equal to the probability of random guessing. The model obviously did not acquire any knowledge from the training samples. However, the classification accuracy of Gaussian Bayes classifier, Bernoulli Bayes classifier, and polynomial Bayes classifier is not much different, and they are all approximately 0.4. It is difficult to determine the pros and cons of the three solely based on classification accuracy.


### 3.3. Unequal Cost Classification Results

In Section3.2, the model validation adopted an intuitive metric, that is, the classification accuracy of the model (correctly classified samples/total samples). However, in some special cases, classification accuracy is not a universal measure for classification task. Table 2 gives the confusion matrix of the rock burst classification task in the case of multiclass classification [21]. The elements on the diagonal represent the correct classification, whose misclassification costs are all zero. Other elements represent the corresponding misclassification cost in the case of misclassification. Obviously, this is an asymmetric matrix. For example, Costvn represents that violent rock bursts were misclassified as no rock burst. This may lead to the lack of timely and effective prevention and control measures at the project field, resulting in huge casualties. Costnv represents that it misclassified no rock burst as violent rock burst. This will lead to excessive precautions in areas where there is no danger of rock bursts, resulting in economic and efficiency losses. The first misclassification cost is much higher than the second. Therefore, the risk assessment of the rock burst is a unequal cost classification.Table 2
Generalized confusion matrix for risk evaluation task of rock burst.
Misclassification costPrediction gradeNo rock burstMedium rock burstStrong rock burstViolent rock burstActual gradeNo rock burst0CostNMCostNSCostNVMedium rock burstCostMN0CostMSCostMVStrong rock burstCostSNCostSM0CostSVviolent rock burstCostVNCostVMCostVS0Consider the following case. There are two classifiers A and B; classifier A misclassified all no rock burst cases as violent rock burst cases while classifier B misclassified all violent rock burst cases as no rock burst. The classification accuracy of the two is equal which is zero. However, the misclassification cost for classifier A is lower, and classifier A is better than classifier B. Therefore, the classification accuracy cannot fully reflect the performance of a classifier, especially when used for unequal cost classification tasks.This article introduces receiver operating characteristic (ROC) curve to measure the performance of Bayesian classifiers [22]. The original intention of ROC curve is designed to serve two classification tasks with unequal costs. Combined with the actual needs of this article, the ROC curve is extended to four classifications. We drew the ROC curve of each category and finally used the average curve to reflect the performance of the entire classifier. According to ROC performance criteria, the classifier with the largest area under the curve (AUC) has the strongest performance. Figure 3 reflects the unequal cost classification performance of the three Bayesian classifiers.Figure 3
ROC curves and corresponding AUCs for three Bayesian models.
(a)(b)(c)


## 3.1. Sample Data Collection

Rock burst is caused by the concentration of on-site stress exceeding the energy storage limit of the rock mass. There are two main control factors of a rock burst, the field stress indexes and the rock mass properties. In order to reflect the nature and characteristics of rock burst comprehensively, this article adopted the tangential stress of surrounding rockσθ (MPa), the ratio of tangential stress to uniaxial tensile strength σθ/σt, the ratio of uniaxial compressive strength to uniaxial tensile strength σc/σt, and elastic energy index WET as the attributes of the data sample to evaluate the risk of rock burst. According to the commonly accepted standards for ranking a rock burst grades, a rock burst can be divided into four grades [17]. That is, no rock burst, moderate rock burst, strong rock burst, and severe rock burst. These four levels are used as data labels to participate in the construction of statistical learning models by digitizing. 0 means no rock burst risk, 1 means medium rock burst risk, 2 means strong rock burst risk, and 3 means violent rock burst risk.This article collected 111 rock burst cases from different underground excavations over the world [4, 11, 12, 18–20]. Each data sample comprises four control factors of the rock burst mentioned above in this research as well as the corresponding rock burst grades, of which the dataset recorded a total of 13 samples labelled 0, 29 samples labelled 1, 55 samples labelled 2, and 15 samples labelled 3. The details of the data samples used to build the Bayes models are shown in Table 1.Table 1
Sample data.
NumberCaseσθσθ/σcσc/σtWETRiskLabel1Yu Zixi Hydropower Station Water Diversion Tunnel, China900.5315.049Strong22Yi Chuan No. 2 Tunnel, China900.4129.737.3Medium13Cavern Highway Tunnel, Iran62.60.3817.539Medium14Takiyami Power Station Headrace tunnel, Japan55.40.3224.119.3Strong25Yutoo Hydropower Station, Norway300.3423.976.6Strong26Yutoo Hydropower Station, Norway48.750.2721.695Strong27Yutoo Hydropower Station, Norway800.4426.875.5Medium18Guanyuk Road Tunnel, Japan890.3828.435Strong29Jin Ping No. 2 Diversion Tunnel, China98.60.8218.463.8Strong210Jin Ping No. 2 Diversion Tunnel, China350.2614.342.9Medium111Tasmanian River Hydropower Station, Australia157.30.5813.186.27Violent312Tasmanian River Hydropower Station, Australia148.40.4517.535.08Medium113Tasmanian River Hydropower Station, Australia132.10.3920.864.63Strong214Tasmanian River Hydropower Station, Australia127.90.2828.93.67Medium115Tasmanian River Hydropower Station, Australia107.50.236.042.29No016Tasmanian River Hydropower Station, Australia96.410.1947.931.87No017Tasmanian River Hydropower Station, Australia167.20.6613.26.83Violent318Tasmanian River Hydropower Station, Australia118.50.2233.752.89Medium119Hui Ze Lead-Zinc Mine, China34.150.634.483.17Medium120Jin Chuan No. 2 Mine, China600.4448.9764.86Medium121Jin Chuan No. 2 Mine, China600.9026.8412.15Medium122Jin Chuan No. 2 Mine, China600.5649.4986.11Medium123Jin Chuan No. 2 Mine, China600.69712.052.85Medium124Jin Chuan No. 2 Mine, China600.40216.043.5Medium125Jin Chuan No. 2 Mine, China600.43913.132.12Medium126Ma Luping lead-zinc mine, China63.80.5824.46.31Strong227Ma Luping lead-zinc mine, China2.60.136.671.39No028Ma Luping lead-zinc mine, China44.40.37245.1Medium129Ma Luping lead-zinc mine, China13.50.4511.22.03Medium130Ma Luping lead-zinc mine, China70.40.6424.46.31Strong231Ma Luping lead-zinc mine, China3.80.196.671.39No032Ma Luping lead-zinc mine, China57.60.48245.1Strong233Ma Luping lead-zinc mine, China19.50.6511.22.03Strong234Ma Luping lead-zinc mine, China81.40.7424.46.31Violent335Ma Luping lead-zinc mine, China4.60.236.671.39No036Ma Luping lead-zinc mine, China73.20.61245.1Strong237Ma Luping lead-zinc mine, China30111.22.03Violent338Birmingham Iron Mine, America15.20.2839.681.92No039Birmingham Iron Mine, America88.90.62710.73.62Violent340Birmingham Iron Mine, America59.820.69711.72.78Strong241Birmingham Iron Mine, America32.30.47910.11.1No042Birmingham Iron Mine, America30.10.3423.976.6Violent343Birmingham Iron Mine, America18.80.1127.227No044Birmingham Iron Mine, America340.2325.257.6Medium145Birmingham Iron Mine, America38.20.7213.591.6No046Birmingham Iron Mine, America11.30.1318.753.6No047Birmingham Iron Mine, America920.3524.588Medium148Birmingham Iron Mine, America62.40.2724.749Violent349Birmingham Iron Mine, America43.40.3218.965.6Violent350Birmingham Iron Mine, America110.121.434.7No051Cheng Chao Iron Mine, China18.70.237.521.5No052Cheng Chao Iron Mine, China28.60.2310.222.5Strong253Cheng Chao Iron Mine, China29.80.2311.524.6Strong254Cheng Chao Iron Mine, China33.60.2214.455.2Strong255Cheng Chao Iron Mine, China26.90.299.83.7Strong256Cheng Chao Iron Mine, China55.90.4420.38.1Violent357Cheng Chao Iron Mine, China59.90.628.261.8Medium158Cheng Chao Iron Mine, China680.6417.517.2Violent359Dong Guashan Copper Mine, China105.50.569.747.27Strong260Dong Guashan Copper Mine, China105.50.6214.055.76Strong261Dong Guashan Copper Mine, China105.50.5511.113.97Strong262Tong Yu Tunnel K21 + 680, China47.560.8116.715Medium163Tong Yu Tunnel K21 + 680, China43.620.5624.416Medium164Da Xiangling Tunnel YK55 + 119, China25.70.4345.91.7No065Da Xiangling Tunnel ZK55 + 154, China26.90.4229.92.4Medium166Da Xiangling Tunnel YK55 + 819, China40.40.5634.31.9Medium167Da Xiangling Tunnel ZK55 + 854, China39.40.628.33.4Strong268Da Xiangling Tunnel YK56 + 080, China38.20.53213.6Strong269Da Xiangling Tunnel YK56 + 109, China45.70.6621.54.1Strong270Da Xiangling Tunnel YK56 + 177, China35.80.5217.84.3Strong271Da Xiangling Tunnel YK56 + 343, China39.40.5725.63.8Strong272Da Xiangling Tunnel ZK56 + 374, China40.60.6125.63.7Strong273Da Xiangling Tunnel YK56 + 421, China390.5629.24.8Strong274Da Xiangling Tunnel YK61 + 305, China57.20.7132.25.5Violent375Da Xiangling Tunnel YK61 + 382, China55.60.4949.54.7Strong276Da Xiangling Tunnel YK61 + 400, China56.90.4645.55.2Strong277Da Xiangling Tunnel ZK61 + 440, China62.10.47555Strong278Da Xiangling Tunnel YK61 + 445, China29.70.2642.93.7Medium179Da Xiangling Tunnel YK61 + 450, China29.10.3136.13.2Medium180Da Xiangling Tunnel YK61 + 493, China27.80.3142.81.8No081Da Xiangling Tunnel YK61 + 827, China30.30.3428.33Medium182Da Xiangling Tunnel YK61 + 382, China55.60.4949.54.7Strong283Da Xiangling Tunnel ZK56 + 451, China41.60.61253.7Strong284Da Xiangling Tunnel ZK56 + 479, China40.10.5531.34.6Strong285Da Xiangling Tunnel ZK61 + 201, China58.20.6932.15.9Violent386Da Xiangling Tunnel ZK61 + 352, China56.80.550.95.2Strong287The Hegula Tunnel, Norway54.40.3421.924.9Strong288West Code Tunnel, Norway47.30.2721.345Strong289Sewage Tunnel, Norway73.50.4221.64.8Strong290La Xiwa Hydropower Station Underground Chamber ,China56.30.3223.788.9Strong291Sheng Qiao Hydropower Station Diversion Tunnel, China29.70.33256.5Strong292Er Tan Hydropower Station Branch tunnel, China88.20.4129.456.9Medium193Jin Ping Hydropower Station Diversion Tunnel, China97.20.8118.754Strong294Yu Zixi Water Hydropower Station Diversion Tunnel, China87.50.5315.148.8Strong295A gold mine roadway, South Africa1080.525.787.1Violent396Hoist underground cavern, South Africa93.90.4725.066.3Strong297Galena Gold Mine, America520.3255.2Strong298Sewage Tunnel, Norway109.40.5831.156.9Strong299Sewage Tunnel, Norway85.70.4624.035.8Medium1100Sewage Tunnel, Norway79.80.4428.574.7Strong2101Long Yangxia Hydropower Station Underground Chamber, China21.60.1231.037.6No0102Tai Pingyi Hydropower Station Underground Chamber, China59.40.3617.748.9Medium1103Li Jiaxia Hydropower Station Underground Chamber, China10.80.0923.535.5No0104Kotoyo Tunnel, Japan790.3626.834.8Strong2105Vietas Water Power Station Diversion Tunnel, Sweden810.4527.275.1Medium1106Raibl metal mine tunnel, Italy103.60.7417.725.2Violent3107Brugge Hydropower Station Underground Chamber, Switzerland35.60.2329.257.9No0108Waterfall Ditch Hydropower Station Underground Chamber, China43.40.3520.55Strong2109Qin Ling Tunnel of Xi Kang Railway, China60.70.5414.26.2Strong2110Bayu Tunnel, China74.20.3921.37.1Violent3111Kamchik Tunnel, Uzbekistan59.60.418.96.9Violent3


## 3.2. Model Training and Verification

This article uses the Scikit-Learn statistical learning platform with Python language to build and train the Naive Bayes model. In order to eliminate the dimensional difference of the sample attributes, the sample attributes are standardized according to formula (7). All attributes are scaled between [0, 1], x is the original value of the attribute, and x′ is the standardized sample attribute value.(7)f:x⟶x′=x−xminxmax−xmin.In order to compare the classification effect of the Naive Bayes models under the four prior distributions, trained models should use unseen samples that have not participated in the training for model verification. In statistical learning, the number of training samples have a great impact on the model’s performance. This article collected 111 learning samples, which are typically small. If a part of samples are retained for model validation, it will further reduce the number of samples participating in model training and harm the model performance. Hence, this paper uses a 10-fold cross-validation method to verify the accuracy of the model. That is, all training samples are randomly divided into 10 disjoint subsets. In each round of training, nine subsets are selected for model training, and the remaining one is used for model validation. The final classification accuracy is the average accuracy of 10 round validations. This strategy ensures that all samples participate in model training and also ensures that the validation samples do not participate in, that is, the independence of verification. Figure1 shows a schematic diagram of the 10-fold cross-validation.Figure 1
Schematic diagram of the 10-fold cross-validation process of the model.Figure2 shows the training process of four Bayesian models with different prior distributions.Figure 2
Model training process.
(a)(b)(c)(d)We can see from Figure2 that as the training samples enter the model, the accuracy of the Gaussian Bayes classifier is steadily increasing. Until all samples are involved in training, the model accuracy stabilizes at 0.4. Complement Naive Bayes classifier has the worst accuracy, only 0.25, which is equal to the probability of random guessing. The model obviously did not acquire any knowledge from the training samples. However, the classification accuracy of Gaussian Bayes classifier, Bernoulli Bayes classifier, and polynomial Bayes classifier is not much different, and they are all approximately 0.4. It is difficult to determine the pros and cons of the three solely based on classification accuracy.


## 3.3. Unequal Cost Classification Results

In Section3.2, the model validation adopted an intuitive metric, that is, the classification accuracy of the model (correctly classified samples/total samples). However, in some special cases, classification accuracy is not a universal measure for classification task. Table 2 gives the confusion matrix of the rock burst classification task in the case of multiclass classification [21]. The elements on the diagonal represent the correct classification, whose misclassification costs are all zero. Other elements represent the corresponding misclassification cost in the case of misclassification. Obviously, this is an asymmetric matrix. For example, Costvn represents that violent rock bursts were misclassified as no rock burst. This may lead to the lack of timely and effective prevention and control measures at the project field, resulting in huge casualties. Costnv represents that it misclassified no rock burst as violent rock burst. This will lead to excessive precautions in areas where there is no danger of rock bursts, resulting in economic and efficiency losses. The first misclassification cost is much higher than the second. Therefore, the risk assessment of the rock burst is a unequal cost classification.Table 2
Generalized confusion matrix for risk evaluation task of rock burst.
Misclassification costPrediction gradeNo rock burstMedium rock burstStrong rock burstViolent rock burstActual gradeNo rock burst0CostNMCostNSCostNVMedium rock burstCostMN0CostMSCostMVStrong rock burstCostSNCostSM0CostSVviolent rock burstCostVNCostVMCostVS0Consider the following case. There are two classifiers A and B; classifier A misclassified all no rock burst cases as violent rock burst cases while classifier B misclassified all violent rock burst cases as no rock burst. The classification accuracy of the two is equal which is zero. However, the misclassification cost for classifier A is lower, and classifier A is better than classifier B. Therefore, the classification accuracy cannot fully reflect the performance of a classifier, especially when used for unequal cost classification tasks.This article introduces receiver operating characteristic (ROC) curve to measure the performance of Bayesian classifiers [22]. The original intention of ROC curve is designed to serve two classification tasks with unequal costs. Combined with the actual needs of this article, the ROC curve is extended to four classifications. We drew the ROC curve of each category and finally used the average curve to reflect the performance of the entire classifier. According to ROC performance criteria, the classifier with the largest area under the curve (AUC) has the strongest performance. Figure 3 reflects the unequal cost classification performance of the three Bayesian classifiers.Figure 3
ROC curves and corresponding AUCs for three Bayesian models.
(a)(b)(c)


## 4. Engineering Case Analysis

The Bayesian classifier based on Gaussian distribution is found to have the best performance among classifiers with four prior distributions, that is, the strongest rock burst risk assessment ability. The trained Gaussian classifier was used to evaluate the rock burst hazard of kimberlite in a diamond mine in Northern Canada.The diamond mine is located in the North Slave Lake region, 250 kilometers north of Yellowknife, the capital of Yukon Province, Canada. The mining area includes three kimberlite pipes with the diameters more than 500 m, all of which are blasted underground. Kimberlite is a kind of igneous rock that usually contains minerals such as diamonds. The lithology of the kimberlite is hard, which is prone to burst. Figure4 shows an aerial view of this diamond mine and a mining kimberlite stope.Figure 4
An aerial view of the diamond mine and a schematic diagram of mining kimberlite pillars [23].
(a)(b)In order to obtain the corresponding risk assessment attributes, laboratory tests of rock properties and project site investigations are implemented. Twelve locations in two kimberlite pipes are selected to investigate the tangential stresses of the surrounding rock mass. The rock samples near the survey locations are extracted to test the rock burst attribute parameters. The distribution of the survey area in kimberlite pipes and the schematic diagram of the lab tests are shown in Figure5.Figure 5
Schematic test of mechanical parameters and field sampling points of kimberlite.
(a)(b)Table3 shows the rock burst attributes of 12 kimberlite samples. The Naive Bayes classifier with Gaussian prior distribution is used to evaluate the rock burst hazard, and the corresponding hazard evaluation results are obtained.Table 3
Sample attributes and corresponding rock burst risk assessment results.
Sampling pointσθ (MPa)σθ/σcσc/σtWETRisk assessment results118.1670.3731.43.3No2210.3518.91.7Medium331.160.3821.22.3Strong446.3760.6225.13.2Medium548.640.6418.62.5Medium622.920.4401.5Medium799.0880.8830.15.2Medium835.1560.4425.62.5Strong915.840.3222.92.8No1013.020.228.51.2Medium1121.120.424.22.3Medium1229.1210.5117.12.2StrongThe hazard evaluation results show that seven out of the 12 sampling areas had medium rock burst hazard, three had strong rock burst hazard, and the remaining two areas had no rock burst hazard. The evaluation results matched mining logs which recorded two rock burst cases in this diamond mine in 2017. The first rock burst case occurred in the haul load at the N9750 level of the A154 kimberlite pipe, near the No. 7 sampling location. The second rock burst case occurred in the haul road at the N9850 level of A418 pipe, near the No. 12 sampling location. The mining log recorded the first rock burst case situation in which the rock mass was peeling off from the roof of the roadway and accompanied by a small-scale ejection, damaging the bolt and wire mesh support of the roadway. The length of the rock burst area is about 3 m, and no rock burst sound signal was recorded at the location before the occurrence. The second rock burst was even more severe. The rock mass ejected from the roof and the roadway support was completely damaged. The length of the rock burst area reached 12 m. According to the rock burst evaluation rules, the first rock burst can be determined as a medium rock burst, and the second one can be determined as a strong rock burst. The rock burst cases verify the effectiveness of our proposed Gaussian distribution Bayesian model.


## 5. Conclusion and Discussion

(1)
This paper uses a statistical learning model to obtain information from the collected data of 111 rock burst cases and trains Naive Bayes classifiers based on different prior distributions. By comparing the classification performance of the four classifiers with different prior distributions, we finally determined that the Naive Bayesian model based on Gaussian distribution has the strongest classification ability. The trained model was used to determine the rock burst hazard of kimberlite in a diamond mine in Canada. The evaluation results obtained from the proposed model match the rock burst observations on the site, which verifies the applicability of the model.(2)
The core part of constructing a Bayesian classifier is to determine the distribution type of training samples. The four control factors of rock burst selected in this paper represent the two necessary conditions inducing a rock burst. The interactions over these control factors are complicated, which is difficult to determine the data distribution of attributes in advance solely through the rock burst mechanism. This paper compares four commonly used prior distributions through model validation and finds that the Gaussian Bayesian classifier works the best. Compared with the other three prior distributions, the selected rock burst attributes follow the general Gaussian distribution.(3)
Although the performance of the model is in line with expectations, it can be seen that the performance of the model still needs to be improved. For example, the general classification accuracy of the model in Figure2 is slightly higher than 0.4, and there is still much room for improvement. Another shortcoming is the “partial discipline” of the model. By analyzing the ROC curve, it can be seen that the model has better performances on the evaluation of no and strong rock burst hazards than on moderate and severe rock burst hazards.(4)
Taking the characteristics of the proposed Bayesian model into account, the following measures are expected to improve model performance. The straightforward one is to more accurately derive the data distribution of the control factors of the rock burst. Combining the physical characteristics of the rock burst, we can explore the mathematical laws of the control factors, so as to establish the prior distribution of the data more accurately and improve the performance of the Bayesian model. The second is to collect more learning samples. If it is difficult to accurately determine the prior distribution for data, collecting more samples can endow the model the better prediction performance. The first method allows the Bayesian model not to be limited to several commonly used distributions. Building accurate prior distribution will cause a qualitative improvement in model performance. However, when the first method fails to realize, the second one can make up for the complement of model performance. However, due to the fixed prior distribution of the model, the second method can only achieve limited performance improvement of the model.(5)
This article discards the conventional strategies based on the mechanism for rock burst hazard evaluation. This paper starts from the data, uses the statistical learning model to automatically obtain data information, and matches the knowledge obtained by the new input to complete the task of rock burst risk assessment. The ideas in this article can provide references for the rock burst prevention and prediction of mine dynamic disasters.



---

*Source: 5576480-2021-04-05.xml*