李 盛，郭民之，康文倩

(云南師范大學(xué)數(shù)學(xué)學(xué)院，昆明650500)

微震災(zāi)害是最難探測和預(yù)測的自然災(zāi)害之一，越來越多的先進(jìn)的微震預(yù)警系統(tǒng)使我們更好的了解巖體運(yùn)動過程和微震災(zāi)害預(yù)測方法的定義，然而到目前為止很多方法的準(zhǔn)確性不盡如人意.在很多高能量(大于10 000 J)微震事件中，由于微震過程的復(fù)雜性和不均衡性導(dǎo)致統(tǒng)計技術(shù)不足以預(yù)測微震災(zāi)害.因此，有必要尋找更好的方法預(yù)測微震災(zāi)害，使用機(jī)器學(xué)習(xí)方法就是一種途徑.Lesniak，Isakow［1］用數(shù)據(jù)聚類技術(shù)評估微震災(zāi)害和Kabiesz［2］使用人工神經(jīng)網(wǎng)絡(luò)對微震災(zāi)害預(yù)測都有一定效果.在許多文獻(xiàn)中，提到的方法以“危險”和“無危險”兩種狀態(tài)形式呈現(xiàn)，正數(shù)(危險狀態(tài))和負(fù)數(shù)(無危險狀態(tài))的不平衡分布是微震災(zāi)害預(yù)測中的一個嚴(yán)重問題，當(dāng)前使用的方法不足以使預(yù)測達(dá)到良好的敏感性和特異性.Kijko［3］使用非參數(shù)方法對礦山微震災(zāi)害進(jìn)行分析，Sikora［4］用歸納和修剪的分類規(guī)則預(yù)測煤礦微震的危害，這兩種方法所取得的效果不是太理想.Bukowska［5］在其論文中提出在發(fā)生微震時能量大于10 000 J情況下，在眾多因數(shù)中，有一些因數(shù)影響微震災(zāi)害的發(fā)生.微震預(yù)測可以用不同的方式來定義，但主要目的都是預(yù)測微震活動可導(dǎo)致巖爆的精確日期和時間.

煤礦礦震危險性預(yù)測的方法主要有:模糊聚類法和神經(jīng)網(wǎng)絡(luò)算法，但模糊聚類法在確定分類指標(biāo)權(quán)值時，存在一定的主觀性，而傳統(tǒng)的神經(jīng)網(wǎng)絡(luò)算法都是基于大樣本數(shù)據(jù)的，推廣性較差.基于隨機(jī)森林分類原理的預(yù)測方法有不需要對數(shù)據(jù)預(yù)處理，能有效的解決不平衡、高維分類問題，能較好的容忍噪聲并且不會過擬合，分類結(jié)果穩(wěn)定等優(yōu)點(diǎn).礦震危險性預(yù)測是一個非線性、高維的多類模式識別問題，隨機(jī)森林方法更適合解決這類模式識別問題.因此，本文提出了一種基于隨機(jī)森林法的煤礦微震危險性預(yù)測方法.

1 隨機(jī)森林法原理及性質(zhì)

1.1 原理

隨機(jī)森林法(random forests)是 Breiman［6］于2001年提出的一種新的組合分類器算法.隨機(jī)森林分類是由很多決策樹分類模型{h(x，θk)，k=1，2，…}組成的組合分類模型，參數(shù)集{θk}是獨(dú)立同分布的隨機(jī)向量，在給定自變量x下，每個決策樹分類模型都由一票投票權(quán)來選擇最優(yōu)的分類結(jié)果.利用Bootstrap抽樣從原始訓(xùn)練集抽取k個樣本，每個樣本的樣本容量均與原始訓(xùn)練集一樣，對k個樣本分別建立k個決策樹模型，得到k種分類結(jié)果.根據(jù)k種分類結(jié)果對每個記錄進(jìn)行投票表決決定其最終分類，如圖1.

隨機(jī)森林法通過構(gòu)造不同的訓(xùn)練集增加分類模型的差異，以提高組合分類模型的外推預(yù)測能力.通過k輪訓(xùn)練，得到分類模型序列{h1(x)，h2(x)，…，hk(x)}，再用它們構(gòu)成多分類模型系統(tǒng)，該系統(tǒng)的最終分類結(jié)果采用簡單多數(shù)投票法，最終的分類決策為:

