基于改進的HHT-SVM微震信號特征提取及分類識別研究*

熊璐偉,鐘曉陽,李庶林,葉龍珍,方 鑫,鄭宗檳

(1.廈門大學(xué) 建筑與土木工程學(xué)院,福建 廈門 361005;2.江西鐵山垅鎢業(yè)有限公司,江西 贛州 341000;3.自然資源部丘陵山地地質(zhì)災(zāi)害防治重點實驗室,福建 福州 350002)

0 引言

隨著淺層礦產(chǎn)資源的日益減少,我國礦山開采逐漸向深部發(fā)展,工程上面臨的沖擊地壓等巖石破裂災(zāi)害日益頻繁,嚴重威脅著礦山的生產(chǎn)安全。相較于傳統(tǒng)的監(jiān)測手段,微震監(jiān)測技術(shù)可監(jiān)測到巖體的微破裂活動并能動態(tài)地揭示巖體的損傷程度[1],且已在礦山、邊坡等巖體破裂災(zāi)害預(yù)警中取得不錯的成果[2-4],保障大量人員的生命財產(chǎn)安全,為各相關(guān)行業(yè)的生產(chǎn)安全提供有力支持。由于礦山工作環(huán)境復(fù)雜,爆破、電磁、機械設(shè)備、車輛等信號不時產(chǎn)生,很大程度上影響對微震信號的識別。其中電磁、機械設(shè)備、車輛信號具有較為固定的頻率、幅值和波形特征,容易剔除,而爆破致使巖石破裂產(chǎn)生的爆破信號在波形上與微震信號有著較強的相似性,難以準確辨識。對有效信號進行提取是微震監(jiān)測預(yù)警的關(guān)鍵一步,因此如何提高對微震信號、爆破信號的識別成功率具有重要意義。

近年來,在微震信號識別與分析領(lǐng)域,已有大量學(xué)者進行研究。裴琳[5]對隧洞開挖時搜集到的信號進行時頻和頻域分析,將信號分為巖石破裂信號、爆破信號及施工干擾信號3類。朱權(quán)潔等[6]對礦山爆破振動、巖石破裂信號進行小波包分解,求取指定特征頻帶上重構(gòu)信號的分形盒維數(shù),建立特征向量,使用SVM實現(xiàn)對信號的識別。陸菜平等[7]使用快速傅里葉變換對微震信號的幅頻特性與功率譜密度進行研究。何正祥等[8]運用梅爾倒譜系數(shù)法實現(xiàn)對微震事件信號的準確識別。Li等[9]使用多分形方法對微震與爆破信號進行處理,指出2個信號在簡單分形維數(shù),奇異指數(shù)范圍與多分形頻譜寬度存在一定差異,可作為信號識別特征。

以上許多學(xué)者從頻率角度入手對微震信號進行分類識別。本文對微震信號、爆破信號采用經(jīng)ICEEMDAN改進的HHT算法進行處理,分析2種信號的瞬時能量分布規(guī)律,尋找2種信號的差異點,并對差異點進行量化處理,提取相應(yīng)特征,構(gòu)建多維特征向量,建立經(jīng)功率譜熵特征加權(quán)與GSWOA算法改進后的SVM支持向量機微震信號識別網(wǎng)絡(luò),研究結(jié)果可為辨識微震信號與爆破信號提供1種新的方法與思路。

1 基本理論

1.1 支持向量機

支持向量機(support vector machine,SVM)算法在解決小樣本、非線性和高維的二分類問題方面有較大優(yōu)勢,已被廣泛應(yīng)用于信號識別、文本分類等領(lǐng)域。其原理是通過使用核函數(shù)將低維特征空間轉(zhuǎn)化為高維度特征空間,尋求能夠?qū)颖具M行分類的最優(yōu)超平面。假設(shè)訓(xùn)練樣本集為(xi,yi),i=1,2,…,l,xi∈Rn,yi∈R,SVM通過映射核函數(shù)φ(x)將xi映射到高維特征空間,超平面表達式如式(1)所示:

ω(x)=wφ(x)+b

(1)

式中:w為權(quán)重向量;φ(x)為SVM的核函數(shù);b代表原點與超平面間的距離。

