高煜，胡賓鑫，朱峰，張華，宋廣東，高國(guó)防，龐江波，鐘國(guó)棟，全霓

(1.齊魯工業(yè)大學(xué)(山東省科學(xué)院) 激光研究所， 山東 濟(jì)南 250300；2.陜西高速集團(tuán)白泉分公司， 陜西 西安 710065)

0 引言

隨著礦山開采深度和強(qiáng)度的增加，巖體失穩(wěn)引發(fā)的沖擊地壓、頂板冒落、突水、瓦斯突出等動(dòng)力災(zāi)害頻發(fā)，嚴(yán)重影響礦山安全生產(chǎn)。微震監(jiān)測(cè)技術(shù)是一種有效的巖體穩(wěn)定性監(jiān)測(cè)預(yù)警手段，廣泛應(yīng)用于煤礦、非煤礦山、隧道、硐室、大壩、邊坡、高鐵路基等工程穩(wěn)定性監(jiān)測(cè)預(yù)警領(lǐng)域。該技術(shù)利用微震監(jiān)測(cè)系統(tǒng)進(jìn)行現(xiàn)場(chǎng)實(shí)時(shí)監(jiān)測(cè)，結(jié)合震源定位技術(shù)確定微震事件的時(shí)空信息及能量，從而對(duì)巖體變形導(dǎo)致的破壞活動(dòng)范圍、穩(wěn)定性及發(fā)展趨勢(shì)等作出科學(xué)評(píng)價(jià)。微震初至波(一般為P波)是檢波器最先接收到的地震波，位于環(huán)境噪聲和有效信號(hào)的分界處，此處伴隨能量急劇變化。準(zhǔn)確地拾取微震初至波到時(shí)是實(shí)現(xiàn)震源定位的前提，是微震監(jiān)測(cè)關(guān)鍵技術(shù)之一[1]。

傳統(tǒng)的微震初至波到時(shí)拾取方法主要有人工拾取和自動(dòng)拾取。人工拾取方法主要依靠肉眼判別微震初至波到時(shí)，效率低，耗時(shí)長(zhǎng)。自動(dòng)拾取方法根據(jù)微震波的振幅、頻率等特征構(gòu)造特征函數(shù)，利用數(shù)學(xué)算法檢測(cè)微震初至波到時(shí)，主要有長(zhǎng)短時(shí)窗能量比值(Short Term Average/Long Term Average，STA/LTA)法[2-4]、赤池信息準(zhǔn)則[5]、高階統(tǒng)計(jì)量方法[6]、小波變換[7-8]、神經(jīng)網(wǎng)絡(luò)法、分形分維法等[9-11]。其中STA/LTA法應(yīng)用較普遍。該方法中STA為短時(shí)窗平均值，代表微震初至波到來時(shí)突變的能量；LTA為長(zhǎng)時(shí)窗平均值，表示噪聲能量。當(dāng)檢波器檢測(cè)到微震初至波時(shí)，相應(yīng)的STA/LTA明顯增大，達(dá)到觸發(fā)閾值時(shí)對(duì)應(yīng)的時(shí)刻認(rèn)為是初至波到時(shí)[12]。STA/LTA法計(jì)算簡(jiǎn)單，適合實(shí)時(shí)處理且穩(wěn)定可靠，但準(zhǔn)確性嚴(yán)重受制于觸發(fā)閾值的設(shè)定、長(zhǎng)短時(shí)窗大小及特征函數(shù)的選取，尤其是對(duì)低信噪比信號(hào)的拾取準(zhǔn)確率較低。

針對(duì)上述問題，以實(shí)際煤礦巖石破裂微震監(jiān)測(cè)為例，提出了一種基于隨機(jī)森林的微震初至波到時(shí)自動(dòng)拾取方法：首先，提取微震數(shù)據(jù)特征并標(biāo)注特征樣本類別；然后，構(gòu)建隨機(jī)森林模型以識(shí)別微震初至波；最后，將微震數(shù)據(jù)測(cè)試集樣本輸入隨機(jī)森林模型，得到每個(gè)樣本屬于某一類別的概率，將概率不小于0.5的第1個(gè)數(shù)據(jù)采樣點(diǎn)判定為微震初至波到時(shí)采樣點(diǎn)。實(shí)驗(yàn)結(jié)果表明該方法的準(zhǔn)確性優(yōu)于STA/LTA法。