其中，h(x)表示組合分類模型，hi是單個決策樹分類模型，y表示輸出變量，I(·)為示性函數(shù).公式(1)表明使用多數(shù)投票決策的方式來確定最終的分類.

1.2 隨機(jī)森林算法

在分類回歸樹(CART)算法中，每個內(nèi)部節(jié)點(diǎn)都是原始數(shù)據(jù)集的子集，根節(jié)點(diǎn)包含了所有的原始數(shù)據(jù).在每個內(nèi)部節(jié)點(diǎn)處，從所有屬性中找出最好的分裂方式進(jìn)行分裂，再對后續(xù)節(jié)點(diǎn)依次進(jìn)行分裂，直到葉節(jié)點(diǎn)，最后通過剪枝使測試誤差最小.而隨機(jī)森林法與算法有所不同，單棵樹的生長可歸納為以下3點(diǎn):

1)用Bagging方法形成新的訓(xùn)練集:在樣本數(shù)為N的原始訓(xùn)練集中，有放回地隨機(jī)選取N個樣本形成一個新的訓(xùn)練集，以此生成一棵分類樹;

2)隨機(jī)選取特征對分類回歸樹的內(nèi)部節(jié)點(diǎn)進(jìn)行分裂:設(shè)共有 M個特征，取任意正整數(shù) m，且m?M.對于每個內(nèi)部節(jié)點(diǎn)，從M個特征中隨機(jī)抽取m個特征作為候選特征，選擇這m個特征上最好的分裂方式對節(jié)點(diǎn)進(jìn)行分裂;

3)每棵樹自由生長，不進(jìn)行剪枝.

1.3 隨機(jī)森林的泛化誤差

泛化誤差的大小可以衡量分類器的性能，泛化誤差越小，則分類器的性能越好，反之則性能越差.給定分類器 h1，h2(x)，…，hk(x)，對于輸入變量 x和輸出變量y，定義邊緣函數(shù)為:其中I(·)為示性函數(shù)，avk(·)為取平均數(shù).邊緣函數(shù)衡量給定的分類器集合將輸入變量x分到正確類別的平均票數(shù)，與分到其他類別的平均票數(shù)的最小差值，因此邊緣函數(shù)值越大，分類的可信度就越高，分類器的泛化誤差表示為 PE*=Px，y(mg(x，y)＜0).

當(dāng)隨機(jī)森林中樹的棵數(shù)足夠多時，根據(jù)大數(shù)定律我們可以得到:隨著森林中樹的棵數(shù)增加，對于序列，幾乎處處收斂于

公式(3)表明當(dāng)隨機(jī)森林中樹的棵數(shù)增加時，不會發(fā)生過擬合，并得到一個有限的泛化誤差值.隨機(jī)森林的泛化誤差上界由給出，其中，s表示單棵樹的分類性能，ˉρ表示樹與樹之間的相關(guān)性.由此可知，單棵樹的分類性能越好，樹與樹之間的相關(guān)度越低，隨機(jī)森林的泛化誤差上界越小.

2 標(biāo)準(zhǔn)化均方誤差(NMSE)

采用五折交叉驗(yàn)證的方法來判斷各種機(jī)器學(xué)習(xí)方法結(jié)果的可靠性.計算中通過隨機(jī)建立的5個訓(xùn)練集建立5個模型，對訓(xùn)練集和測試集分別得到5個標(biāo)準(zhǔn)化均方誤差(NMSE)，再得出5次平均的NMSE.令為因變量均值，為從訓(xùn)練集得到的模型對一個數(shù)據(jù)集(可能是訓(xùn)練集本身也可能是測試集)的預(yù)測值，這里的NMSE定義為:

顯然，如果什么模型都不用，只用均值來做預(yù)測，那么 NMSE等于1.所以，如果在回歸時得到NMSE大于1，這個模型就很糟糕了，還不如沒有模型.僅僅對于訓(xùn)練集來說，其NMSE等于1－R2，R2為回歸系數(shù).但是對于測試集來說，其NMSE與測試集的R2沒有什么關(guān)系.

3 算例

3.1 樣本選取