上述求解超平面問題可轉(zhuǎn)化為求解式(1)的最小值問題,如式(2)所示:

(2)

式中:C為懲罰因子;l為樣本集數(shù)量;Lε為損失函數(shù)。

為了考慮擬合誤差,引入松弛變量,使得優(yōu)化目標變?yōu)槿缡?3)所示:

(3)

式中:ξi為松弛變量。

為了求解式(3)的最優(yōu)解,引入拉格朗日函數(shù)將上述問題轉(zhuǎn)換成對偶問題,則可得到最終的超平面表達式,如式(4)所示:

(4)

本文選用徑向基核函數(shù)(RBF),SVM的核心是求解1個二次規(guī)劃問題,使得目標函數(shù)f(xi)最小,如式(5)所示:

(5)

式中:f(xi)為目標函數(shù);標簽yi=±1。

1.2 功率譜熵特征加權(quán)算法

支持向量機在對樣本進行分類時,將樣本的每一維特征都賦予相同的權(quán)重,但事實上,樣本的各維度特征對分類效果的貢獻大小不一定一致[10]。若對分類效果貢獻較低的特征參數(shù)賦予與對分類效果貢獻較高的特征參數(shù)相同的權(quán)重,則必然會影響最終的分類效果。因此本文使用信息學(xué)中熵的概念對特征向量賦予權(quán)重。在信息學(xué)中,熵能有效度量隨機變量的不確定性程度,能反映隨機變量包含的信息量。本文采用功率譜熵對樣本特征向量賦予權(quán)重,若樣本特征的功率譜熵越大,則說明這一維度特征的不確定性越高,對分類效果的貢獻較小,應(yīng)分配較小的權(quán)重;反之,則應(yīng)分配較大的權(quán)重。

若訓(xùn)練樣本集數(shù)量為N,特征維度為d,則特征向量矩陣為Φ=[Φ1,Φ2,…,Φd],其具體表達式如式(6)所示:

(6)

樣本集的第k維度特征可表達為Φk,對樣本的第k維度特征的功率譜熵進行計算,計算公式如式(7)所示:

(7)

式中:H(k)為第k維度特征的功率譜熵;pl為離散序列信號的第l頻段的能量占總能量的比重大小。

H=[H(1),H(2),…,H(d)]表示d維特征的功率譜熵,第k個維度的特征權(quán)重計算公式如式(8)所示:

(8)

式中:ωq(k)為第k個維度的特征權(quán)重。

對樣本特征向量加權(quán)后,SVM的優(yōu)化問題轉(zhuǎn)化為如式(9)所示:

(9)

式中:ωq為特征權(quán)重向量。

1.3 全局搜索策略的鯨魚算法

根據(jù)鯨魚的捕食習(xí)慣,文獻[11]中提出WOA算法,主要包括包圍獵物、旋轉(zhuǎn)搜尋及隨機搜尋3部分。但由于此算法存在收斂速度慢且易陷入局部最優(yōu)解的問題,本文采用劉磊等[12]提出的全局搜索策略的鯨魚算法(global search whale optimization algorthm,GSWOA),其原理如下:

1)包圍獵物

鯨魚將目前搜捕獵物的信息進行共享,然后鯨魚向離獵物最近的鯨魚靠近,逐漸縮小包圍圈,實現(xiàn)對獵物的包圍,此時鯨魚位置更新表達式如式(10)所示:

(10)

式中:t為迭代搜尋次數(shù);X為鯨魚位置;X*是全局最優(yōu)位置;A和B為系數(shù)矩陣;D為權(quán)重系數(shù)。

2)自定義權(quán)重與變螺旋搜尋

在對鯨魚位置進行更新時,增加自適應(yīng)權(quán)重來減少處于最佳位置鯨魚對其鯨魚位置更新的影響,減少算法陷入局部最優(yōu)的風險;同時為了使鯨魚能以多種方式搜尋獵物,引入變螺旋搜尋方法,自適應(yīng)權(quán)重及位置更新公式如式(11)～(12)所示:

(11)

(12)

