基于改進VMD和多尺度排列熵的混合聲音事件特征提取

劉 爽， 沈希忠

（上海應用技術大學 電氣與電子工程學院，上海 201418）

聲音事件的識別和分析，在各個領域都體現出它的重要研究價值。通過識別某些動物的鳴叫，了解它們的習性，為進一步了解自然提高效率[1]；通過對風雨、雷電等聲音識別分析，提高對自然界環境氣候的認識，為天氣預報提供有效準確的信息[2-3]；由于不受光照和遮擋物的影響，聲音識別也廣泛應用于安防和軍事等領域[4-5]。特征提取是聲音事件識別系統中最重要的部分。

經驗模態分解（empirical mode decomposition，EMD）算法是Huang等[6]在1998年針對故障診斷中很難準確提取含有故障的特征這一現狀提出的，在很多方面都體現了其優越性，但在結果分析時發現模態會混疊。Wu等[7-8]提出了集合經驗模態分解（ensemble empirical mode decomposition，EEMD）方法，該方法可在一定程度上抑制模態混疊現象，但同樣仍然有一定程度的分解誤差。Dragomiretskiy等[9]于2014 年提出了變分模態分解（variational modal decomposition，VMD）方法，相比于EMD及EEMD等遞歸式的分解模式，VMD是少層數、高精度的非遞歸、變分形式的分解模態[10]。VMD方法大多應用于故障診斷領域[11-13]，隨著VMD與其他算法的結合，也逐漸應用于運行可靠性評估[14]、風功率預測[15]等研究。

利用VMD和EMD算法對帶有發動機信號、蟬鳴信號的混合聲音事件信號進行分解，得到本征模態函數（intrinsic mode function，IMF）分量，通過相關系數對各分量合成和拼接，并確定重構信號所屬類型，再做多尺度排列熵（multiscale permutation entropy，MPE）的分析，最終完成對混合聲音事件特征提取的目的，并對比VMDMPE和EMD-MPE的結果。

1 基本原理

1.1 EMD算法

EMD算法將非平穩信號按照其自身的尺度特征變為平穩信號，得到具有不同特征尺度的IMF分量。

對于任一信號X(t)，EMD算法分解流程[16]，具體可以描述為：

(1) 用3次樣條曲線分別連接X(t)上所有極大、極小值點，形成上、下包絡線。將信號X(t)上下包絡線均值m1的差記為

h1接 替X(t)成 為新的信號，重復上述操作，直到hi滿足IMF的2個條件，成功篩選出IMF1，記作C1。

(2) 將C1分 離出來，得到新信號r1，即

令r1重 復(1)，直到第n階的殘余信號rn成為單調遞增或者單調遞減函數，即

(3) 原始信號X(t)可以表示為

式中：rn(t)為 殘量，代表信號中的平均趨勢；Cj(t)代表信號從高到低不同頻率段的成分。能量在整個分解過程中沒有損失，可利用分解得到的各IMF分量和殘量來精準重構原始信號。

1.2 MPE

MPE是由Aziz等[17]在改進排列熵時提出，用于分析信號在不同粒度采樣（即多尺度）下的狀態，能夠在一定程度上實現抗干擾和抗噪。基本思想是多尺度粗粒化不同時間序列，然后計算它們的的排列熵，即MPE[18]。其計算流程如圖1所示。

圖1 MPE計算流程圖Fig. 1 Flow chart of MPE

（1） 粗粒化處理時間序列X={xi,i=1,2,···,N}， 得到粗粒化序列yj(s)

式中：s為尺度因子；N為序列長度；[N/s]表 示對N/s取整。

（2） 對進行時間重構得到

式中：m為嵌入維度；τ為延遲時間；l為 第l個重構分量，l=1,2,···,N?(m?1)τ。

（3） 將時間重構序列按升序排列，可得到符號序列S(r)=(l1,l1,···,lm)， 其中，r=1,2,···,R，且R?m!， 計算每一種符號序列出現的概率Pr。

（4）計算每個粗粒化序列的排列熵

當Pr=1/m!時 ，Hp(m)達 到 最 大值 ln(m!)。歸一化處理MPE值Hp(m)，得到多MPE值，即

式中，Hp(m)為歸一化處理后的排列熵值。

