宋英華，龐昭勝，李墨瀟，江 晨，齊 石

(1.武漢理工大學(xué) 中國應(yīng)急管理研究中心，湖北 武漢 430070；2.武漢理工大學(xué) 安全科學(xué)與應(yīng)急管理學(xué)院，湖北 武漢 430070)

0 引言

巖爆是堅硬巖石開采和土木建筑中最常見的由連續(xù)巖體應(yīng)力過大引起的破壞之一，其發(fā)生總是伴隨著巖塊開裂、剝落、拋擲、巨響等現(xiàn)象[1]。巖爆具有突發(fā)性、破壞性強(qiáng)、難控制等特點，極易對現(xiàn)場施工人員及工程設(shè)備造成損傷[2]。隨著地下巖體開采深度和挖掘需求的增加，巖爆事故的發(fā)生頻率與影響程度也越來越嚴(yán)重，巖爆已成為愈發(fā)被重視的工程問題。

由于巖爆機(jī)理的復(fù)雜性，且各指標(biāo)數(shù)據(jù)測取具有較大的不穩(wěn)定性，巖爆分類預(yù)測在世界范圍內(nèi)都一直是迫切需要解決的難題，巖爆烈度等級評價研究主要有3大類[3]：1)第1類是基于巖爆機(jī)理判據(jù)，如應(yīng)力強(qiáng)度Russense判據(jù)、Barton判據(jù)，能量理論中能量比指標(biāo)判據(jù)等直接進(jìn)行評價；2)第2類則需要對巖爆現(xiàn)場進(jìn)行實測，諸如微震法、聲發(fā)射法等；3)第3類是基于巖爆影響因素的綜合預(yù)測方法，是目前巖爆預(yù)測研究的熱點。第3類方法又有2種不同的判別方式，一是基于巖爆工程實例數(shù)據(jù)，如采用XGBoost[4]、神經(jīng)網(wǎng)絡(luò)[5]、隨機(jī)森林[6]等機(jī)器學(xué)習(xí)模型進(jìn)行評判；二是基于巖爆指標(biāo)判據(jù)的預(yù)測方法，主要運(yùn)用模糊綜合評判模型[7]，理想點模型[8]，云模型[9]等進(jìn)行預(yù)測。

針對巖爆的隨機(jī)性和模糊性特點，云模型的使用在巖爆預(yù)測中成為當(dāng)下研究熱點，已被證明具有一定的可靠性。Liu等[10]將云模型與粗糙集理論進(jìn)行結(jié)合，根據(jù)沖擊地壓標(biāo)準(zhǔn)生成正態(tài)云模型，運(yùn)用粗糙集理論確定權(quán)重，最后利用最大隸屬度原理確定巖爆等級；周東良等[11]引進(jìn)改進(jìn)的AHP、熵權(quán)法、博弈論和模糊熵理論組合賦權(quán)法與二維云模型綜合評判；劉曉悅等[12]將指標(biāo)的權(quán)重融合到多維云模型中，生成多指標(biāo)的等級綜合云進(jìn)行預(yù)測并驗證了模型的可靠性及實用性。目前，云模型在巖爆預(yù)測方面的研究中，大多學(xué)者都采用正態(tài)云模型的正向云發(fā)生器算法，而正向云發(fā)生器算法的數(shù)字特征一般是經(jīng)驗值，往往具有較強(qiáng)的主觀性。逆向云發(fā)生器算法是基于數(shù)據(jù)生成數(shù)字特征，具有較強(qiáng)的客觀性，在巖爆預(yù)測研究中的使用較少。由于巖爆災(zāi)害等級確定受多指標(biāo)綜合影響，多維云模型也應(yīng)運(yùn)而生。多維云模型的權(quán)重分配研究目前大多是依靠主觀賦權(quán)與客觀賦權(quán)相結(jié)合的方式，但此類方式往往面臨著主觀權(quán)重的精度以及權(quán)重結(jié)合方法的科學(xué)性不足等問題，使問題復(fù)雜化。

本文采用逆向云發(fā)生器MBCT-SR算法[13]，解決確定云模型數(shù)字特征及權(quán)重時，由于主觀干擾導(dǎo)致的預(yù)測結(jié)果與實際情況存在偏差等問題，該算法基于一定數(shù)量的樣本實例，計算客觀云數(shù)字特征，并結(jié)合多維云模型建立動態(tài)適應(yīng)度函數(shù)，通過改進(jìn)的遺傳算法反求最優(yōu)權(quán)重，從而建立完整的巖爆預(yù)測云模型，為巖爆預(yù)測研究提供更貼合實際案例數(shù)據(jù)的評價方法，使評價結(jié)果更加客觀準(zhǔn)確。

