基于主成分分析和改進(jìn)Bayes判別的巖爆等級預(yù)測

劉 劍，周宗紅，劉 軍，洪貞群

( 1. 昆明理工大學(xué) 國土資源工程學(xué)院，云南 昆明 650093；2. 元陽縣華西黃金有限公司，云南 紅河 661400 )

巖爆是深埋地下工程的地質(zhì)災(zāi)害之一，對人員、設(shè)備和國家財產(chǎn)造成巨大的威脅。近年來，我國深部地下工程逐年增多，巖爆災(zāi)害也呈頻發(fā)趨勢[1-3]。

目前，相關(guān)學(xué)者從不同的方向?qū)r爆進(jìn)行預(yù)測研究，提出了許多對巖爆烈度的預(yù)測方法。蔡美峰[4]等基于開挖擾動能量積聚理論，對三山島金礦未來深部開采過程中可能誘發(fā)巖爆的地點和等級進(jìn)行預(yù)測；陳炳瑞[5]等研發(fā)了傳感—采集—傳輸一體化集成、32位A/D與元器件聯(lián)合降噪和微震信號遞歸神經(jīng)網(wǎng)絡(luò)識別等技術(shù)，使微震監(jiān)測技術(shù)朝著巖石工程災(zāi)害自動監(jiān)測、分析與智能預(yù)警方向快速發(fā)展；李克鋼[6]等提出一種基于改進(jìn)綜合賦權(quán)的巖爆傾向性評價方法，選取應(yīng)力條件和圍巖對應(yīng)的15個因素作為巖爆傾向性指標(biāo)，判斷巖爆傾向性等級與可靠性；李寧[7]等建立了粗糙集理論和粒子群支持向量機(jī)( RS-PSOSVM )的巖爆預(yù)測模型，并將該模型用于大相嶺隧道巖爆的預(yù)測；吳順川[8]等基于PCAPNN原理建立巖爆預(yù)測模型，在保證預(yù)測精度的前提下提高收斂速度；湯志立[9]等基于9種經(jīng)典機(jī)器學(xué)習(xí)算法建立巖爆預(yù)測模型，并應(yīng)用于西藏多雄拉隧道的巖爆預(yù)測；謝學(xué)斌[10]等基于CRITIC-XGB算法建立巖爆傾向性預(yù)測模型，為巖爆傾向性等級預(yù)測提供一種新的可靠方法；周科平[11]等采用熵權(quán)法和云模型判定巖爆等級，建立熵權(quán)-正態(tài)云模型并應(yīng)用于公路隧道工程中；王佳信[12]等基于PCA-DDA建立沖擊地壓的判別模型，對復(fù)雜的地質(zhì)災(zāi)害進(jìn)行預(yù)測。

雖然相關(guān)學(xué)者對于巖爆進(jìn)行了大量研究并取得了一定成果，但是巖爆的機(jī)理十分復(fù)雜，采用數(shù)學(xué)或者力學(xué)的方法很難達(dá)到準(zhǔn)確的預(yù)測精度。

近年來Bayes優(yōu)異分類性能受到廣泛關(guān)注，付玉華[13]等基于Bayes判別選取4個指標(biāo)對巖爆進(jìn)行了預(yù)測；宮鳳強(qiáng)[14]等提出了地下工程巖爆發(fā)生及烈度分級預(yù)測的Bayes判別分析方法。上述巖爆預(yù)測理論都取得了較好的預(yù)測結(jié)果，但考慮到巖爆樣本所屬類別分布不均衡，導(dǎo)致了Bayes判別在決策邊界附近出現(xiàn)誤判，降低了判別的準(zhǔn)確率。筆者提出了基于閾值改進(jìn)的Bayes判別模型，在考慮到數(shù)據(jù)樣本的先驗概率和錯誤損失的前提下，不斷更新先驗概率的值，進(jìn)而求得后驗概率。設(shè)置后驗概率的閾值來消除決策邊界附近分類不確定性的風(fēng)險。分類邊界處通過最大后驗概率判別存在較大誤判風(fēng)險，基于閾值改進(jìn)的Bayes判別可以消除誤判風(fēng)險，提高判別的正判率。