1.3 粒子群算法優化VMD參數

粒 子 群 算 法（particle swarm optimization，PSO）根據鳥群覓食特性建立模型，主要作用就是尋找得到最優的結果，已經應用于很多算法的優化處理中[19-21]。VMD算法中懲罰因子α和分解層數K影響著VMD的分解效果，但往往就只判斷K這一個參數，就會只得到一種最優值。因此，本文利用粒子群算法求解2個參數的最優值，其中選用的適應度函數的包絡熵平均值，PSO優化VMD參數過程如圖2。過程如下：

圖2 PSO優化VMD參數流程圖Fig. 2 Flow chart of PSO optimization VMD parameters

（1） 初始化PSO的參數，并選擇包絡熵平均值作為適應度函數。

（2） 以參數組合[α,K]作為PSO中的粒子，確定粒子的最初位置及移動速度。

（3） 在粒子不同的情況下做VMD分解，計算相應的適應度值，通過對比找到個體最優和群體最優值。

（4） 更新粒子的速度和位置。

（5） 計算最新粒子的適應度值，與原適應度值比較，將更好的結果替換成為新的個體最優和群體最優值，繼續下一時刻的運算。

（6） 當迭代次數達到最大設定時就得到最佳適應度值及粒子的位置，否則循環至步驟（3）迭代。

2 基本實驗分析

2.1 基于EMD算法分析

類似如上VMD的流程，得到經EMD分解cos信號的時域圖，如圖3所示。

圖 3 仿真信號及其對應EMD處理的時域圖Fig. 3 Simulation signal and its corresponding EMD processing time-domain diagram

將得到的分量與3種cos信號關聯，計算相關系數，如表1所示。由表1可知，與這3種cos信號的相關系數最高的是IMF5、IMF4、IMF1與VMD算法的相關系數相比，在經過處理得到分量與原信號的相關性方面，EMD算法略遜于VMD算法。

表1 EMD得到的各分量與3個cos信號的相關系數Tab. 1 The correlation coefficients of each component obtained by EMD and the three cosine signals

2.2 利用VMD-MPE對實際信號分析

2.2.1 模態參數K的確定

實際錄制一段包含有混合聲音事件的文件，截取長度為1.04 s的信號，作為特征提取的樣本信號，利用VMD算法對其進行分解，聲音事件信號包括發動機聲音、蟬鳴聲音信號以及周圍環境聲。由于VMD算法中模態分量個數的多少都會影響到后面的結果，太多，會產生模態混疊，導致特征模態函數沒有任何物理意義；太少，會丟掉原本有用的信息。因此，在對樣本信號進行VMD分解前，要先確定模態分量個數K。本文采用計算各模態分量的中心頻率來確定K值。不同K值下各模態分量的中心頻率如表2所示。可以看出：當K=4、5、6時，中心頻率存在突然減小的情況，屬于過度分解狀態，所以選擇K=3。

表2 不同K值對應的中心頻率Tab. 2 The center frequency corresponding to different K values

2.2.2 利用VMD算法分解實際錄制語音

利用2.2.1確定模態參數K=3，用VMD算法對實際錄制的這段車載語音文件進行分解，得到圖4。將得到的各分量與發動機信號、蟬鳴信號互相關，得到各相關系數，如表3所示。

表3 VMD處理的各分量與2種信號相關系數Tab. 3 The correlation coefficients between the VMD components and the two signals

圖4 VMD分解圖Fig. 4 The VMD decomposition diagram

由表3可見，與發動機信號相關系數最高的是IMF1，與蟬鳴信號相關系數最高的是IMF2，因此，表示經VMD算法分解得到的3個分量依次為發動機信號、蟬鳴信號以及周圍環境產生的信號。

2.2.3 MPE計算

利用MPE對上述得到的分量進行處理，經過VMD分解和MPE的特征提取后，此特征向量可以在一定程度上區分信號。本文取嵌入維數m=4，提取的尺度因子s=20，結果如圖5所示。

圖5 嵌入維度m =4， 尺度因子s =20時MPE分布Fig. 5 MPE distribution when embedding dimension m=4 and scale factor s=20

