孟子軒 程 巍 張?zhí)煊?呂 君 滕鵬曉
(1 中國科學(xué)院聲學(xué)研究所 北京 100190)
(2 中國科學(xué)院噪聲與振動重點(diǎn)實(shí)驗(yàn)室 北京 100190)
(3 中國科學(xué)院大學(xué) 北京 100049)
次聲是指頻率低于20 Hz 的聲信號。自然環(huán)境與人類社會中廣泛存在著次聲信號,許多物理現(xiàn)象在發(fā)生和發(fā)展過程中都會伴隨低頻次聲信號的產(chǎn)生,例如自然活動中的地震、臺風(fēng)、閃電、火山噴發(fā)和海嘯以及人類活動中的核和化學(xué)爆炸、火箭發(fā)射、飛機(jī)起飛等事件[1]。從災(zāi)害預(yù)防的角度看,對次聲信號的識別分類能夠起到預(yù)警的作用;在軍事對抗領(lǐng)域,通過次聲信號來識別敵對方的軍事活動對國防安全具有十分重要的意義。
自1996 年全面禁止核試驗(yàn)條約組織(Comprehensive nuclear-test-ban treaty organization,CTBTO)成立以來,次聲成為國際監(jiān)測系統(tǒng)(International monitoring system,IMS)所使用的4種主要監(jiān)測核爆的手段之一,利用機(jī)器學(xué)習(xí)方法對次聲信號進(jìn)行識別分類的研究也由此展開[2]。對于次聲信號的識別分類任務(wù)而言,由于樣本數(shù)量較少,因此研究的關(guān)鍵問題在于信號的特征提取[3?4]。吳涢暉[5]采用了8 種不同的特征提取方式對化學(xué)爆炸、閃電和臺風(fēng)3 類事件進(jìn)行處理,對比了支持向量機(jī)(Support vector machines,SVM)、BP神經(jīng)網(wǎng)絡(luò)、長短時(shí)神經(jīng)網(wǎng)絡(luò)和卷積神經(jīng)網(wǎng)絡(luò)(Convolutional neural network,CNN) 4 種分類器的分類性能,結(jié)果表明SVM的識別性能最好。同時(shí)指出,其構(gòu)建的分類流程對人工設(shè)計(jì)的特征有較高的要求,需要研究人員對各種信號的特征進(jìn)行深入的研究和挑選,以找出區(qū)分度較大的特征。
譚笑楓等[6]以數(shù)據(jù)驅(qū)動為出發(fā)點(diǎn),采用CNN進(jìn)行特征提取、模型訓(xùn)練和識別分類,以簡化特征設(shè)計(jì)過程,其提出的方法在CTBTO 提供的化爆與地震兩類次聲數(shù)據(jù)集上達(dá)到了82.72%的準(zhǔn)確率,但該方法對數(shù)據(jù)量的要求較高,難以應(yīng)用于小樣本數(shù)據(jù)集的情況。本文從矩陣分解的角度,考慮采用淺層的模型對次聲信號進(jìn)行特征提取,以適用于小樣本場景下的次聲信號分類任務(wù)。
戴翊靖等[7]考慮到次聲信號特性與樣本量小的特點(diǎn),采用非負(fù)矩陣分解(Non-negative matrix factorization,NMF)進(jìn)行次聲降噪。該方法的研究對象為一段混合的次聲信號,將預(yù)訓(xùn)練的平穩(wěn)噪聲作為監(jiān)督項(xiàng),通過NMF 方法將信號中的目標(biāo)部分與噪聲部分分離,然后恢復(fù)目標(biāo)信號。其目的在于獲得較高信噪比的目標(biāo)信號,以便于后續(xù)的次聲監(jiān)測任務(wù)。本文聚焦于分類的具體應(yīng)用,對于大氣次聲波的產(chǎn)生和傳播機(jī)理研究以及次聲事件的影響和防治策略都有重要意義。特別是本研究使用的NMF 方法處理結(jié)果是在可靠的先驗(yàn)信息標(biāo)注的多條數(shù)據(jù)庫上獲得的,其目的在于提取目標(biāo)信號的基本組成部分,分別作為不同類信號的特征。特征提取后使用SVM 等分類器完成識別分類任務(wù),是次聲監(jiān)測的重要環(huán)節(jié),也是大氣次聲學(xué)的一項(xiàng)基礎(chǔ)性研究。
NMF 是一種廣泛應(yīng)用于圖像識別、語聲增強(qiáng)以及聲事件識別等領(lǐng)域的算法[8?10],其基本框架由Lee 等[11]提出。定義非負(fù)矩陣V ∈Rf×t,此處非負(fù)的含義是指V中的任一元素Vij≥0。NMF 的目的是希望得到兩個(gè)非負(fù)矩陣W與H,同時(shí)保證二者乘積與一個(gè)V盡可能接近:
在本文中,V為原信號經(jīng)短時(shí)傅里葉變換(Short-time Fourier transform,STFT)后得到的時(shí)頻圖,而W ∈Rf×d與H ∈Rd×t則表示經(jīng)訓(xùn)練后從數(shù)據(jù)中學(xué)習(xí)到的特征,其具體含義往往與實(shí)際問題相關(guān)。W可以認(rèn)為是V中的基本組成部分,被稱為字典矩陣,其列向量被稱為字典原子;而H被稱為激活矩陣,表示在V中這些部分在相應(yīng)時(shí)間點(diǎn)上的線性組合計(jì)權(quán),其行向量被稱為字典原子的激活系數(shù);參數(shù)d則表示字典矩陣中字典原子的個(gè)數(shù)[12?13],此時(shí)有
式(2)中,wi與hi分別為字典原子與相應(yīng)的激活系數(shù)。由于通常設(shè)定d ?min(f,t),即只用很少的字典原子來描述原信號,因此只有在W包含了原信號中最主要的組成時(shí),才會使得式(1)成立[12],本文將這些表示原信號基本部分的字典原子作為特征輸入,采用SVM等分類器進(jìn)行識別分類任務(wù)。
W和H可通過最小化V與WH之間的距離度量函數(shù)得到,即求解如下的優(yōu)化問題[11]:
式(3)中,D表示V與WH之間的距離度量函數(shù),其定義為
式(5)中,i=1,2,···,f;j=1,2,···,t。應(yīng)用梯度下降法對式(5)進(jìn)行迭代求解,J(W,H)分別對Wil與Hlj的偏導(dǎo)數(shù)為
其中,l=1,2,···,d,T表示矩陣轉(zhuǎn)置。則基于梯度下降法的迭代規(guī)則,由式(6)和式(7)可得
其中,λil和μlj是迭代步長。乘性更新法則對步長進(jìn)行限制,規(guī)定其取值為
將式(10)與式(11)分別帶入式(8)與式(9)中,得到歐幾里得距離度量下的乘性更新法則:
式(12)與式(13)的非增性在文獻(xiàn)[15]中得到了證明??梢钥闯?,采用式(12)與式(13)迭代求解式(3)時(shí)只進(jìn)行乘法和加法運(yùn)算,從而保證了迭代過程及結(jié)果的非負(fù)性。
對于不同的距離度量函數(shù)而言,其乘性更新法則可采用統(tǒng)一的形式給出[13]:
式(14)中,β的不同取值代表采用不同的距離度量函數(shù),?表示對矩陣元素的乘法,同時(shí)指數(shù)運(yùn)算和除法運(yùn)算均表示對元素的運(yùn)算。
一個(gè)完整的識別分類系統(tǒng)如圖1 所示。分類模型的選擇主要包括傳統(tǒng)機(jī)器學(xué)習(xí)算法與深度神經(jīng)網(wǎng)絡(luò)兩類。SVM 是一種廣泛應(yīng)用于次聲信號識別分類領(lǐng)域的機(jī)器學(xué)習(xí)模型,由于其分類面完全由支持向量確定,因此在小樣本場景下有較好的表現(xiàn)性能[16]。
圖1 識別分類系統(tǒng)框圖Fig.1 Block diagram of the classification system
SVM 的基本模型是定義在特征空間上的間隔最大的二分類線性分類器[14],即求解如下的優(yōu)化問題:
式(15)中,w為分類超平面的法向量,xi為第i個(gè)輸入,yi為輸入對應(yīng)的標(biāo)簽,b則表示超平面的偏置量。式(15)是在數(shù)據(jù)集線性可分的情況下得到的,又被稱為線性可分SVM。對于一般的線性不可分?jǐn)?shù)據(jù)集,通常是對式(15)添加懲罰項(xiàng),允許部分樣本不滿足約束,則式(15)可改寫為
式(16)中,?(·)表示損失函數(shù)。若采用合頁損失函數(shù),同時(shí)引入松弛變量ξi≥0,可以得到線性SVM模型,如下:
采用拉格朗日乘子法,可以得到式(17)的對偶問題,如下:
式(18)中,αi表示拉格朗日乘子。式(15)與式(18)到的SVM 分類面均為線性分類面,通過采用核技巧,可以得到更為一般的SVM形式,如下:
式(19)中,κ(·)表示核函數(shù)。若選用線性核函數(shù),則式(19)退化為式(18)。幾種常用的核函數(shù)包括線性核、高斯核、多項(xiàng)式核等,需在實(shí)際使用時(shí)進(jìn)行選擇。以上提出的SVM 均為二分類模型,對于K種類別的多分類問題,目前最常用的方法是分別構(gòu)建K個(gè)獨(dú)立的SVM[17]:當(dāng)訓(xùn)練第k個(gè)模型時(shí),使用當(dāng)前類別yk作為正樣本數(shù)據(jù),而將其余的K?1 個(gè)類別作為負(fù)樣本。
由于所提特征是一種二維數(shù)據(jù),因此同時(shí)使用CNN 進(jìn)行分類。LeNet-5 是一種經(jīng)典的CNN 結(jié)構(gòu)[18?19],其網(wǎng)絡(luò)結(jié)構(gòu)如圖2所示。
圖2 CNN 結(jié)構(gòu)圖Fig.2 CNN structure diagram
LeNet-5 是一種較為淺層的CNN,其輸入經(jīng)過兩層卷積層、兩層池化層后與全連接神經(jīng)網(wǎng)絡(luò)相連,最后輸出其類別概率。為保證泛化能力,每層網(wǎng)絡(luò)均添加正則化項(xiàng),并采用Dropout進(jìn)行隨機(jī)失活。
中國科學(xué)院聲學(xué)研究所計(jì)劃在全國范圍內(nèi)建設(shè)廣域多臺陣次聲監(jiān)測網(wǎng)絡(luò),現(xiàn)已在新疆、遼寧、云南等地建設(shè)了固定式次聲臺陣,并通過機(jī)動式次聲探測站對酒泉、文昌等發(fā)射基地的次聲信號進(jìn)行收集。通過逐步多通道相關(guān)方法等陣列處理算法得到每個(gè)信號的相速度、入射角和聲源位置等聲學(xué)參數(shù),并結(jié)合聲源信息的驗(yàn)證,對于信號的類別可以進(jìn)行準(zhǔn)確的標(biāo)注,經(jīng)過準(zhǔn)確標(biāo)注的實(shí)測數(shù)據(jù)集為本研究奠定了重要的基礎(chǔ),有利于推動我國次聲監(jiān)測技術(shù)的發(fā)展,對于自然災(zāi)害預(yù)警以及國防建設(shè)等領(lǐng)域具有重大意義。
本文使用的數(shù)據(jù)為通過實(shí)地布陣自行采集的次聲時(shí)域信號,數(shù)據(jù)集大小為105 條,包括4 種類型的信號,其中爆炸信號17 條,地震信號22 條,閃電信號41條,再入信號25條。本文采用原信號進(jìn)行過STFT 后的時(shí)頻圖作為輸入。4 類信號經(jīng)STFT 后的時(shí)頻圖如圖3所示。
圖3 4 類信號時(shí)頻圖Fig.3 The spectrograms of four kinds of signals
圖3 所有標(biāo)注出的紅色方框區(qū)域均為信號部分。就信號的頻率維特征而言,爆炸事件與地震事件的主要頻率分布在5 Hz 以下,其中地震事件在該頻帶中的分布比較均勻,而爆炸事件則分布在2~4 Hz,相對集中;閃電事件的頻帶較寬,在20 Hz以上仍有部分信號,主要部分明顯地集中在5 Hz 和15 Hz 左右,分成了兩個(gè)部分;再入事件的頻率則分布在15 Hz 以下,主要部分在5 Hz 以下,持續(xù)時(shí)間較短。從時(shí)頻圖上可以看出各信號之間存在明顯差異。由于各個(gè)信號均為實(shí)際采集信號,則不可避免地混有了噪聲成分,這干擾了識別分類任務(wù)。為了找出目標(biāo)信號的主要部分,同時(shí)對原始數(shù)據(jù)進(jìn)行降維,需要進(jìn)行特征提取工作。圖3中呈現(xiàn)的4類次聲波頻譜隨時(shí)間的變化,和持續(xù)時(shí)間的長短都有顯著的差異,其形成機(jī)制復(fù)雜,本研究尚未涉及,暫時(shí)擱置,留待進(jìn)一步的研究。
本文采用NMF 對圖3 所示的4 類信號時(shí)頻圖進(jìn)行特征提取。在NMF 中,字典原子個(gè)數(shù)d是需要預(yù)先設(shè)置的參數(shù),可采用經(jīng)驗(yàn)的參考取值[20],但在實(shí)際應(yīng)用時(shí)仍需進(jìn)行調(diào)整。首先對圖3 所示的4 類信號進(jìn)行NMF 分解,采用歐幾里得距離作為代價(jià)函數(shù),迭代次數(shù)設(shè)置為200 次,在d=8 時(shí)得到的爆炸信號收斂曲線如圖4所示。
圖4 爆炸信號收斂曲線Fig.4 The convergence curves of four kinds of explosion signals
圖4 中黑實(shí)線表示全部事件的度量函數(shù)值平均值,上下兩條虛線則分別表示全部事件度量函數(shù)值的上下限。為了保證特征能夠充分反映原始信號,要求選取的迭代次數(shù)使得信號的距離度量函數(shù)充分收斂。從圖4 中可以看出,經(jīng)過50 次迭代后度量函數(shù)已趨于收斂,同時(shí)經(jīng)實(shí)驗(yàn)驗(yàn)證,繼續(xù)增大迭代次數(shù)對于分類準(zhǔn)確率的影響不大,綜合計(jì)算效率與分類準(zhǔn)確率,本文選取迭代次數(shù)為200。4 類信號經(jīng)NMF分解得到的W矩陣與H矩陣如圖5所示。
圖5 4 類事件NMF 分解結(jié)果Fig.5 NMF decomposition results of four kinds of events
圖5 中左側(cè)為信號的字典矩陣,右側(cè)為激活矩陣。結(jié)合圖3 與圖5,可以看出字典矩陣是對信號時(shí)頻譜的一種降維表示,不同類之間的字典矩陣有明顯差異,因此本文選取W作為信號的特征向量進(jìn)行分類實(shí)驗(yàn)。
本文使用的數(shù)據(jù)集中訓(xùn)練集與測試集的比例為7 : 3,由于數(shù)據(jù)量有限,因此不再設(shè)計(jì)驗(yàn)證集,而是在訓(xùn)練時(shí)采用四折交叉驗(yàn)證的方式進(jìn)行模型選擇。識別分類實(shí)驗(yàn)包括兩部分,分別是基于經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical mode decomposition,EMD)和NMF 的識別分類實(shí)驗(yàn)。分類實(shí)驗(yàn)在AMD?4600H平臺上進(jìn)行,操作系統(tǒng)為Windows10,所用軟件為Python3.6.8,CNN 的開發(fā)框架為Tensorflow2.4.0,SVM模型由Sklearn模塊提供。
2.3.1 基于EMD的分類過程
基于EMD 的特征提取主要是對分解得到的各階本征模態(tài)函數(shù)(Intrinsic mode function,IMF)分量進(jìn)行處理[3,5],可選擇的處理方式包括計(jì)算分量各階矩、能量、信息熵、波形特征、分量比等,經(jīng)處理后的各分量仍可繼續(xù)提取其能量、信息熵或波形特征等特征。本文采用時(shí)域能量、時(shí)域熵、EMD能量、EMD 熵、EMD 奇異值以及希爾伯特邊際譜(Hilbert marginal spectrum,HMS)[21]4 種特征作為對比。以HMS 為例,對識別分類過程進(jìn)行說明。提取HMS 時(shí)首先對時(shí)域信號進(jìn)行EMD,將得到的IMF進(jìn)行希爾伯特變換,構(gòu)造出原信號的解析信號,從而得到原信號的希爾伯特譜,HMS 即為希爾伯特譜的時(shí)間維積分結(jié)果,反映了瞬時(shí)頻率的時(shí)域幅值累加。對得到的HMS進(jìn)行主成分分析(Principal component analysis,PCA),取前30維作為SVM輸入進(jìn)行分類。本文采取隨機(jī)優(yōu)化的方式對SVM 進(jìn)行參數(shù)選擇,以相同條件進(jìn)行100次分類實(shí)驗(yàn),在測試集上的分類結(jié)果如圖6所示。
圖6 HMS-SVM 測試集分類結(jié)果Fig.6 Test set classification results of HMS-SVM
圖6中橫坐標(biāo)為HMS-SVM在測試集上的分類準(zhǔn)確率,縱坐標(biāo)則表示在100 次實(shí)驗(yàn)中不同準(zhǔn)確率出現(xiàn)的次數(shù),藍(lán)色柱狀圖為實(shí)驗(yàn)數(shù)據(jù)記錄,紅色曲線為使用正態(tài)分布擬合的結(jié)果。從圖6 中可以看出,HMS-SVM 在測試集上的分類準(zhǔn)確率最大值為87.5%,平均值為69.76%,方差為7.65%,多數(shù)分類結(jié)果集中在60%~80%之間。本文還提取了時(shí)域能量、時(shí)域熵、EMD能量、EMD熵等4 種特征,在上述實(shí)驗(yàn)條件下進(jìn)行分類實(shí)驗(yàn),分別使用SVM 與一維CNN作為分類器,其結(jié)果見表1。
表1 5 種特征分類結(jié)果Table 1 Five feature classification results
從表1 中可以看出,在本文所使用的數(shù)據(jù)集中,5 種特征中最高的準(zhǔn)確率為EMD 熵特征,達(dá)到了68.71%。進(jìn)行EMD 后使用熵作為特征的準(zhǔn)確率提升較高,而使用能量作為特征的準(zhǔn)確率則有所下降。在進(jìn)行實(shí)驗(yàn)前無法確定最優(yōu)的特征提取方式,需要根據(jù)分類結(jié)果對特征向量進(jìn)行設(shè)計(jì)以獲取最佳特征。
2.3.2 基于NMF的分類過程
使用NMF 進(jìn)行特征提取時(shí),除字典原子個(gè)數(shù)外,還需要考慮距離度量函數(shù)的選取。本文比較了歐幾里得距離(β=2)、廣義K-L 散度(β=1)和IS散度(β=0)三種距離度量函數(shù)下不同字典原子個(gè)數(shù)的分類結(jié)果。分別采用SVM與CNN作為分類器,在相同條件下進(jìn)行100 次分類實(shí)驗(yàn),取分類準(zhǔn)確率的平均值作為最終結(jié)果,實(shí)驗(yàn)結(jié)果如圖7 所示。3 種距離度量函數(shù)下,在驗(yàn)證集上平均準(zhǔn)確率的最大值見表2。
表2 不同距離度量函數(shù)下的分類結(jié)果Table 2 Classification results under different distance measurement functions(單位: %)
圖7 不同度量函數(shù)下的分類結(jié)果Fig.7 Classification results under different metric functions
圖7 中綠色虛線表示訓(xùn)練集上的分類結(jié)果,藍(lán)色點(diǎn)劃線表示測試集上的分類結(jié)果,左右兩側(cè)分別為使用SVM 和CNN 的分類結(jié)果。表2 表示采用SVM 與CNN 作為分類器時(shí),不同距離度量函數(shù)下的平均準(zhǔn)確率最大值以及相應(yīng)的字典原子數(shù)。
結(jié)合圖7 與表2,在同樣使用NMF 方法作為特征輸入的情況下,使用SVM 作為分類器的效果好于LeNet 的分類結(jié)果,其平均準(zhǔn)確率的最大值在β=1、d=1 時(shí)取得,即距離度量函數(shù)選取為廣義K-L散度、字典原子數(shù)選取為1的情形。在不同的度量函數(shù)下,訓(xùn)練集上的分類準(zhǔn)確率隨字典原子個(gè)數(shù)的增加而增加,但是在測試集上,分類準(zhǔn)確率隨字典原子個(gè)數(shù)增加呈下降趨勢。這是由于原始信號中實(shí)際上包含了目標(biāo)次聲信號與噪聲成分,而次聲信號的頻率低、頻帶窄,目標(biāo)信號與噪聲信號的頻帶混合嚴(yán)重[7]。當(dāng)字典原子數(shù)較低時(shí),其取值主要由目標(biāo)信號決定;而當(dāng)字典原子數(shù)較高時(shí),其取值同時(shí)受到噪聲影響。因此當(dāng)較高維度的特征向量輸入分類器后,產(chǎn)生了過擬合現(xiàn)象,識別準(zhǔn)確率下降。對于CNN 來說,其訓(xùn)練集上的準(zhǔn)確率較高,但在測試集上效果并不好,平均準(zhǔn)確率均在70%以下,這說明對于小樣本的次聲信號識別分類任務(wù)而言,采用較深度的復(fù)雜模型并不是一個(gè)理想的選擇,其效果差于SVM分類器。
2.3.3 實(shí)驗(yàn)結(jié)果分析
兩組實(shí)驗(yàn)在測試集上的分類結(jié)果對比如圖8所示,圖中增加了采用一維CNN 對HMS 特征進(jìn)行分類的結(jié)果。從圖8 中可以看出,基于HMS 的方法雖然提高了時(shí)頻分析的分辨能力,但其特征提取的方式較為直接,采用幅值相加的方式將全部分量值降為一維信號,對于原始信號來說有較多的信息損失,因此其識別的準(zhǔn)確率并不高。在測試集上,采用NMF-SVM 進(jìn)行分類的準(zhǔn)確率最高,并且在字典原子數(shù)較低時(shí)即取得了較好的結(jié)果。
圖8 兩組實(shí)驗(yàn)分類結(jié)果對比Fig.8 Comparison of classification results of two groups of experiments
結(jié)果表明,對于次聲分類這一小樣本分類任務(wù)而言,采用過高維度的特征向量或過深的網(wǎng)絡(luò)模型進(jìn)行分類效果欠佳,很容易受到訓(xùn)練模型過擬合的影響,因此應(yīng)該采用較低維的特征向量與淺層的分類模型,以提高泛化能力,獲得更好的分類性能。
對于小樣本的識別分類任務(wù),其關(guān)鍵問題在于特征向量的設(shè)計(jì)。本文以嚴(yán)格標(biāo)定的數(shù)據(jù)集為基礎(chǔ),針對爆炸、地震、閃電和再入這4類次聲信號的識別分類問題,引入了NMF 對次聲信號進(jìn)行特征提取,并分別使用SVM 與CNN 作為訓(xùn)練模型進(jìn)行分類。實(shí)驗(yàn)結(jié)果表明,NMF-SVM 在測試集上的平均準(zhǔn)確率可以達(dá)到83.13%,在當(dāng)前數(shù)據(jù)集上獲得了最佳的性能,是一種適用于次聲信號識別分類的方法。