1 主成分分析和閾值改進(jìn)Bayes理論

1.1 主成分分析法理論

主成分分析是一種通過降維技術(shù)把多個變量轉(zhuǎn)化成幾個主成分的多元統(tǒng)計學(xué)方法，得到的主成分能夠反映原始變量的大部分信息，通常表示為原始變量的組合。為了達(dá)到主成分所包含的信息不重疊，要求各成分之間互不相關(guān)[15]。設(shè)有p個變量第i次試 驗 中的 取值 為得到原始矩陣為

步驟1：計算xij*標(biāo)準(zhǔn)差。

主成分分析的結(jié)果易受量綱影響，為了消除量綱對變量的影響，先把各變量的單位標(biāo)準(zhǔn)化得到矩陣x*。xij*為矩陣x*的一個數(shù)，計算公式為

式中，為變量xj的觀測值的平均值；sjj為變量xj的觀測值的方差；為標(biāo)準(zhǔn)差。

步驟2：計算相關(guān)系數(shù)矩陣R。

和rpo為對稱矩陣，指標(biāo)間相關(guān)系數(shù)rop計算公式為

步驟3：確定主成分個數(shù)。

通常特征值大于1，累計方差貢獻(xiàn)率大于80%，特征值記為λ1,λ2,…,λn，方差貢獻(xiàn)率vs和累計方差貢獻(xiàn)率vt計算公式分別為

步驟4：計算提取主成分的對應(yīng)得分。

其中主成分系數(shù)矩陣為U=(p1,p2,…,pn)，若提取k個主成分得分

1.2 改進(jìn)的Bayes判別理論

k組p維數(shù)據(jù)沿某個方向投影，投影需滿足組與組之間盡可能分開的條件[16]。

先驗概率ip計算公式為

式中，nω為類別屬于ω的樣本數(shù)量為樣本總體。

式中，iμ為第i個總體的均值；iΣ為協(xié)方差矩陣。

根據(jù)Bayes理論，樣品X來自總體Gi的后驗概率，計算公式為

在不考慮誤判代價的情況下，判別規(guī)則為

Bayes判別分析模型是根據(jù)已知分類標(biāo)記的訓(xùn)練樣本集構(gòu)建一個分類器，然后根據(jù)構(gòu)建好的分類器對已知屬性和未知標(biāo)記的數(shù)據(jù)進(jìn)行分類。根據(jù)最大后驗概率原則，確定最終的分類。Bayes判別分析模型在兩個類別決策邊界附近，屬于該類的決策有很大的不確定性，易出現(xiàn)誤判現(xiàn)象。本文通過設(shè)置警戒值，在分類邊界處尋找一個閾值，作為Bayes判別分析模型有效性的標(biāo)準(zhǔn)，當(dāng)后驗概率大于閾值時，可以認(rèn)為最大后驗概率分類有效；當(dāng)后驗概率小于閾值時拒絕分類，表明該樣本后驗概率位于分類邊界附近，分類具有較大的風(fēng)險。

考慮到Bayes判別先驗概率差距較大時，對后驗概率分類造成較大影響，根據(jù)先驗概率差距大造成誤判的現(xiàn)象，提出一種按先驗概率比例來確定閾值的方法，先驗概率比例系數(shù)kj的計算公式為

式中，i=1，2，3，4；mp和np分別為1級、2級、3級和4級先驗概率；maxmp為最大先驗概率；minnp為最小先驗概率。

當(dāng)樣本數(shù)據(jù)各類別為相等數(shù)量時，各類別先驗概率相等，kj=1。

后驗概率比例系數(shù)計算公式為

式中，Ω為閾值，即后驗概率臨界值。

kj＞ki拒絕判別；k j＜ki接受判別；當(dāng)時，拒絕判別；當(dāng)時，接受判別。

當(dāng)設(shè)Ω=1時任何決策都會被拒絕，當(dāng)Ω過小( 樣本個數(shù)的倒數(shù) )時，所有樣本的決策都不會被拒絕，無法提高模型的泛化能力。