圖5中，“—”代表IMF1（即發動機信號），“”代表IMF2（即蟬鳴信號），“”代表IMF3（即周圍環境產生的信號）。可以明顯看出，嵌入維數m=4時，發動機信號被完全分離出來，蟬鳴信號和周圍環境產生的信號的MPE值有部分交疊現象，此時若根據MPE值進行分類識別不能較好地區分各個變量。

2.3 利用EMD-MPE對實際信號分析

2.3.1 利用EMD算法分解實際錄制語音

利用EMD算法對2.2.1中錄制的1.04 s的信號進行處理，結果如圖6所示，IMF分量個數是16，將得到的各分量與發動機信號、蟬鳴信號互相關，得到各相關系數，如表4所示。

圖6 EMD分解圖Fig. 6 The EMD decomposition diagram

由表4可見，與發動機信號相關系數較高的是IMF6～IMF11，與蟬鳴信號相關系數最高的是IMF1～IMF3，將IMF1～IMF3合成一個信號，記為信號1，IMF6～IMF11合成一個信號，記為信號2，剩余的分量合成一個信號，記為信號3，再將得到的3個信號與發動機信號、蟬鳴信號互相關，得到相關系數，如表5所示。

表4 EMD處理的各分量與2種信號相關系數Tab. 4 The correlation coefficients between the EMD components and the two signals

表5 重構后分量與2種信號的相關系數Tab. 5 The correlation coefficient between the reconstructed components and the two signals

由表5可見，與發動機信號相關系數較高的是信號2，與蟬鳴信號相關系數最高的是信號1，因此，表示經EMD算法分解并合理合成得到的3個信號依次為蟬鳴信號、發動機信號以及周圍環境產生的信號，將合成的3個信號表示如圖7所示。

圖7 EMD重構的3個信號Fig. 7 3 components reconstructed by EMD

2.3.2 MPE的計算

利用MPE對EMD算法分解并合理合成得到的信號進行處理，取嵌入維數m=4，提取的尺度因子s=20，結果如圖8所示。圖8中，“—”代表IMF1（即蟬鳴信號），“”代表IMF2（即發動機信號），“”代表IMF3（即周圍環境產生的信號）。可以明顯看出，嵌入維數m=4時，發動機信號被完全分離出來，而在尺度因子s=6之后蟬鳴信號和周圍環境產生的信號的MPE值就已經有部分交疊現象。因此，與VMD-MPE相比，若根據EMD-MPE進行分類識別會更有難度。

圖8 嵌入維度m =4， 尺度因子s =20時MPE分布Fig. 8 MPE distribution when embedding dimension m=4 and scale factor s=20

3 復雜實驗分析

該節的實驗收集了稍復雜的混合聲音信號，包括了汽車發動機聲音信號、警報聲音信號以及蟋蟀聲音信號等混合信號，首先利用EMD、VMD以及PSO算法算法改進VMD算法中的兩個參數后，對聲音信號進行分解，通過相關系數分析確定每個分量的類型，然后結合MPE對分解得到的分量進行特征提取，得到MPE分布圖，最后通過對比得出實驗結果。

3.1 分解實驗

3.1.1 EMD分解

利用EMD算法分解該節中的混合聲音，結果如圖9所示，其相關系數，如表6所示。

圖9 EMD分解圖Fig. 9 The EMD decomposition diagram

表6 EMD處理的各分量與3種信號相關系數Tab. 6 The correlation coefficients between the EMD components and the three signals

根據相關系數重構信號，可以看出，與發動機信號相關最高的是IMF4～IMF10，重組成為信號1，與警報聲音信號相關最高的是IMF1～IMF3，重組成為信號2，而蟋蟀聲音信號沒有被分解出來，剩余的分量重組成為信號3，可以認為是周圍環境聲音信號，重組后的結果如圖10所示，相關系數如表7所示。

表7 重構后的信號與3種信號的相關系數Tab. 7 The correlation coefficient between the reconstructed signal and the three signals

圖10 EMD重構的3個信號Fig. 10 3 signals reconstructed by EMD

可以看出，由EMD分解并重構后的3個信號依次為汽車發動機聲音信號、警報聲音信號以及周圍環境信號，而在混合聲音信號中的蟋蟀聲音信號并沒有被EMD算法分解出來。

3.1.2 VMD分解