1 隨機(jī)森林算法

隨機(jī)森林是一種有監(jiān)督的機(jī)器學(xué)習(xí)算法[13-15]，主要思想是將多個(gè)獨(dú)立的弱分類器即多棵決策樹進(jìn)行組合，形成一個(gè)強(qiáng)分類器。隨機(jī)森林算法原理：首先，從訓(xùn)練集中有放回地隨機(jī)抽取樣本，生成單棵決策樹，當(dāng)大量決策樹生成時(shí)，對(duì)其進(jìn)行組合構(gòu)成隨機(jī)森林；然后，通過驗(yàn)證集進(jìn)行多次驗(yàn)證，得到訓(xùn)練好的隨機(jī)森林模型；最后，將測(cè)試集樣本輸入訓(xùn)練好的隨機(jī)森林模型中，每一棵決策樹對(duì)其進(jìn)行判斷和分類，分別得到一個(gè)分類結(jié)果，采用最大投票法則確定待測(cè)樣本類型。

隨機(jī)森林的構(gòu)建過程如下。

(1) 隨機(jī)且有放回地抽樣。設(shè)訓(xùn)練集中樣本數(shù)量為N，從中隨機(jī)抽取1個(gè)樣本并放回，重復(fù)N次，得到的N個(gè)樣本組成決策樹的根節(jié)點(diǎn)。

(2) 決策樹節(jié)點(diǎn)分裂。為了降低每個(gè)特征之間的相關(guān)性，決策樹在節(jié)點(diǎn)分裂時(shí)從所有特征中隨機(jī)選取一部分，依次計(jì)算其基尼系數(shù)，將基尼系數(shù)較小的特征作為該節(jié)點(diǎn)的分裂準(zhǔn)則，根據(jù)分裂函數(shù)將當(dāng)前決策樹在該節(jié)點(diǎn)分為左子樹和右子樹。不斷重復(fù)上述過程，直到該節(jié)點(diǎn)所有樣本屬于同一類，停止分裂。

(3) 構(gòu)建隨機(jī)森林。執(zhí)行步驟(1)、(2)多次后產(chǎn)生大量決策樹，將其集成、組合構(gòu)成隨機(jī)森林。

2 微震初至波到時(shí)自動(dòng)拾取方法

基于隨機(jī)森林的微震初至波到時(shí)自動(dòng)拾取方法步驟如下。

(1) 微震數(shù)據(jù)特征提取及處理。當(dāng)微震事件發(fā)生時(shí)，微震記錄的主要特征會(huì)發(fā)生變化，尤其是初至波到來時(shí)刻，微震波振幅和能量(振幅的平方)明顯增大，相鄰時(shí)刻振幅比(后一時(shí)刻振幅絕對(duì)值與前一時(shí)刻振幅絕對(duì)值的比值)相應(yīng)增大，因此提取微震波振幅、能量及相鄰時(shí)刻振幅比作為特征，并對(duì)每個(gè)樣本進(jìn)行特征類別標(biāo)記。

(2) 構(gòu)建微震初至波到時(shí)自動(dòng)拾取隨機(jī)森林模型。定義提取的微震數(shù)據(jù)特征向量總數(shù)為T，第θ(θ=1,2，…，T)個(gè)特征向量為Uθ，Uθ對(duì)應(yīng)標(biāo)簽Lθ，微震數(shù)據(jù)特征向量與標(biāo)簽的集合為S={(Uθ,Lθ)}。將S隨機(jī)劃分為訓(xùn)練集H與驗(yàn)證集V，訓(xùn)練集和驗(yàn)證集樣本個(gè)數(shù)分別為α，β，1<β<α

從訓(xùn)練集H中隨機(jī)且有放回地抽取α次，生成α個(gè)樣本。從微震數(shù)據(jù)特征中隨機(jī)選取一部分，計(jì)算每個(gè)特征的基尼系數(shù)。設(shè)訓(xùn)練集H中有n個(gè)類別，第i(i=1,2，…，n)個(gè)類別所占比例為pi，則H對(duì)應(yīng)的基尼系數(shù)為