2 基于閾值改進(jìn)的Bayes判別巖爆等級綜合預(yù)測模型

2.1 閾值改進(jìn)Bayes判別評價步驟

閾值改進(jìn)Bayes判別預(yù)測模型是一個綜合評判巖爆等級的模型，通過主成分分析消除信息重疊對Bayes判別分析的影響，建立改進(jìn)的Bayes判別模型。通過閾值判斷后驗概率的有效性，相較于傳統(tǒng)的Bayes判別分析，拒絕了巖爆決策邊界的附近分類( 分類決策有很大不確定性 )，提高了模型的正判率?；陂撝蹈倪M(jìn)的Bayes判別綜合預(yù)測模型流程如圖1所示。

圖1 閾值改進(jìn)Bayes判別綜合預(yù)測模型流程Fig. 1 Flow chart of Bayes discriminant synthesis prediction model is improved by threshold value

2.2 建立巖爆傾向性指標(biāo)體系

巖爆發(fā)生的機(jī)理十分復(fù)雜。影響因素較多，主要包括巖體條件、應(yīng)力水平、開挖方法、工程埋深及巖體的儲能情況等。相關(guān)學(xué)者[2,17]研究發(fā)現(xiàn)，巖爆易發(fā)生在應(yīng)力集中程度較高、高儲能和完整的脆性硬巖中，主要破壞形式為張拉破壞、伴隨剪切破壞( 抗剪記錄較少，無法進(jìn)行分析 )。其中圍巖最大切應(yīng)力( Maximum Tangential Stress，MTS )和應(yīng)力集中系數(shù)( Stress Concentration Factor，SCF )能夠反映圍巖應(yīng)力集中程度；單軸抗拉強(qiáng)度( Uniaxial Tensile Strength，UTS )能夠反映巖爆斷面破壞形式；單軸抗壓強(qiáng)度( Uniaxial Compressive Strength，UCS )能夠反映巖石的堅硬程度和完整性；脆性系數(shù)( Brittleness Index，BI )能夠反映巖石的脆性大小，彈性能量指數(shù)( Elastic Energy Index，EEI )能夠反映巖體儲存和釋放能量的性能。故本文選取以下6個指標(biāo)：圍巖最大切應(yīng)力( MTS )、單軸抗拉強(qiáng)度( UTS )、單軸抗壓強(qiáng)度( UCS )、應(yīng)力集中系數(shù)( SCF )、脆性系數(shù)( BI )和彈性能量指數(shù)( EEI )等作為巖爆等級預(yù)測的評價指標(biāo)，并按巖爆烈度由弱到強(qiáng)分為4個等級[17]，分別為無巖爆( 1級 )、弱巖爆( 2級 )、中等巖爆( 3級 )和強(qiáng)巖爆( 4級 )。依據(jù)所選指標(biāo)整理出典型巖爆工程案例數(shù)據(jù)[18-24]，見表1。

表1 實測巖爆實例數(shù)據(jù)Table 1 Actual data of rockburst

2.3 主成分分析法結(jié)果分析

為了消除指標(biāo)量綱對模型判別精度的影響，對表1所列44組數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理。各指標(biāo)的相關(guān)系數(shù)見表2。

表2 各指標(biāo)相關(guān)系數(shù)Table 2 Correlation coefficient of each index

由表2可知，部分指標(biāo)的相關(guān)系數(shù)大于0.5，從而使得指標(biāo)之間信息重疊，有較為明顯的相關(guān)性，所以采取主成分分析法消除指標(biāo)之間的相關(guān)性，消除非獨立性對Bayes判別模型影響，提高巖爆預(yù)測模型的準(zhǔn)確性。

將標(biāo)準(zhǔn)化后的數(shù)據(jù)進(jìn)行主成分分析，其累計方差貢獻(xiàn)率見表3。

表3 累計方差貢獻(xiàn)率Table 3 Cumulative variance contribution rate