1 巖爆預(yù)測多維云模型

1.1 多維正態(tài)云模型

正態(tài)云模型[14](圖1)可實現(xiàn)定性概念與定量表示的相互轉(zhuǎn)換，反映定性概念的不確定性，其性質(zhì)通過期望Ex、熵En以及超熵He3個數(shù)字特征表示。其中期望Ex為云的重心，最能代表定性概念；熵En表示定性概念的離散程度；超熵He則表示熵的離散程度，是熵的不確定性度量。

圖1 云的數(shù)字特征Fig.1 Digital features of cloud model

在一維正態(tài)云模型的基礎(chǔ)上，引入精確數(shù)值表示的n維集合U={X1,X2,X3,…,Xn}，其中U為n維定量論域，C為U上一定性概念。對于U中任意值X=(x1,x2,x3,…,xn)(X∈U)，都存在X對U有一定約束的隨機(jī)數(shù)μ(X(x1,x2,x3,…,xn))，稱為確定度。多維正態(tài)云模型確定度式[15]為式(1)：

(1)

式中：i表示n維數(shù)據(jù)中第i維數(shù)據(jù)(i=1,2,3,…,n)；yi為服從以熵Eni為均值，超熵的平方Hei2為方差的正態(tài)分布隨機(jī)數(shù)；xi滿足以期望Exi為均值；yi為方差的正態(tài)分布規(guī)律。

多維云模型包含2種云發(fā)生器，一種是正向云發(fā)生器(CGn)，另一種為逆向云發(fā)生器(CGn-1)，2種發(fā)生器可進(jìn)行雙向轉(zhuǎn)換，如圖2所示。

圖2 多維云發(fā)生器Fig.2 Multi-dimensional cloud generator

n維正向云發(fā)生器是指將定性概念通過n組數(shù)字特征N(Exn,Enn,Hen)表示，并生成一定數(shù)量的云滴P(X(x1,x2,x3,…,xn))，μj的過程，公式(1)中用j來表示巖爆等級序列(j=1,2,3,4)；逆向云發(fā)生器是指基于一定數(shù)量的云滴，計算云數(shù)字特征的過程，在云滴較少的情況下數(shù)字特征一般為估計值。

1.2 MBCT-SR算法及綜合確定度

通過式(1)可知數(shù)字特征決定確定度的取值及其波動范圍，傳統(tǒng)云模型數(shù)值特征一般采用經(jīng)驗式[16]進(jìn)行計算，即式(2)：

(2)

式中：各指標(biāo)上下限Cmax、Cmin和固定值K均為經(jīng)驗常數(shù)，這使得云模型規(guī)避了巖爆預(yù)測的模糊性特征，導(dǎo)致主觀性較強(qiáng)。云模型的逆向云發(fā)生器算法基于大量云滴數(shù)據(jù)，客觀計算云模型數(shù)字特征，減小經(jīng)驗常數(shù)的主觀性。但傳統(tǒng)的SBCT-1stM逆向算法常伴隨計算結(jié)果漂移及不穩(wěn)定等現(xiàn)象，因此，本文采用多步還原的逆向云變換MBCT-SR算法，提高結(jié)果的準(zhǔn)確度和穩(wěn)定性。具體算法步驟如下：

Step 1：輸入m個樣本點xk，k表示m個樣本點中所取樣本序列(k=1,2,3,…,m)。

Step 2：計算樣本均值，得到期望Ex的估計值，即

(3)

(4)

Step 4：將Step 3中的結(jié)果代入式(5)計算出En2與He2的估計值，即式(5)：

(5)

多維模型的建立需結(jié)合各指標(biāo)權(quán)重，由式(6)計算綜合確定度Ω[17]，其中i為n維各項數(shù)據(jù)序列(i=1,2,3,…,n)，最終達(dá)到預(yù)測的目的。

(6)

2 指標(biāo)權(quán)重反分析