式中:w(t)為自適應(yīng)權(quán)重;Xrand(t)為隨機鯨魚位置;p為[0,1]均勻分布隨機數(shù)。

3)最優(yōu)鄰域擾動

為了提高收斂速度以及避免陷入局部最優(yōu)解,引入了最優(yōu)鄰域擾動策略,在當前最優(yōu)值附近進行隨機搜尋,尋找更佳的全局值,其定義如式(13)所示:

(13)

按照式(10)～(13)完成種群位置更新,在預(yù)定迭代次數(shù)下求解最佳適應(yīng)度對應(yīng)的參數(shù)值。

1.4 改進的希爾伯特-黃變換

希爾伯特-黃變換(hilbert-huang transform,HHT)是目前分析爆破、振動信號的主要方法。其主要由經(jīng)驗?zāi)B(tài)分解(empirical mode decomposition,EMD)與Hilbert變換組成。其中EMD算法分解的本征模態(tài)分量(intrinsic mode function,IMF)存在模態(tài)混疊現(xiàn)象,導(dǎo)致分解效果不理想。改進自適應(yīng)完備經(jīng)驗?zāi)Ｊ椒纸?improved complete ensemble empirical mode decomposition with adaptive noise,ICEEMDAN)算法通過改進自適應(yīng)分解過程,大幅減少了信號分解后產(chǎn)生的殘余分量,獲得更加良好的重構(gòu)信號,一定程度上改善了模態(tài)混疊問題。因此本文采用ICEEMDAN算法[13]對HHT算法進行改進。ICEEMDAN算法的計算流程可分為以下5個步驟:

1)基于原始信號S構(gòu)造N個含可控噪聲信號,如式(14)所示:

S(i)=S+β0E1(ω(i))i=1,2,…,N

(14)

式中:S(i)表示第i個構(gòu)造信號;β0為第1次分解時信號的噪聲標準偏差;β0的取值范圍一般為[0.1,0.3],本文取0.2;ω(i)表示第i個零均值單位方差白噪聲;E1(·)表示第1個IMF算子。

2)計算S(i)的局部均值并求平均,得到第1個殘余分量,計算公式如式(15)所示:

(15)

式中:M(·)為局部均值函數(shù)。

3)將原信號S減去第1個殘余分量r1得到第1個模態(tài)(k=1),如式(16)所示:

(16)

4)計算第k個模態(tài)(k≥2),計算公式如式(17)～(18)所示:

(17)

(18)

5)返回4),計算第k+1個模態(tài),直到殘余分量不可再分或迭代次數(shù)滿足終止條件,本文將迭代次數(shù)設(shè)置為100次。

使用互信息法對分解后的模態(tài)進行篩選,計算各IMF分量與原信號的互信息量MIi,并定義平均互信息量MIm,計算公式如式(19)所示:

(19)

式中:K為ICEEMDAN算法的分解個數(shù)。

若MIi>MIm,則認為該分量為有效分量。

希爾伯特變換主要原理為先將信號使用ICEEMDAN算法分解成多個IMF分量,使用互信息法進行篩選,對有效分量進行Hilbert變換得到希爾伯特譜,其表達式如式(20)所示:

(20)

式中:Re表示取實部進行運算;F為分解信號頻率;T為各信號采樣點對應(yīng)時間;n為分解信號個數(shù);ai(T)表示瞬時幅值;φi(T)表示瞬時頻率。

將Hilbert譜的平方對頻率進行積分,可得到相應(yīng)的Hilbert瞬時能量譜,如式(21)所示:

(21)

式中:IE(T)為希爾伯特譜瞬時能量。

2 基于改進HHT微震信號瞬時能量特征分析

2.1 工程背景簡介