首先利用VMD算法預分解，然后根據中心頻率的情況確定分解層數K值，對于該節中的混合聲音信號，中心頻率如表8所示，由于K=4時，IMF2與IMF3出現模態混疊現象，因此選用K=3。

表8 不同K值對應的中心頻率Tab. 8 The center frequency corresponding to different K values

選用K=3，此時利用VMD算法分解，結果如圖11所示。

圖11 VMD分解圖Fig. 11 The VMD decomposition diagram

將得到的各分量與汽車發動機信號、警報聲音信號、蟋蟀聲音信號互相關，得到各相關系數，如表9所示。由表9中相關系數的高低可以判斷信號屬于哪種類型，經VMD算法分解得到的3個分量，依次代表汽車發動機信號、警報聲音信號、蟋蟀聲音信號。

表9 VMD處理各分量與3種信號相關系數Tab. 9 The correlation coefficients between the VMD components and the three signals

3.1.3 PSO-VMD分解

利用PSO算法優化VMD的懲罰因子和分解層數2個參數，以VMD分量包絡熵平均值作為適應度函數，如圖12所示，當滿足停止條件時得到迭代結果，懲罰因子α=962， 分解層數K=3。將其帶入VMD算法，結果如圖13所示，各相關系數如表10所示。

表10 PSO-VMD處理的各分量與3種信號相關系數Tab. 10 Correlation coefficients between each componentprocessed by PSO-VMD and the three signals

圖12 PSO算法優化VMD參數圖Fig. 12 PSO algorithm optimization VMD parameter diagram

圖13 PSO-VMD分解圖Fig. 13 PSO-VMD decomposition diagram

因此，可以判斷，經PSO優化VMD算法參數，分解得到的3個分量依次代表汽車發動機信號、警報聲音信號、蟋蟀聲音信號。

3.2 MPE特征提取

3.2.1 參數選擇

首先在VMD算法分解信號的情況下，假設尺度因子取s=20的 情況下，比較嵌入維度m的情況，如圖14所示，可以很明顯的看出在m=4時，效果是最好的，3個分量很容易就能夠區分出來，因此確定嵌入維度m=4。

嵌入維度確定m=4后 ，確定尺度因子s，如圖15所示，對比后發現，尺度因子s=20是最合適的值。

圖15 MPE在不同尺度因子下的分布圖Fig. 15 Distribution of MPE at different scale factors

3.2.2 實驗對比

將EMD、VMD、PSO-VMD分解后得到的分量經MPE提取特征，對比得到的MPE分布圖，如圖16所示。

圖16（a）中的3個分量很好區分，但是EMD分解時并沒有將蟋蟀聲音信號分解出來，從某種程度上來說，可以認定該信息被丟失，而圖16（b）、16（c）中的3個分量分別代表了汽車發動機聲音信號、警報聲音信號、蟋蟀聲音信號，并且對于區分各分量一目了然。因此，VMD及PSO-VMD算法優于EMD算法。

圖16 不同算法下的MPE分布圖Fig. 16 MPE distribution diagram under different algorithms

4 結 語

本文利用PSO算法改進VMD算法，并且將PSO-VMD算法、VMD算法、EMD算法以及MPE應用于混合聲音事件的特征提取。

在基本實驗分析中，利用VMD算法和EMD算法對2種簡單的混合信號進行分解，一種信號是存在3個cos信號的混合信號，另一種信號是錄制的一段帶有發動機信號、蟬鳴信號等信號的混合聲音事件信號，結合MPE計算EMD-MPE和VMDMPE的值，并得到分布圖，通過對比可以得到VMD算法的結果優于EMD算法。

此外，還利用PSO優化VMD算法參數，利用3個方法對復雜的混合聲音信號做分解，混合聲音信號包含有汽車發動機聲音信號、警報聲音信號、蟋蟀聲音信號等信號。利用EMD算法分解該信號并重組信號，發現只能得到汽車發動機聲音信號和警報聲音信號，而蟋蟀聲音信號沒有被分解出來，VMD算法和PSO-VMD算法分解穩定，3種聲音信號都被成功分解出來。在結合MPE后，可以通過圖示輕松的看出VMD算法和PSO-VMD算法優于EMD算法。