本文數(shù)據(jù)來自UCI機(jī)器學(xué)習(xí)數(shù)據(jù)庫中的seismic－bumps數(shù)據(jù)集，數(shù)據(jù)來自位于波蘭的采用長壁開采法的煤礦每8 h監(jiān)測一次的實(shí)時數(shù)據(jù).數(shù)據(jù)集描述在高能量(大于104J)情況下關(guān)于煤礦開采微震災(zāi)害的預(yù)測問題，數(shù)據(jù)集包括2 584個樣本，9個自變量，其中包括2種微震危害評估方法［7］x1、x2(a表示不危險，b表示低危險，c表示高危險，d表示危險狀態(tài))，x3(W表示繼續(xù)工作，N表示準(zhǔn)備轉(zhuǎn)移)，x4為由微震檢波器檢測到的能量，x5為檢測到的脈沖，x6為當(dāng)前記錄能量與先前8次記錄能量的平均的偏差，x7為當(dāng)前記錄脈沖與先前8次記錄脈沖的平均的偏差，x8為上一次轉(zhuǎn)移登記能量值的總和，x9為上一次轉(zhuǎn)移登記的最大能量值.一個因變量(y)，即危險級別，1表示下一次轉(zhuǎn)移是高能量顫動(危險狀態(tài))，0表示下一次轉(zhuǎn)移不是高能量顫動(無危險狀態(tài))，其中有167個正數(shù)1.使用R軟件對數(shù)據(jù)集進(jìn)行分析［8－9］.

表1 數(shù)據(jù)集變量

3.2 數(shù)據(jù)分析

數(shù)據(jù)集描述在高能量(大于104J)情況下關(guān)于煤礦開采微震災(zāi)害的預(yù)測問題，通過已有數(shù)據(jù)自身特點(diǎn)的分類分析，預(yù)測下一次高能量顫動采取的策略:“1”表示將有高能量震動(危險狀態(tài))，“0”表示接下來無高能量震動(無危險狀態(tài))，其中有167個正數(shù)1，占總數(shù)據(jù)的6.5%.數(shù)據(jù)集中所測能量大于104J的樣本有2 035個，占總數(shù)據(jù)的78.8%，在高能量顫動發(fā)生時采取轉(zhuǎn)移的樣本占高能量樣本的8.2%，數(shù)據(jù)的部分特征見表2.

表2 微震能量大于104J的數(shù)據(jù)總結(jié)

上表可以看出是否轉(zhuǎn)移與高能量顫動的發(fā)生是密切相關(guān)，但不是只要是高能量發(fā)生就必須要轉(zhuǎn)移.因?yàn)槲⒄疬^程既包含線性的關(guān)系，也包含了復(fù)雜的非線性關(guān)系，所以分析數(shù)據(jù)集各變量之間是否具有線性關(guān)系.由圖2可以看出因變量與各自變量間不存在明顯的線性相關(guān)性，線性模型在分析多樣本、高維度的高能量煤礦礦山微震數(shù)據(jù)效果不好.

下面對危險評估方法進(jìn)行分析，選取危險級別為“1”的樣本集和緊接危險級為“1”后的一個樣本構(gòu)成的樣本集為分析對象，通過對微震危害評估方法分析，發(fā)現(xiàn)seismic評估方法對緊接危險級“1”后的一個樣本評估比seismoacoustic評估方法好點(diǎn)，但都不能準(zhǔn)確地預(yù)測微震災(zāi)害.

表3 危險評估方法與危險級別的比例分析

用傳統(tǒng)的線性模型很難分析多變量、不平衡分布的微震數(shù)據(jù)，故提出機(jī)器學(xué)習(xí)方法，機(jī)器學(xué)習(xí)方法對數(shù)據(jù)沒有任何假定，產(chǎn)生的結(jié)果用交叉驗(yàn)證的方法來判斷，脫離了假定分布?假設(shè)檢驗(yàn)?p值的經(jīng)典統(tǒng)計過程。這種基于算法或程序的模型預(yù)測效果相當(dāng)好，而且交叉驗(yàn)證的結(jié)果也容易被廣大實(shí)際工作者所理解和接受.

表4給出了隨機(jī)森林方法對數(shù)據(jù)集的NMES分析結(jié)果，并與決策樹、Bagging算法、支持向量機(jī)(SVM)和最近鄰法所得到的NMSE結(jié)果進(jìn)行比較，發(fā)現(xiàn)這幾種機(jī)器學(xué)習(xí)方法對數(shù)據(jù)集分析效果都較好，其中隨機(jī)森林方法對數(shù)據(jù)集構(gòu)造的訓(xùn)練集和測試集的為9.67×10－5和0.067 33，是上述方法中結(jié)果最理想的，說明隨機(jī)森林方法在處理高能量煤礦礦山微震數(shù)據(jù)時能很好的控制誤差.