黃沙鎢礦位于贛州市于都縣,地處南嶺東西構(gòu)造帶大余——會昌隆起帶與新華復(fù)系于山隆起帶的反接復(fù)合區(qū),礦區(qū)構(gòu)造復(fù)雜,裂隙發(fā)育。由于黃沙礦區(qū)開采時間長,礦脈密集,采深大,坑場體積大,特別是礦脈密集的采空區(qū)或大采區(qū),有利于形成直通地表的大型倒楔形不穩(wěn)固結(jié)構(gòu)體,使通過采空區(qū)的斷層不穩(wěn)定性加大,導(dǎo)致地壓的增強,致使局部區(qū)域出現(xiàn)冒頂、片幫等問題。針對日益嚴重的地壓問題,于2023年引入具有大范圍、立體監(jiān)測特點的實時在線多通道微震監(jiān)測系統(tǒng)。微震系統(tǒng)采用加速度型傳感器,其響應(yīng)頻率范圍為50～5 000 Hz,數(shù)據(jù)采集儀的最高采集頻率為30 kHz。從微震系統(tǒng)數(shù)據(jù)庫內(nèi),選取典型的爆破信號與微震信號進行分析。

2.2 信號瞬時能量變化規(guī)律分析

通過改進的HHT變換對拾取波形進行處理分析,揭示信號瞬時能量的變化規(guī)律特征,典型微震、爆破信號的時域及瞬時能量圖如圖1～2所示。

圖1 典型微震信號時域及瞬時能量圖Fig.1 Time domain and transient energy diagram of typical microseismic signal

通過分析峰前能量的變化規(guī)律,爆破信號的瞬時能量在爆破事件開始后迅速達到能量峰值,而由于微震信號來源于巖體內(nèi)部微破裂的產(chǎn)生、發(fā)展、貫通過程中產(chǎn)生的應(yīng)力波[14],導(dǎo)致微震信號從事件開始發(fā)生至能量達到峰值的過程相對緩慢,即事件孕育過程速度較慢。分析峰后瞬時能量的變化規(guī)律,從峰值能量點后瞬時能量線的下降斜率可以明顯看出爆破信號在達到能量峰值后衰減的速度超過微震信號的峰后衰減速度。從圖2可以看出爆破信號的后期部分衰減波動較大,這是由于爆破信號往往會伴隨著2～3個地震子波,且爆破信號的地震子波能量較高,對爆破信號后期的幅值衰減有著較大的影響,導(dǎo)致爆破信號具有持續(xù)時間相對較長的特點[15]。

圖2 典型爆破信號時域及瞬時能量圖Fig.2 Time domain and transient energy diagram of typical blasting signal

2.3 信號瞬時能量變化規(guī)律研究及特征提取

當信號瞬時能量初次達到1%峰值能量時標志著事件開始,當信號能量最后一次衰減至1%峰值能量時標志著事件的終止,如圖1(b)所示,A點為事件的起始點,C點為事件的終止點。通過峰值能量點B將整個事件(AC段)劃分為事件孕育階段(AO段)與事件衰減階段(OC段)。

通過將信號事件按瞬時能量劃分為3個階段進行分析,發(fā)現(xiàn)微震信號、爆破信號在孕育速度、持續(xù)時間、峰后衰減速度具有一定差異。對差異點進行定量分析,將事件孕育階段時間與整個事件持續(xù)時間之比記作量化結(jié)果1,對事件的孕育速度進行描述;將整個事件持續(xù)時間記作量化結(jié)果2;將事件衰減階段按時間分為20個區(qū)域,并統(tǒng)計每個區(qū)域內(nèi)的累計能量值即瞬時能量對該區(qū)域時間進行積分,將峰值能量后的6個區(qū)域的累計能量值進行記錄,并將6個區(qū)域間的累計能量的衰減率記作量化結(jié)果3～7,對信號峰后衰減速度進行描述,衰減率計算公式如式(22)所示:

(22)

式中:Ei表示第i個部分的區(qū)間累計能量。

從微震監(jiān)測系統(tǒng)中隨機抽取微震信號與爆破信號各200組,對這400組信號采用上述量化方法進行處理后,得到各量化結(jié)果概率密度分布如圖3所示。

3 改進的SVM信號識別網(wǎng)絡(luò)

3.1 改進的SVM信號識別網(wǎng)絡(luò)建立

基于功率譜熵特征加權(quán)與GSWOA算法改進的SVM信號識別網(wǎng)絡(luò)的具體實現(xiàn)步驟如下:

