宋英華,龐昭勝,李墨瀟,江 晨,齊 石
(1.武漢理工大學(xué) 中國(guó)應(yīng)急管理研究中心,湖北 武漢 430070;2.武漢理工大學(xué) 安全科學(xué)與應(yīng)急管理學(xué)院,湖北 武漢 430070)
巖爆是堅(jiān)硬巖石開采和土木建筑中最常見的由連續(xù)巖體應(yīng)力過大引起的破壞之一,其發(fā)生總是伴隨著巖塊開裂、剝落、拋擲、巨響等現(xiàn)象[1]。巖爆具有突發(fā)性、破壞性強(qiáng)、難控制等特點(diǎn),極易對(duì)現(xiàn)場(chǎng)施工人員及工程設(shè)備造成損傷[2]。隨著地下巖體開采深度和挖掘需求的增加,巖爆事故的發(fā)生頻率與影響程度也越來越嚴(yán)重,巖爆已成為愈發(fā)被重視的工程問題。
由于巖爆機(jī)理的復(fù)雜性,且各指標(biāo)數(shù)據(jù)測(cè)取具有較大的不穩(wěn)定性,巖爆分類預(yù)測(cè)在世界范圍內(nèi)都一直是迫切需要解決的難題,巖爆烈度等級(jí)評(píng)價(jià)研究主要有3大類[3]:1)第1類是基于巖爆機(jī)理判據(jù),如應(yīng)力強(qiáng)度Russense判據(jù)、Barton判據(jù),能量理論中能量比指標(biāo)判據(jù)等直接進(jìn)行評(píng)價(jià);2)第2類則需要對(duì)巖爆現(xiàn)場(chǎng)進(jìn)行實(shí)測(cè),諸如微震法、聲發(fā)射法等;3)第3類是基于巖爆影響因素的綜合預(yù)測(cè)方法,是目前巖爆預(yù)測(cè)研究的熱點(diǎn)。第3類方法又有2種不同的判別方式,一是基于巖爆工程實(shí)例數(shù)據(jù),如采用XGBoost[4]、神經(jīng)網(wǎng)絡(luò)[5]、隨機(jī)森林[6]等機(jī)器學(xué)習(xí)模型進(jìn)行評(píng)判;二是基于巖爆指標(biāo)判據(jù)的預(yù)測(cè)方法,主要運(yùn)用模糊綜合評(píng)判模型[7],理想點(diǎn)模型[8],云模型[9]等進(jìn)行預(yù)測(cè)。
針對(duì)巖爆的隨機(jī)性和模糊性特點(diǎn),云模型的使用在巖爆預(yù)測(cè)中成為當(dāng)下研究熱點(diǎn),已被證明具有一定的可靠性。Liu等[10]將云模型與粗糙集理論進(jìn)行結(jié)合,根據(jù)沖擊地壓標(biāo)準(zhǔn)生成正態(tài)云模型,運(yùn)用粗糙集理論確定權(quán)重,最后利用最大隸屬度原理確定巖爆等級(jí);周東良等[11]引進(jìn)改進(jìn)的AHP、熵權(quán)法、博弈論和模糊熵理論組合賦權(quán)法與二維云模型綜合評(píng)判;劉曉悅等[12]將指標(biāo)的權(quán)重融合到多維云模型中,生成多指標(biāo)的等級(jí)綜合云進(jìn)行預(yù)測(cè)并驗(yàn)證了模型的可靠性及實(shí)用性。目前,云模型在巖爆預(yù)測(cè)方面的研究中,大多學(xué)者都采用正態(tài)云模型的正向云發(fā)生器算法,而正向云發(fā)生器算法的數(shù)字特征一般是經(jīng)驗(yàn)值,往往具有較強(qiáng)的主觀性。逆向云發(fā)生器算法是基于數(shù)據(jù)生成數(shù)字特征,具有較強(qiáng)的客觀性,在巖爆預(yù)測(cè)研究中的使用較少。由于巖爆災(zāi)害等級(jí)確定受多指標(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)重時(shí),由于主觀干擾導(dǎo)致的預(yù)測(cè)結(jié)果與實(shí)際情況存在偏差等問題,該算法基于一定數(shù)量的樣本實(shí)例,計(jì)算客觀云數(shù)字特征,并結(jié)合多維云模型建立動(dòng)態(tài)適應(yīng)度函數(shù),通過改進(jìn)的遺傳算法反求最優(yōu)權(quán)重,從而建立完整的巖爆預(yù)測(cè)云模型,為巖爆預(yù)測(cè)研究提供更貼合實(shí)際案例數(shù)據(jù)的評(píng)價(jià)方法,使評(píng)價(jià)結(jié)果更加客觀準(zhǔn)確。
正態(tài)云模型[14](圖1)可實(shí)現(xiàn)定性概念與定量表示的相互轉(zhuǎn)換,反映定性概念的不確定性,其性質(zhì)通過期望Ex、熵En以及超熵He3個(gè)數(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上一定性概念。對(duì)于U中任意值X=(x1,x2,x3,…,xn)(X∈U),都存在X對(duì)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í)序列(j=1,2,3,4);逆向云發(fā)生器是指基于一定數(shù)量的云滴,計(jì)算云數(shù)字特征的過程,在云滴較少的情況下數(shù)字特征一般為估計(jì)值。
通過式(1)可知數(shù)字特征決定確定度的取值及其波動(dòng)范圍,傳統(tǒng)云模型數(shù)值特征一般采用經(jīng)驗(yàn)式[16]進(jìn)行計(jì)算,即式(2):
(2)
式中:各指標(biāo)上下限Cmax、Cmin和固定值K均為經(jīng)驗(yàn)常數(shù),這使得云模型規(guī)避了巖爆預(yù)測(cè)的模糊性特征,導(dǎo)致主觀性較強(qiáng)。云模型的逆向云發(fā)生器算法基于大量云滴數(shù)據(jù),客觀計(jì)算云模型數(shù)字特征,減小經(jīng)驗(yàn)常數(shù)的主觀性。但傳統(tǒng)的SBCT-1stM逆向算法常伴隨計(jì)算結(jié)果漂移及不穩(wěn)定等現(xiàn)象,因此,本文采用多步還原的逆向云變換MBCT-SR算法,提高結(jié)果的準(zhǔn)確度和穩(wěn)定性。具體算法步驟如下:
Step 1:輸入m個(gè)樣本點(diǎn)xk,k表示m個(gè)樣本點(diǎn)中所取樣本序列(k=1,2,3,…,m)。
Step 2:計(jì)算樣本均值,得到期望Ex的估計(jì)值,即
(3)
(4)
Step 4:將Step 3中的結(jié)果代入式(5)計(jì)算出En2與He2的估計(jì)值,即式(5):
(5)
多維模型的建立需結(jié)合各指標(biāo)權(quán)重,由式(6)計(jì)算綜合確定度Ω[17],其中i為n維各項(xiàng)數(shù)據(jù)序列(i=1,2,3,…,n),最終達(dá)到預(yù)測(cè)的目的。
(6)
依據(jù)文獻(xiàn)收集[3,6,11,12,16-18]的真實(shí)樣本數(shù)據(jù)192組,采用遺傳算法對(duì)指標(biāo)權(quán)重進(jìn)行反分析。設(shè)指標(biāo)權(quán)重為未知變量,以最大化滿足樣本數(shù)據(jù)為結(jié)果,通過優(yōu)化算法求全局最優(yōu)解。具體流程見圖3。此方法降低權(quán)重確定主觀性,其分析步驟為:①確定評(píng)價(jià)指標(biāo);②收集樣本實(shí)例并進(jìn)行數(shù)據(jù)預(yù)處理;③計(jì)算數(shù)字特征,構(gòu)建多維云模型;④由式(6),建立優(yōu)化適應(yīng)度函數(shù);⑤利用遺傳算法,尋求全局最優(yōu)解;⑥回代,驗(yàn)證完整多維云模型的有效性及可行性。
圖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ù)時(shí),由式(6)可得綜合確定度中的參數(shù)yi為服從以熵Eni為均值,超熵的平方Hei2為方差的正態(tài)分布隨機(jī)數(shù)(i=1,2,3,…,n),故yi需由MBCT-SR算法計(jì)算獲得其正態(tài)分布均值與方差,再進(jìn)行隨機(jī)取值,最終建立動(dòng)態(tài)適應(yīng)度函數(shù)fitness[18]如式(7):
(7)
式中:j為巖爆等級(jí)序列;i為評(píng)價(jià)指標(biāo)序列;k為樣本序號(hào);Ωj為第j個(gè)巖爆等級(jí)的綜合確定度;n,p,m分別為指標(biāo)總數(shù)、巖爆等級(jí)總數(shù)和樣本總數(shù);qk與Qk分別是第k個(gè)樣本中的預(yù)測(cè)等級(jí)和實(shí)際等級(jí)。
由式(7)可知,優(yōu)化算法以預(yù)測(cè)結(jié)果滿足實(shí)例結(jié)果最大化為目標(biāo),規(guī)避主觀干擾,建立適應(yīng)度函數(shù),求出最優(yōu)權(quán)重。求得權(quán)重代入式(6),根據(jù)綜合確定度獲得最終評(píng)價(jià)結(jié)果。
巖爆等級(jí)受多因素綜合影響,本文基于真實(shí)數(shù)據(jù)分析,規(guī)避主觀因素,對(duì)數(shù)據(jù)樣本具有較高的要求,數(shù)據(jù)精度直接影響預(yù)測(cè)結(jié)果。結(jié)合以往研究經(jīng)驗(yàn),通過文獻(xiàn)收集[3,6,11,12,16-21]的方法,最終確定數(shù)據(jù)樣本較多、影響程度較大的3項(xiàng)指標(biāo):應(yīng)力比Ts=σθ/σt、巖石脆性指數(shù)B=σc/σt以及彈性應(yīng)變儲(chǔ)能指數(shù)Wet作為巖爆等級(jí)預(yù)測(cè)評(píng)價(jià)指標(biāo),參考王元漢等[22]的相關(guān)研究,巖爆等級(jí)依托于各指標(biāo)分為無巖爆(Ⅰ)、輕微巖爆(Ⅱ)、中等巖爆(Ⅲ)、強(qiáng)巖爆(Ⅳ)4個(gè)等級(jí)評(píng)價(jià)。
本文研究數(shù)據(jù)樣本較全,不存在數(shù)據(jù)缺失的情況,但數(shù)據(jù)收集時(shí)無法避免存在噪點(diǎn)的情況出現(xiàn)。因此采用偏差分析結(jié)合箱線圖進(jìn)行預(yù)處理,采用aσ(a=2或3)規(guī)則[23]進(jìn)行修剪,其中σ表示原始數(shù)據(jù)方差,a表示以a倍方差為臨界值對(duì)離散數(shù)據(jù)進(jìn)行截取。最終選取192 組國(guó)內(nèi)外實(shí)際案例進(jìn)行預(yù)測(cè),其分布規(guī)律如圖4所示。
圖4 原始數(shù)據(jù)樣本箱線Fig.4 Box line diagram of original data sample
從圖4所知,由于各指標(biāo)量綱不同,各指標(biāo)數(shù)值差異較大,不利于進(jìn)行科學(xué)計(jì)算及綜合分析,故按式(8)對(duì)各指標(biāo)進(jìn)行歸一化處理,即得
(8)
式中:xik為第i個(gè)指標(biāo)中第k個(gè)原始樣本數(shù)據(jù);ximin,ximax分別是第i個(gè)指標(biāo)中原始樣本的最小值和最大值。
式(8)中原始數(shù)據(jù)歸一化后其分布結(jié)果如圖5所示。由圖5可得,預(yù)測(cè)指標(biāo)巖石脆性指數(shù)B=σc/σt依舊存在一定數(shù)量的離散點(diǎn),數(shù)據(jù)中心偏移較為嚴(yán)重,而應(yīng)力比Ts=σθ/σt以及彈性應(yīng)變儲(chǔ)能指數(shù)Wet分布較為均勻。指標(biāo)上下四分位數(shù)跨度較小,其中Ts下四分位偏低,中位數(shù)靠近0.5,無離散點(diǎn)。彈性應(yīng)變儲(chǔ)能指數(shù)Wet上四分位數(shù)偏高,中位數(shù)稍低,但總體分布比較合理。由歸一化數(shù)據(jù)分析初步猜測(cè),本數(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組國(guó)內(nèi)外實(shí)例數(shù)據(jù)做測(cè)試集,列舉部分如表1所示,以剩余163組數(shù)據(jù)為訓(xùn)練樣本,并與其他云模型預(yù)測(cè)方法對(duì)比,驗(yàn)證其準(zhǔn)確性。
表1 國(guó)內(nèi)外巖爆實(shí)例數(shù)據(jù)Table 1 Data of rock burst cases at home and abroad
將163 組巖爆實(shí)例數(shù)據(jù)按級(jí)分類,并采用MBCT-SR算法對(duì)各級(jí)各指標(biāo)巖爆數(shù)據(jù)分別計(jì)算,由于算法計(jì)算結(jié)果并不唯一,任選其一組將各數(shù)字特征示于表2 。
表2 數(shù)字特征計(jì)算結(jié)果Table 2 Digital feature calculation results
根據(jù)表2中各數(shù)字特征繪制多維云模型,如圖6所示。從圖6所知,“·”描繪數(shù)據(jù)代表無巖爆(Ⅰ級(jí)),“+”描繪數(shù)據(jù)點(diǎn)代表弱巖爆(Ⅱ級(jí)),“○”點(diǎn)代表中巖爆(Ⅲ級(jí)),“*”點(diǎn)表示強(qiáng)巖爆(Ⅳ級(jí))。指標(biāo)Ts與實(shí)際巖爆順序相同,相對(duì)于其他2個(gè)指標(biāo)分布較為均勻,各等級(jí)界限較為明顯。指標(biāo)Wet的分布順序與巖爆順序也相同,但在無巖爆與弱巖爆中數(shù)據(jù)跨度較大,其分布特征相比于指標(biāo)Ts較不穩(wěn)定。B維度數(shù)據(jù)點(diǎn)分布順序較亂,跨度較大,很明顯相對(duì)于其他2指標(biāo)規(guī)律性較差,最不適用于巖爆評(píng)價(jià)。此次MBCT-SR算法基于大量樣本數(shù)據(jù)計(jì)算數(shù)字特征并生成三維云圖的分析中,對(duì)第3.1節(jié)中假設(shè)進(jìn)行了初步驗(yàn)證,并針對(duì)圖6 做進(jìn)一步假設(shè),求得最終權(quán)重順序應(yīng)滿足ωTs>ωWet>ωB。
圖6 巖爆實(shí)例數(shù)據(jù)聚類云模型Fig.6 Clustering cloud model of rock burst case data
將各數(shù)字特征代入式(7)構(gòu)建適應(yīng)度函數(shù)反求權(quán)重,在遺傳算法中,適應(yīng)度函數(shù)的每次使用將重新計(jì)算數(shù)字特征,確保結(jié)果真實(shí)性。
常規(guī)遺傳算法迭代曲線保留原始樣本最優(yōu)解,尋求全局最優(yōu)解,迭代曲線呈現(xiàn)單調(diào)上升趨勢(shì)。本文由式(7)使用動(dòng)態(tài)適應(yīng)度函數(shù),同一最優(yōu)解會(huì)出現(xiàn)不同結(jié)果,滿足巖爆評(píng)價(jià)實(shí)際情況。為降低時(shí)間成本,減少迭代次數(shù),最快獲取較為穩(wěn)定全局最優(yōu)解。本文增加種群數(shù)量,設(shè)置初始化種群數(shù)量N=1 000,迭代次數(shù)f=100,采用錦標(biāo)賽選擇和精英保留方法,增加圖3中4個(gè)權(quán)重約束條件進(jìn)行約束,運(yùn)行Matlab最終獲得迭代曲線如圖7所示。
圖7 迭代曲線Fig.7 Iteration curve
由圖7得迭代48 次后滿足條件樣本數(shù)量保持穩(wěn)定,此時(shí)各指標(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)重比例相符,滿足樣本實(shí)例。
將所求結(jié)果通過前文表1中未參與反求權(quán)重的29組實(shí)例樣本進(jìn)行評(píng)價(jià)。由于巖爆事故本身具有較強(qiáng)隨機(jī)性,應(yīng)將評(píng)價(jià)結(jié)果的合理波動(dòng)考慮進(jìn)來,增大預(yù)測(cè)范圍,分析評(píng)價(jià)結(jié)果更多的可能性。本文對(duì)樣本數(shù)據(jù)進(jìn)行客觀評(píng)價(jià),采取保守預(yù)測(cè)方法,將確定度差值不足0.1的巖爆等級(jí)進(jìn)行傾向性預(yù)測(cè),其預(yù)測(cè)結(jié)果進(jìn)行跨區(qū)間表示,例如:Ⅱ~Ⅲ級(jí)巖爆,此方法滿足實(shí)際巖爆預(yù)測(cè)要求,使預(yù)測(cè)更加嚴(yán)謹(jǐn),并與其他云模型評(píng)價(jià)方法進(jìn)行結(jié)果對(duì)比,進(jìn)一步驗(yàn)證本文方法的可靠性與有效性,具體情況見表3。
表3 巖爆傾向性評(píng)價(jià)結(jié)果Table 3 Evaluation index judgment interval
結(jié)果顯示,預(yù)測(cè)結(jié)果與實(shí)際基本相符,與其他模型的預(yù)測(cè)結(jié)果基本吻合。本文與反賦權(quán)一維云模型權(quán)重確定方法相同。由表3得,反賦權(quán)一維云模型傾向性預(yù)測(cè)準(zhǔn)確率為72%。組合賦權(quán)多維云模型預(yù)測(cè)結(jié)果準(zhǔn)確率為86%,比一維正向云發(fā)生器的準(zhǔn)確率更高。在逆向云發(fā)生器中,MBCT-SR算法與優(yōu)化算法結(jié)合評(píng)價(jià)精準(zhǔn)預(yù)測(cè)準(zhǔn)確率可達(dá)89%,傾向性預(yù)測(cè)準(zhǔn)確率可達(dá)100%。預(yù)測(cè)結(jié)果對(duì)比表明,本文評(píng)價(jià)模型更加合理且有效,算法充分考慮實(shí)際巖爆的隨機(jī)性與模糊性,更貼合實(shí)際情況。
1)結(jié)合當(dāng)前機(jī)器學(xué)習(xí)的熱門方向,依據(jù)國(guó)內(nèi)外192 組巖爆實(shí)例數(shù)據(jù),選取應(yīng)力比Ts=σθ/σt、巖石脆性指數(shù)B=σc/σt以及彈性應(yīng)變儲(chǔ)能指數(shù)Wet作為評(píng)價(jià)指標(biāo),將逆向云發(fā)生器與優(yōu)化算法相結(jié)合建立綜合評(píng)價(jià)模型。本模型引進(jìn)MBCT-SR算法,規(guī)避常規(guī)正向云發(fā)生器主觀性過強(qiáng)、絕對(duì)性較大的缺點(diǎn),建立更加客觀的多維巖爆預(yù)測(cè)云模型。
2)指標(biāo)權(quán)重優(yōu)化結(jié)果來自于真實(shí)數(shù)據(jù)與多維云模型。本文依據(jù)圖5、6中的數(shù)據(jù)可視化分析,對(duì)權(quán)重分布規(guī)律進(jìn)行了滿足ωTs>ωWet>ωB的初步假設(shè),并由最終計(jì)算結(jié)果求得最優(yōu)權(quán)重為ωTs=0.632 1、ωWet=0.257 7和ωB=0.110 2,對(duì)假設(shè)進(jìn)行了科學(xué)驗(yàn)證。使得本指標(biāo)權(quán)重確定方法更加準(zhǔn)確嚴(yán)謹(jǐn),最終結(jié)果與主觀認(rèn)知相符。
3)相比其他預(yù)測(cè)模型可得,本模型最大化降低主觀因素干擾,在選取應(yīng)力比Ts=σθ/σt、巖石脆性指數(shù)B=σc/σt以及彈性應(yīng)變儲(chǔ)能指數(shù)Wet作為評(píng)價(jià)指標(biāo)進(jìn)行預(yù)測(cè)時(shí),精準(zhǔn)預(yù)測(cè)率達(dá)89%,傾向性預(yù)測(cè)結(jié)果準(zhǔn)確率可達(dá)到100%,其評(píng)價(jià)結(jié)果更加準(zhǔn)確、真實(shí)。
中國(guó)安全生產(chǎn)科學(xué)技術(shù)2022年3期