(1)

取基尼系數(shù)較小的特征作為決策樹當(dāng)前節(jié)點(diǎn)進(jìn)行分裂。當(dāng)計(jì)算得到的基尼系數(shù)為最小值0，即落在葉子節(jié)點(diǎn)中的所有微震數(shù)據(jù)特征屬于同一類別時(shí)，決策樹停止分裂。設(shè)訓(xùn)練集H被劃分為左子樹集合HL和右子樹集合HR，其基尼系數(shù)分別為GHL和GHR，則H的基尼系數(shù)可表示為

(2)

在決策樹的每個(gè)節(jié)點(diǎn)上隨機(jī)選取若干微震數(shù)據(jù)特征進(jìn)行測(cè)試，選取基尼系數(shù)最小的特征作為最優(yōu)分裂屬性劃分微震數(shù)據(jù)樣本。

將上述過程執(zhí)行γ次，生成γ棵決策樹，并對(duì)其進(jìn)行組合，生成隨機(jī)森林模型。將驗(yàn)證集V輸入隨機(jī)森林模型中，直到輸出結(jié)果滿足預(yù)設(shè)精度要求，得到最優(yōu)隨機(jī)森林模型。在構(gòu)建隨機(jī)森林模型過程中，振幅、能量和相鄰時(shí)刻振幅比3種特征的重要度分別為0.285 7,0.419 4,0.294 9。

(3) 預(yù)測(cè)類別。將微震測(cè)試集輸入最優(yōu)隨機(jī)森林模型，得到每個(gè)待測(cè)樣本屬于某一類別的概率。微震波判定為二分類問題，判定公式為

(3)

式中：L為微震特征向量對(duì)應(yīng)的標(biāo)簽；γ1,γ2分別為被分為正樣本決策樹、負(fù)樣本決策樹的數(shù)量。

若預(yù)測(cè)得到的概率不小于0.5，則被分為標(biāo)簽1是微震波，被分為標(biāo)簽0不是微震波，且預(yù)測(cè)概率不小于0.5的第1個(gè)采樣點(diǎn)判定為微震初至波到時(shí)采樣點(diǎn)。

基于隨機(jī)森林的微震初至波到時(shí)自動(dòng)拾取流程如圖1所示。

圖1 基于隨機(jī)森林的微震初至波到時(shí)自動(dòng)拾取流程Fig.1 Automatic picking process of microseismic first arrival based on random forest

3 實(shí)驗(yàn)及結(jié)果分析

3.1 實(shí)驗(yàn)數(shù)據(jù)

選取在煤礦井下巷道深孔中布置的微震傳感器采集的微震數(shù)據(jù)來驗(yàn)證基于隨機(jī)森林的微震初至波到時(shí)自動(dòng)拾取流程方法的有效性。提取每個(gè)微震數(shù)據(jù)的振幅、能量及相鄰時(shí)刻振幅比作為特征組成微震數(shù)據(jù)樣本，并對(duì)每個(gè)樣本進(jìn)行特征類別標(biāo)記。

若直接對(duì)微震初至波到時(shí)所對(duì)應(yīng)的樣本標(biāo)記為正樣本，其余標(biāo)記為負(fù)樣本，則會(huì)出現(xiàn)正負(fù)樣本分布不均衡的現(xiàn)象，導(dǎo)致訓(xùn)練效果不理想，因此需要擴(kuò)大正樣本數(shù)量。本文對(duì)微震初至波到時(shí)至微震波結(jié)束時(shí)刻對(duì)應(yīng)的所有樣本均標(biāo)記為正樣本，其余均標(biāo)記為負(fù)樣本，形成一個(gè)微震樣本數(shù)據(jù)集并存儲(chǔ)。該數(shù)據(jù)集中樣本數(shù)為8 000，將其隨機(jī)劃分為訓(xùn)練集與驗(yàn)證集，所含樣本數(shù)分別為5 600，2 400。

測(cè)試集來源于同一煤礦井下巷道深孔中的相同微震傳感器在不同時(shí)刻采集的8通道微震數(shù)據(jù)S1—S8，如圖2所示。采樣頻率為2 kHz。S1，S2，S7，S8為低信噪比微震數(shù)據(jù)，S3—S6為高信噪比微震數(shù)據(jù)。