表4 幾種機(jī)器學(xué)習(xí)方法對數(shù)據(jù)集的NMSE分析

表5 幾種機(jī)器學(xué)習(xí)方法對數(shù)據(jù)集的預(yù)測精度分析

對數(shù)據(jù)集做預(yù)測精度分析，先把數(shù)據(jù)分成E≥105J數(shù)據(jù)、E≥104J數(shù)據(jù)、E＜104J數(shù)據(jù)和整個數(shù)據(jù)4個樣本集，用表5中的機(jī)器學(xué)習(xí)方法對各個樣本集做預(yù)測精度分析，上述方法對危險級為“0”的樣本的預(yù)測都較理想，預(yù)測精度都在90%以上，對危險級為“1”的樣本的預(yù)測較差，可能由其樣本數(shù)在總樣本數(shù)中占的比例較小的緣故.其中，E≥105J的樣本的預(yù)測結(jié)果較其他樣本集的結(jié)果好，說明高能量的檢測與微震的發(fā)生是緊密相關(guān)的.隨機(jī)森林方法對數(shù)據(jù)的預(yù)測分析效果最理想.

4 結(jié)語

在煤礦開采中微震引發(fā)的一系列監(jiān)測數(shù)據(jù)可以描述為時間系列，其中既包含線性的關(guān)系，也包含了復(fù)雜的非線性關(guān)系，微震過程的復(fù)雜性和不均衡性

導(dǎo)致線性模型不足以預(yù)測微震災(zāi)害.本文注重采用機(jī)器學(xué)習(xí)方法分析在高能量(大于)情況下關(guān)于煤礦開采微震災(zāi)害的預(yù)測問題，發(fā)現(xiàn)隨機(jī)森林法、決策樹、Bagging算法、隨機(jī)森林、支持向量機(jī)和最近鄰方法對處理高能量微震數(shù)據(jù)都具有較好的誤差容忍性，其中隨機(jī)森林法的五折交叉驗(yàn)證的值都較低，能很好的控制誤差，在對高能量煤礦礦山微震預(yù)測分析時，隨機(jī)森林法效果最理想.本文的不足之處是沒能結(jié)合國內(nèi)的煤礦礦山微震數(shù)據(jù)進(jìn)行比較分析.

［1］LE S＇NIAK A，ISAKOW Z.Space-time clustering of seismic events and hazard assessment in the Zabrze-Bielszowice coal mine，Poland［J］.International Journal of Rock Mechanics and Mining Sciences，2009，46(5):918－928.

［2］KABIESZ J.Effect of the form of data on the quality ofmine tremors hazard forecasting using neural networks［J］.Geotechnical and Geological Engineering，2006，24(5):1131－1147.

［3］KIJKO A，LASOCKI S，GRAHAM G.Non－parametric seismic hazard in mines［J］.Pure and Applied Geophysics，2001，158(9/10):1655 －1675.

［4］SIKORA M.Induction and pruning of classification rules for prediction of microseismic hazards in coal mines［J］.Expert Systems with Applications，2011，38(6):6748－6758.

［5］BUKOWSKA M.The probability of rockburst occurrence in the Upper Silesian Coal Basin area dependent on natural mining conditions［J］.Journal of Mining Science，2006，42(6):570－577.

［6］BREIMAN L.Random forests［J］.Machine Learning，2001，45(1):5－32.

［7］BARANSKI A，DRZEWIECKI J，KABIESZ J，et al.Rules of application of the comprehensive and detailed rockburst hazard assessment methods in hard － coal mines［J］.Expert Systems with Applications:An International Journal，2011，38(6):6748 －6758.

［8］吳喜之.復(fù)制數(shù)據(jù)統(tǒng)計方法——基于R的應(yīng)用［M］.北京:中國人民大學(xué)出版社，2012.

［9］吳喜之.統(tǒng)計學(xué):從數(shù)據(jù)到結(jié)論［M］.北京:中國統(tǒng)計出版社，2013.