由表3可知，前3個主成分累計方差貢獻(xiàn)率為86.719%。國內(nèi)外學(xué)者認(rèn)為累計方差貢獻(xiàn)率大于80%和特征值不小于1( 圖2 )的主成分指標(biāo)能充分取代其他指標(biāo)，可以將計算得出的3個主成分Y1，Y2和Y3作為巖爆綜合預(yù)測模型的指標(biāo)。

圖2 主成分分析碎石Fig. 2 Principal component analysis of lithotripsy

主成分Y1，Y2和Y3與6個巖爆指標(biāo)關(guān)系式為

由式( 15 )可以看出，3個主成分可以充分代表圍巖最大切應(yīng)力( MTS )、單軸抗壓強(qiáng)度( UCS )、單軸抗拉強(qiáng)度( UTS )、應(yīng)力集中系數(shù)( SCF )、脆性系數(shù)( BI )和彈性能量指數(shù)( EEI )6個指標(biāo)的絕大部分信息。

2.4 改進(jìn)Bayes判別分析

將主成分分析的巖爆指標(biāo)數(shù)據(jù)代入Bayes判別，進(jìn)行巖爆傾向性等級評價。但考慮到Bayes判別模型在分類邊界上受先驗概率影響，位于分類邊界附近，易出現(xiàn)分類錯誤，如圖3所示拒絕區(qū)域內(nèi)樣本最大后驗概率差距較小，位于分類邊界附近。筆者提出一種通過閾值改進(jìn)Bayes判別方法，提高Bayes判別的準(zhǔn)確性和不同工況條件下巖爆預(yù)測的泛化能力，克服傳統(tǒng)Bayes決策邊界決策不確定性的缺點。閾值的確定會根據(jù)樣本數(shù)量做出動態(tài)調(diào)整。

圖3 Bayes判別的巖爆預(yù)測部分?jǐn)?shù)據(jù)Fig. 3 Partial prediction data of rockburst based on Bayes discrimination

根據(jù)式( 8 )計算得出1級巖爆、2級巖爆、3級巖爆、4級巖爆的先驗概率分別為0.205，0.114，0.500，0.182，代入式( 12 )計算出kj為4.4。Bayes計算出44組巖爆數(shù)據(jù)，將該數(shù)據(jù)分別代入式( 13 )計算出ik，按照判別規(guī)則式( 14 )篩選出后驗概率位于分類邊界數(shù)據(jù)，見表4。由表4分析發(fā)現(xiàn)，當(dāng)kj＞ki時拒絕判別，造成巖爆分類樣本部分?jǐn)?shù)據(jù)損失，如表4所示數(shù)據(jù)為決策邊界附近舍棄樣本。巖爆數(shù)據(jù)集各類別的樣本數(shù)通常是不均衡的，造成先驗概率差距較大，一部分樣本后驗概率處在分類邊界處。

表4 決策邊界附近的樣本Table 4 Samples near the decision boundary

如圖4所示，巖爆傾向性等級評價模型準(zhǔn)確率為93.18%，考慮到一些樣本數(shù)據(jù)類別后驗概率差距較小，最大后驗概率方法劃分類別存在不足，即使在本次巖爆分級試驗中分類正確，當(dāng)調(diào)整巖爆數(shù)據(jù)集各類別的樣本數(shù)時，出現(xiàn)在分類邊界上樣本易發(fā)生誤判。因此筆者提出閾值改進(jìn)Bayes判別，對于不滿足閾值樣本拒絕分類，閾值改進(jìn)Bayes判別預(yù)測結(jié)果如圖5所示，表明巖爆傾向性等級評價模型預(yù)測結(jié)果與實際相符。極大提高了巖爆傾向性分級的準(zhǔn)確性。

圖4 Bayes巖爆預(yù)測等級Fig. 4 Bayes rockburst prediction grade

圖5 改進(jìn)Bayes巖爆預(yù)測等級Fig. 5 Improved Bayes rockburst prediction grade

2.5 巖爆預(yù)測模型敏感性分析