依據(jù)文獻(xiàn)收集[3,6,11,12,16-18]的真實樣本數(shù)據(jù)192組，采用遺傳算法對指標(biāo)權(quán)重進(jìn)行反分析。設(shè)指標(biāo)權(quán)重為未知變量，以最大化滿足樣本數(shù)據(jù)為結(jié)果，通過優(yōu)化算法求全局最優(yōu)解。具體流程見圖3。此方法降低權(quán)重確定主觀性，其分析步驟為：①確定評價指標(biāo)；②收集樣本實例并進(jìn)行數(shù)據(jù)預(yù)處理；③計算數(shù)字特征，構(gòu)建多維云模型；④由式(6)，建立優(yōu)化適應(yīng)度函數(shù)；⑤利用遺傳算法，尋求全局最優(yōu)解；⑥回代，驗證完整多維云模型的有效性及可行性。

圖3 指標(biāo)權(quán)重反分析方法流程Fig.3 Flow chart of index weight back analysis method

遺傳算法只能從適應(yīng)度函數(shù)中獲取信息，故步驟④為反求權(quán)重的核心內(nèi)容，在求取適應(yīng)度函數(shù)時，由式(6)可得綜合確定度中的參數(shù)yi為服從以熵Eni為均值，超熵的平方Hei2為方差的正態(tài)分布隨機(jī)數(shù)(i=1,2,3,…,n)，故yi需由MBCT-SR算法計算獲得其正態(tài)分布均值與方差，再進(jìn)行隨機(jī)取值，最終建立動態(tài)適應(yīng)度函數(shù)fitness[18]如式(7)：

(7)

式中：j為巖爆等級序列；i為評價指標(biāo)序列；k為樣本序號；Ωj為第j個巖爆等級的綜合確定度；n，p，m分別為指標(biāo)總數(shù)、巖爆等級總數(shù)和樣本總數(shù)；qk與Qk分別是第k個樣本中的預(yù)測等級和實際等級。

由式(7)可知，優(yōu)化算法以預(yù)測結(jié)果滿足實例結(jié)果最大化為目標(biāo)，規(guī)避主觀干擾，建立適應(yīng)度函數(shù)，求出最優(yōu)權(quán)重。求得權(quán)重代入式(6)，根據(jù)綜合確定度獲得最終評價結(jié)果。

3 模型實例驗證

3.1 評價指標(biāo)選取及數(shù)據(jù)預(yù)處理

巖爆等級受多因素綜合影響，本文基于真實數(shù)據(jù)分析，規(guī)避主觀因素，對數(shù)據(jù)樣本具有較高的要求，數(shù)據(jù)精度直接影響預(yù)測結(jié)果。結(jié)合以往研究經(jīng)驗，通過文獻(xiàn)收集[3,6,11,12,16-21]的方法，最終確定數(shù)據(jù)樣本較多、影響程度較大的3項指標(biāo)：應(yīng)力比Ts=σθ/σt、巖石脆性指數(shù)B=σc/σt以及彈性應(yīng)變儲能指數(shù)Wet作為巖爆等級預(yù)測評價指標(biāo)，參考王元漢等[22]的相關(guān)研究，巖爆等級依托于各指標(biāo)分為無巖爆(Ⅰ)、輕微巖爆(Ⅱ)、中等巖爆(Ⅲ)、強(qiáng)巖爆(Ⅳ)4個等級評價。

本文研究數(shù)據(jù)樣本較全，不存在數(shù)據(jù)缺失的情況，但數(shù)據(jù)收集時無法避免存在噪點的情況出現(xiàn)。因此采用偏差分析結(jié)合箱線圖進(jìn)行預(yù)處理，采用aσ(a=2或3)規(guī)則[23]進(jìn)行修剪，其中σ表示原始數(shù)據(jù)方差，a表示以a倍方差為臨界值對離散數(shù)據(jù)進(jìn)行截取。最終選取192 組國內(nèi)外實際案例進(jìn)行預(yù)測，其分布規(guī)律如圖4所示。

圖4 原始數(shù)據(jù)樣本箱線Fig.4 Box line diagram of original data sample

從圖4所知，由于各指標(biāo)量綱不同，各指標(biāo)數(shù)值差異較大，不利于進(jìn)行科學(xué)計算及綜合分析，故按式(8)對各指標(biāo)進(jìn)行歸一化處理，即得

(8)

式中：xik為第i個指標(biāo)中第k個原始樣本數(shù)據(jù)；ximin，ximax分別是第i個指標(biāo)中原始樣本的最小值和最大值。