1)構(gòu)造特征向量,建立訓(xùn)練、測試數(shù)據(jù)集,設(shè)置信號類型標簽。選取量化結(jié)果1～7構(gòu)建7維特征向量。從信號數(shù)據(jù)庫中隨機抽取200組微震信號樣本、200組爆破信號樣本,將60%樣本作為訓(xùn)練集,40%樣本作為測試數(shù)據(jù)集。將微震信號類型標簽設(shè)置為1,爆破信號類型標簽設(shè)置為2。

2)數(shù)據(jù)預(yù)處理。對訓(xùn)練、測試集數(shù)據(jù)使用[0,1]區(qū)間歸一化方法進行預(yù)處理,計算訓(xùn)練樣本集特征向量的功率譜熵及權(quán)重,計算結(jié)果如表1所示。

表1 特征參數(shù)功率譜熵值及權(quán)重Table 1 Entropy values and weights of power spectrum of feature parameters

3)SVM網(wǎng)絡(luò)參數(shù)選取。SVM的分類器選用徑向基核函數(shù)RBF,經(jīng)GSWOA算法對RBF中懲罰系數(shù)C,核函數(shù)參數(shù)g進行優(yōu)化。設(shè)置鯨魚數(shù)量為30,迭代次數(shù)100,懲罰系數(shù)C與核函數(shù)參數(shù)g的臨界區(qū)域分別設(shè)置為[0,100]、[0,2.00]。通過100次迭代計算后,得到優(yōu)化參數(shù)組合[54,0.21]。

3.2 分類識別結(jié)果

將經(jīng)改進HHT提取的訓(xùn)練集特征向量及其權(quán)重與經(jīng)GSWOA算法尋優(yōu)后得到的C與g的優(yōu)化參數(shù)組合[54,0.21]輸入SVM識別網(wǎng)絡(luò)中,隨后輸入測試集數(shù)據(jù)并進行分類識別,最終的綜合識別成功率為96.250%,微震、爆破信號具體識別情況如圖4所示。

圖4 信號識別效果圖Fig.4 Signal identification effect

為了驗證本文提出的基于功率譜熵特征加權(quán)算法與GSWOA算法改進后的SVM信號識別網(wǎng)絡(luò)的可行性與優(yōu)越性,將其與傳統(tǒng)SVM、GSWOA-SVM、基于功率譜熵特征加權(quán)-SVM方法對比,未進行參數(shù)優(yōu)化的SVM中C與g參數(shù)組合取[10,0.10],分類識別結(jié)果如表2所示。結(jié)果表明改進后的HHT的識別結(jié)果整體優(yōu)于傳統(tǒng)HHT算法,說明改進HHT算法能更加充分地提取微震、爆破信號特征。對于同一特征向量,經(jīng)功率譜熵特征加權(quán)算法及GSWOA算法改進的SVM與改進HHT相結(jié)合的方法對的識別成功率達到96.250%,驗證了本文建立的信號識別網(wǎng)絡(luò)的優(yōu)越性與可行性,表明本文所提方法能有效識別微震與爆破信號。

表2 分類識別結(jié)果Table 2 Classification and identification results

4 結(jié)論

1)使用ICEEMDAN改進HHT算法對微震信號與爆破信號進行處理,分析信號的瞬時能量變化規(guī)律,并對400組信號進行定量研究,提取信號事件孕育速度、事件持續(xù)時間及峰后衰減速度的特征并構(gòu)建多維特征向量。

2)通過定量研究分析,構(gòu)建多維特征向量的各量化結(jié)果分類識別效果具有差異,使用功率譜熵特征加權(quán)算法對多維特征向量賦予權(quán)重,能有效提高對微震信號的識別成功率。

3)改進的HHT分類效果優(yōu)于傳統(tǒng)HHT算法;使用改進HHT對信號進行特征提取并構(gòu)建特征向量,輸入至經(jīng)功率譜熵特征加權(quán)與GSWOA算法改進的SVM信號識別網(wǎng)絡(luò)中,最終測試集的綜合識別成功率為96.250%,表明本文提出的改進算法是有效可行的,能夠?qū)ξ⒄稹⒈菩盘栠M行有效分類識別,為后續(xù)利用微震信號進行監(jiān)測預(yù)警奠定基礎(chǔ)。