為了驗證待測樣本所屬類別分布是否均衡，避免存在先驗概率差距較大，從而造成巖爆預(yù)測模型準(zhǔn)確率下降；同時避免樣本和樣本數(shù)量選擇主觀性的影響，本文分別根據(jù)表1隨機(jī)選取樣本，所屬類別最少樣本數(shù)量與所屬類別最多樣本數(shù)量之比按照1∶1，1∶2，1∶3，1∶4，1∶5和1∶6的比例代入巖爆預(yù)測模型和閾值改進(jìn)巖爆預(yù)測模型。

圖6為樣本數(shù)量對預(yù)測模型的影響，由圖6可知，隨著樣本比例增大，越來越多的樣本后驗概率處于分類邊界附近，巖爆預(yù)測模型準(zhǔn)確率隨之下降，閾值改進(jìn)后巖爆預(yù)測模型準(zhǔn)確率明顯提高。樣本比較小時閾值改進(jìn)模型改善效果不明顯。

圖6 樣本數(shù)量對預(yù)測模型影響Fig. 6 Influence of sample size on prediction model

3 工程實例分析

為了驗證主成分分析和改進(jìn)Bayes判別模型的準(zhǔn)確率和泛化能力，本文選取國內(nèi)隧道、國內(nèi)礦山和國內(nèi)外礦山隧道3個領(lǐng)域的巖爆實例進(jìn)行驗證。

3.1 大相領(lǐng)隧道

選取北京至昆明高速公路中的大相領(lǐng)隧道進(jìn)行研究。該隧道地處四川盆地，地形陡峭、埋深大、地應(yīng)力強(qiáng)，具有典型巖爆現(xiàn)象，選取YK55，YK56，YK61等位置的22組巖爆數(shù)據(jù)進(jìn)行分析[25]。大相領(lǐng)隧道巖爆原始數(shù)據(jù)見表5。將表5數(shù)據(jù)代入式( 8 )，計算每個工程1，2，3，4等級的先驗概率，將計算結(jié)果代入式( 12 )，計算kj為6.5。

表5 大相領(lǐng)隧道巖爆原始數(shù)據(jù)Table 5 Original data of rockburst in Daxiangling tunnel

圖7為大相領(lǐng)隧道Bayes判別巖爆數(shù)據(jù)。

圖7 大相領(lǐng)隧道改進(jìn)Bayes判別巖爆數(shù)據(jù)Fig. 7 Improved Bayes discrimination of rockburst data in Daxiangling tunnel

根據(jù)圖7將1，2，3，4等級的后驗概率分別代入式( 13 )計算出每組數(shù)據(jù)ik，得出ki均大于kj，接受判別，巖爆數(shù)據(jù)的樣本均不在分類邊界附近，故認(rèn)為經(jīng)主成分分析法和Bayes判別的巖爆等級判定是可信的；由圖8可知，大相領(lǐng)隧道預(yù)測結(jié)果與實際結(jié)果相符。主成分分析和改進(jìn)Bayes判別模型應(yīng)用于隧道中巖爆災(zāi)害評價準(zhǔn)確率較高。

圖8 大相領(lǐng)隧道改進(jìn)Bayes判別巖爆預(yù)測等級Fig. 8 Improved Bayes discrimination of rockburst prediction grade in Daxiangling tunnel

3.2 馬路坪礦

馬路坪礦位于洋水背斜東翼北段，開采深度為600～750 m，隨著開采深度和采區(qū)地應(yīng)力增加、圍巖表現(xiàn)出明顯脆性，部分巷道出現(xiàn)巖爆現(xiàn)象[26]。

表6為馬路坪礦巖爆原始數(shù)據(jù)。將表6數(shù)據(jù)代入式( 8 )，計算出每個工程1，2，3，4等級的先驗概率，將計算結(jié)果代入式( 12 )，計算kj為2.5。圖9為馬路坪礦改進(jìn)Bayes判別巖爆數(shù)據(jù)。

表6 馬路坪礦巖爆原始數(shù)據(jù)Table 6 Original data of rockburst in Maluping Mine

圖9 馬路坪礦改進(jìn)Bayes判別巖爆數(shù)據(jù)Fig. 9 Improved Bayes discrimination of distinguish rockburst data in Maluping Mine