式(8)中原始數(shù)據(jù)歸一化后其分布結(jié)果如圖5所示。由圖5可得，預(yù)測指標(biāo)巖石脆性指數(shù)B=σc/σt依舊存在一定數(shù)量的離散點，數(shù)據(jù)中心偏移較為嚴(yán)重，而應(yīng)力比Ts=σθ/σt以及彈性應(yīng)變儲能指數(shù)Wet分布較為均勻。指標(biāo)上下四分位數(shù)跨度較小，其中Ts下四分位偏低，中位數(shù)靠近0.5，無離散點。彈性應(yīng)變儲能指數(shù)Wet上四分位數(shù)偏高，中位數(shù)稍低，但總體分布比較合理。由歸一化數(shù)據(jù)分析初步猜測，本數(shù)據(jù)樣本指標(biāo)B穩(wěn)定性較差，權(quán)重應(yīng)最小，這也滿足巖爆各指標(biāo)權(quán)重分配的主觀認(rèn)知。

圖5 歸一化后樣本箱線Fig.5 Normalized sample box plot

選取29組國內(nèi)外實例數(shù)據(jù)做測試集，列舉部分如表1所示，以剩余163組數(shù)據(jù)為訓(xùn)練樣本，并與其他云模型預(yù)測方法對比，驗證其準(zhǔn)確性。

表1 國內(nèi)外巖爆實例數(shù)據(jù)Table 1 Data of rock burst cases at home and abroad

3.2 計算數(shù)值特征

將163 組巖爆實例數(shù)據(jù)按級分類，并采用MBCT-SR算法對各級各指標(biāo)巖爆數(shù)據(jù)分別計算，由于算法計算結(jié)果并不唯一，任選其一組將各數(shù)字特征示于表2 。

表2 數(shù)字特征計算結(jié)果Table 2 Digital feature calculation results

根據(jù)表2中各數(shù)字特征繪制多維云模型，如圖6所示。從圖6所知，“·”描繪數(shù)據(jù)代表無巖爆(Ⅰ級)，“+”描繪數(shù)據(jù)點代表弱巖爆(Ⅱ級)，“○”點代表中巖爆(Ⅲ級)，“*”點表示強(qiáng)巖爆(Ⅳ級)。指標(biāo)Ts與實際巖爆順序相同，相對于其他2個指標(biāo)分布較為均勻，各等級界限較為明顯。指標(biāo)Wet的分布順序與巖爆順序也相同，但在無巖爆與弱巖爆中數(shù)據(jù)跨度較大，其分布特征相比于指標(biāo)Ts較不穩(wěn)定。B維度數(shù)據(jù)點分布順序較亂，跨度較大，很明顯相對于其他2指標(biāo)規(guī)律性較差，最不適用于巖爆評價。此次MBCT-SR算法基于大量樣本數(shù)據(jù)計算數(shù)字特征并生成三維云圖的分析中，對第3.1節(jié)中假設(shè)進(jìn)行了初步驗證，并針對圖6 做進(jìn)一步假設(shè)，求得最終權(quán)重順序應(yīng)滿足ωTs>ωWet>ωB。

圖6 巖爆實例數(shù)據(jù)聚類云模型Fig.6 Clustering cloud model of rock burst case data

3.3 求取最優(yōu)權(quán)重

將各數(shù)字特征代入式(7)構(gòu)建適應(yīng)度函數(shù)反求權(quán)重，在遺傳算法中，適應(yīng)度函數(shù)的每次使用將重新計算數(shù)字特征，確保結(jié)果真實性。

常規(guī)遺傳算法迭代曲線保留原始樣本最優(yōu)解，尋求全局最優(yōu)解，迭代曲線呈現(xiàn)單調(diào)上升趨勢。本文由式(7)使用動態(tài)適應(yīng)度函數(shù)，同一最優(yōu)解會出現(xiàn)不同結(jié)果，滿足巖爆評價實際情況。為降低時間成本，減少迭代次數(shù)，最快獲取較為穩(wěn)定全局最優(yōu)解。本文增加種群數(shù)量，設(shè)置初始化種群數(shù)量N=1 000，迭代次數(shù)f=100，采用錦標(biāo)賽選擇和精英保留方法，增加圖3中4個權(quán)重約束條件進(jìn)行約束，運(yùn)行Matlab最終獲得迭代曲線如圖7所示。