圖2 8通道微震波形Fig.2 Eight-channel microseismic waveforms

3.2 實(shí)驗(yàn)設(shè)置

在Windows 10系統(tǒng)、基于Jupyter Notebook的Python3.7.1集成環(huán)境下，結(jié)合Scikit learn模塊進(jìn)行實(shí)驗(yàn)。設(shè)置STA/LTA法的短時(shí)窗為200 ms，長(zhǎng)時(shí)窗為800 ms，觸發(fā)閾值為2。隨機(jī)森林模型主要參數(shù)為決策樹的數(shù)量和最大深度，其對(duì)分類準(zhǔn)確率的影響如圖3所示。經(jīng)反復(fù)調(diào)整，選取2個(gè)參數(shù)分別為137，6。此時(shí)經(jīng)驗(yàn)證集驗(yàn)證，微震樣本分類準(zhǔn)確率達(dá)98.5%。

(a) 決策樹數(shù)量

(b) 決策樹最大深度

3.3 實(shí)驗(yàn)結(jié)果與分析

選取高信噪比微震數(shù)據(jù)S3和低信噪比微震數(shù)據(jù)S7作為測(cè)試數(shù)據(jù)。對(duì)S3拾取的微震初至波到時(shí)如圖4所示。可看出人工拾取方法、本文方法、STA/LTA法拾取的微震初至波到時(shí)分別為2 334.6，2 330.0，2 317.7 ms。以人工拾取結(jié)果為標(biāo)準(zhǔn)，本文方法與STA/LTA法的拾取誤差分別為4.6，16.9 ms。對(duì)S7拾取的微震初至波到時(shí)如圖5所示。可看出人工拾取方法、本文方法、STA/LTA法拾取的微震初至波到時(shí)分別為2 318.5，2 335.0，2 344.5 ms，本文方法與STA/LTA法的拾取誤差分別為16.5，26.0 ms?？梢姳疚姆椒ū萐TA/LTA法對(duì)微震初至波到時(shí)的拾取精度高。

(a) 人工拾取方法

(b) 本文方法

(c) STA/LTA法

(a) 人工拾取方法

(b) 本文方法

(c) STA/LTA法

對(duì)8個(gè)通道的微震數(shù)據(jù)分別采用人工拾取方法、STA/LTA法和本文方法拾取初至波到時(shí)，結(jié)果見表1。

表1 8通道微震初至波到時(shí)拾取結(jié)果統(tǒng)計(jì)Table 1 Result statistics of eight-channel microseismic first arrival picking

以人工拾取初至波到時(shí)為標(biāo)準(zhǔn)，由表1可知：本文方法和STA/LTA法對(duì)高信噪比微震數(shù)據(jù)的拾取準(zhǔn)確率均較高；對(duì)于低信噪比微震數(shù)據(jù)，本文方法的拾取準(zhǔn)確率較STA/LTA法高；對(duì)于8個(gè)通道的微震數(shù)據(jù)，STA/LTA法和本文方法對(duì)微震初至波到時(shí)拾取的平均誤差分別為39.5，23.1 ms，本文方法整體優(yōu)于STA/LTA法。

4 結(jié)語

提出了一種基于隨機(jī)森林的微震初至波到時(shí)自動(dòng)拾取方法，并成功應(yīng)用于煤礦巖石破裂微震信號(hào)自動(dòng)處理中。該方法與傳統(tǒng)STA/LTA法相比無需設(shè)定閾值，只需提取每個(gè)微震數(shù)據(jù)的振幅、能量及相鄰時(shí)刻振幅比這3個(gè)特征，并對(duì)每個(gè)樣本進(jìn)行特征類別標(biāo)注，減少了閾值由經(jīng)驗(yàn)判斷造成的誤差；該方法易于實(shí)現(xiàn)，無需設(shè)置特征函數(shù)與時(shí)窗長(zhǎng)度，消除了不同時(shí)窗長(zhǎng)度對(duì)微震初至波到時(shí)自動(dòng)拾取的影響，拾取精度明顯提升。