根據(jù)圖9將1，2，3，4等級的后驗概率分別代入式( 13 )計算出每個工程數(shù)據(jù)ik，得出ik均大于kj，接受判別。認(rèn)為所有巖爆數(shù)據(jù)樣本在主成分分析法和Bayes判別模型分類中均有效。由圖10可知，馬路坪礦巖爆預(yù)測情況與實際相符，說明主成分分析和改進(jìn)Bayes判別應(yīng)用于國內(nèi)礦山巖爆災(zāi)害評價依然有效。

圖10 馬路坪礦改進(jìn)Bayes判別巖爆預(yù)測等級Fig. 10 Improved Bayes discrimination of rockburst prediction grade in Maluping Mine

3.3 其他國內(nèi)外礦山隧道

筆者選取其他國內(nèi)外10個典型的巖爆工程，包括南非金礦、美國CAD-B礦山、前蘇聯(lián)X礦山和瑞士布魯格水電站地下硐室等，進(jìn)行巖爆災(zāi)害分析[27]，巖爆原始數(shù)據(jù)見表7。將表7數(shù)據(jù)代入式( 8 )計算每個工程1，2，3，4等級的先驗概率，將計算結(jié)果代入式( 12 )，計算kj為1.5。圖11為國內(nèi)外典型工程改進(jìn)Bayes判別巖爆數(shù)據(jù)，根據(jù)圖11將1，2，3，4等級的后驗概率分別代入式( 13 )計算出每組數(shù)據(jù)ik。發(fā)現(xiàn)其中4k，7k，k10不大于kj，所以拒絕判別，剔出4k，7k，k10三個樣本，將剩余樣本再次進(jìn)行Bayes判別運算，判別結(jié)果如圖12和13所示。

圖11 國內(nèi)外礦山隧道Bayes判別巖爆數(shù)據(jù)Fig. 11 Bayes discriminant rockburst data of mine tunnels at home and abroad

圖12 國內(nèi)外礦山隧道巖爆預(yù)測等級Fig. 12 Prediction grade of rockburst in mine tunnels at home and abroad

表7 國內(nèi)外礦山隧道巖爆原始數(shù)據(jù)Table 7 Original data of rockburst at home and abroad

Bayes判別出現(xiàn)2個誤判樣本( 均為拒絕判別樣本 )，閾值修正后預(yù)測結(jié)果與實際巖爆等級相符，大大提高巖爆預(yù)測模型的準(zhǔn)確率。

圖13 國內(nèi)外礦山隧道閾值改進(jìn)Bayes判別巖爆數(shù)據(jù)Fig. 13 Threshold improved Bayes discrimination of rockburst data of mine tunnels at home and abroad

4 結(jié) 論

( 1 ) 提出一種閾值改進(jìn)Bayes判別的巖爆綜合預(yù)測模型，模型克服了傳統(tǒng)的判別方法受原始數(shù)據(jù)代表性影響的缺點，通過閾值修正傳統(tǒng)模型決策邊界分類不確定性的缺點，后驗概率小于閾值拒絕判別，提高了巖爆分級的預(yù)測精度。該模型方法簡單、準(zhǔn)確率高，為巖爆分級預(yù)測提供借鑒。

( 2 ) 在不同的地質(zhì)條件下影響巖爆的因素也存在差異，本文中6個影響巖爆的指標(biāo)中存在參數(shù)信息重疊的問題，故采用主成分分析法用3個新指標(biāo)代替原來指標(biāo)，且盡可能保留原始變量的信息，為判別分析消除影響，提高Bayes判別的正判率。

( 3 ) 采取綜合評判模型對44個樣本進(jìn)行分析，驗證主成分分析和閾值改進(jìn)Bayes判別方法的可行性，并將該方法應(yīng)用到大相領(lǐng)隧道、馬路坪礦和其他國內(nèi)外相關(guān)工程，預(yù)測結(jié)果與礦山實際情況相符，但也存在部分樣本丟失的缺點。丟失樣本問題，有待于更深入探索。