圖7 迭代曲線Fig.7 Iteration curve

由圖7得迭代48 次后滿足條件樣本數(shù)量保持穩(wěn)定，此時各指標(biāo)最優(yōu)權(quán)重為ωTs=0.632 1、ωWet=0.257 7和ωB=0.110 2，最優(yōu)權(quán)重與假設(shè)權(quán)重比例相符，滿足樣本實例。

3.4 實例驗證與對比

將所求結(jié)果通過前文表1中未參與反求權(quán)重的29組實例樣本進(jìn)行評價。由于巖爆事故本身具有較強(qiáng)隨機(jī)性，應(yīng)將評價結(jié)果的合理波動考慮進(jìn)來，增大預(yù)測范圍，分析評價結(jié)果更多的可能性。本文對樣本數(shù)據(jù)進(jìn)行客觀評價，采取保守預(yù)測方法，將確定度差值不足0.1的巖爆等級進(jìn)行傾向性預(yù)測，其預(yù)測結(jié)果進(jìn)行跨區(qū)間表示，例如：Ⅱ～Ⅲ級巖爆，此方法滿足實際巖爆預(yù)測要求，使預(yù)測更加嚴(yán)謹(jǐn)，并與其他云模型評價方法進(jìn)行結(jié)果對比，進(jìn)一步驗證本文方法的可靠性與有效性，具體情況見表3。

表3 巖爆傾向性評價結(jié)果Table 3 Evaluation index judgment interval

結(jié)果顯示，預(yù)測結(jié)果與實際基本相符，與其他模型的預(yù)測結(jié)果基本吻合。本文與反賦權(quán)一維云模型權(quán)重確定方法相同。由表3得，反賦權(quán)一維云模型傾向性預(yù)測準(zhǔn)確率為72%。組合賦權(quán)多維云模型預(yù)測結(jié)果準(zhǔn)確率為86%，比一維正向云發(fā)生器的準(zhǔn)確率更高。在逆向云發(fā)生器中，MBCT-SR算法與優(yōu)化算法結(jié)合評價精準(zhǔn)預(yù)測準(zhǔn)確率可達(dá)89%，傾向性預(yù)測準(zhǔn)確率可達(dá)100%。預(yù)測結(jié)果對比表明，本文評價模型更加合理且有效，算法充分考慮實際巖爆的隨機(jī)性與模糊性，更貼合實際情況。

4 結(jié)論

1)結(jié)合當(dāng)前機(jī)器學(xué)習(xí)的熱門方向，依據(jù)國內(nèi)外192 組巖爆實例數(shù)據(jù)，選取應(yīng)力比Ts=σθ/σt、巖石脆性指數(shù)B=σc/σt以及彈性應(yīng)變儲能指數(shù)Wet作為評價指標(biāo)，將逆向云發(fā)生器與優(yōu)化算法相結(jié)合建立綜合評價模型。本模型引進(jìn)MBCT-SR算法，規(guī)避常規(guī)正向云發(fā)生器主觀性過強(qiáng)、絕對性較大的缺點，建立更加客觀的多維巖爆預(yù)測云模型。

2)指標(biāo)權(quán)重優(yōu)化結(jié)果來自于真實數(shù)據(jù)與多維云模型。本文依據(jù)圖5、6中的數(shù)據(jù)可視化分析，對權(quán)重分布規(guī)律進(jìn)行了滿足ωTs>ωWet>ωB的初步假設(shè)，并由最終計算結(jié)果求得最優(yōu)權(quán)重為ωTs=0.632 1、ωWet=0.257 7和ωB=0.110 2，對假設(shè)進(jìn)行了科學(xué)驗證。使得本指標(biāo)權(quán)重確定方法更加準(zhǔn)確嚴(yán)謹(jǐn)，最終結(jié)果與主觀認(rèn)知相符。

3)相比其他預(yù)測模型可得，本模型最大化降低主觀因素干擾，在選取應(yīng)力比Ts=σθ/σt、巖石脆性指數(shù)B=σc/σt以及彈性應(yīng)變儲能指數(shù)Wet作為評價指標(biāo)進(jìn)行預(yù)測時，精準(zhǔn)預(yù)測率達(dá)89%，傾向性預(yù)測結(jié)果準(zhǔn)確率可達(dá)到100%，其評價結(jié)果更加準(zhǔn)確、真實。