基于壓縮感知和小波信息熵的滾動軸承特征提取方法

李中， 王星， 盧春華

(1. 華北電力大學電子與通信工程系， 保定 071003; 2. 華北電力大學河北省電力物聯網技術重點實驗室， 保定 071003)

滾動軸承在工業領域應用廣泛，在線監測診斷分析其運行狀態十分重要。基于 Nyquist-Shannon 采樣定理的傳統采樣方法,采樣的頻率要求遠遠大于故障狀態下信號的最高頻率兩倍,而故障狀態下的滾動軸承振動信號的高頻分量會增加，這給在線監測和故障診斷過程中的信號采集、傳輸、存儲和處理提出了更高的要求。Donoho[1]提出了壓縮感知技術，采集少量信號即可包含傳統采樣信號的大部分信息，為解決這一問題提供了新的思路。

中外學者就振動法分析診斷軸承運行狀態做了大量研究工作。Min等[2]根據小波包分解后的能量分布構造特征向量，確定滾動軸承故障的類型和嚴重程度。陳宗祥等[3]分析比較了不同狀態下經小波包分解后各子頻帶信號能量分布并結合小波包重構和Hilbert解調技術判斷滾動軸承的故障位置。陳湘中等[4]以小波包分解的各節點能量作為特征參量，提出基于蟻群優化的K均值聚類算法識別故障類型。蘇文勝等[5]利用小波包分解求取最大能量自帶的樣本熵，將此樣本熵作為故障診斷的特征參量。張龍等[6]利用加權能量算子對小波包分解最優系數的節點解調分析信號的瞬時能量，最后分析包絡譜頻率成分判斷故障類型。陳慧等[7]應用局部均值算法對最優小波系數加權運算增強滾動軸承的沖擊部分。Wang等[8]采用滑動均值濾波算法對原始數據作預處理，將平滑后的數據作為故障診斷的輸入向量。馮桓榰等[9]應用小波包分解提取振動信號小波包熵作為軸承特征參量。Verstockt[10]基于卷積神經網絡的狀態監測特征學習模型從數據本身中自主提取軸承故障特征，然后在Hilbert解調譜中提取特征頻率判斷軸承的運行狀態。這些成果基于各種信號處理和機器學習等技術，較好解決了滾動軸承振動信號的特征提取和故障診斷問題，但沒有對振動信號進行壓縮處理，其監測和診斷過程存在成本高、耗時長等不足。

針對振動信號監測分析中應用壓縮感知技術的研究主要包括：一是振動信號的壓縮與重構；二是基于重構信號提取特征進行故障診斷。劉暢等[11]以幅值百分比作為閾值對稀疏系數進行稀疏處理，采用高斯隨機測量矩陣實現滾動軸承振動數據的壓縮。唐其琯等[12]根據信號頻域稀疏特性設計一種確定性隨機觀測矩陣，采集的滾動軸承振動信號能以較高精度恢復重構。何天遠等[13]在信號的稀疏分解中基于K-奇異值分解生成字典，增加了重構信號和原始信號的相似度。王海明等[14]用快速譜峭度確定最優濾波器參數，在傅里葉稀疏基下對去除噪聲的信號包絡信號進行重構，分析重構信號頻譜判斷軸承狀態。王江萍等[15]通過高斯隨機測量矩陣采集齒輪振動信號，經正交匹配追蹤(orthogonal matching pursuit，OMP)算法重構信號，從重構信號的Hilbert解調譜中提取特征頻率判斷齒輪的運行狀態。

綜上所述，現有應用壓縮感知理論采集滾動軸承振動信號實現在線監測和故障診斷有關研究中，存在重構計算復雜且難免產生新的誤差，影響故障診斷精度的問題。考慮大多狀態下正常信號高精度重構對監測和故障診斷意義不大，現提出直接提取壓縮后的滾動軸承振動信號特征進行故障診斷，提取壓縮信號的歸一化小波信息熵作為故障診斷的依據。該方法首先利用部分哈達瑪(PartHadamard)測量矩陣完成對振動信號的壓縮，應用小波包分解壓縮信號得到各子頻帶的小波信息熵，將各子頻帶歸一化小波信息熵作為故障診斷的特征參量。在標準數據集上，應用粒子群優化的支持向量機(particle swarm optimization-support vector machine, PSO-SVM)、K近鄰(K-nearest neighbor, KNN)、決策樹(decision tree, DT)、概率神經網絡(probabilistic neural network, PNN)4種分類算法完成了滾動軸承故障診斷實驗。

1 壓縮感知和小波包分解

1.1 壓縮感知

壓縮感知的數學模型為

y=ΦMlength×NΨθ

(1)

式(1)中：y為采集的壓縮信號；Φ為Mlength×N的測量矩陣；Ψ為原始信號的變換基；θ為信號在變換基下稀疏分解的系數；N為原始信號長度；Mlength為壓縮信號的長度，原始信號本身稀疏或在變換基Ψ下可以稀疏表示。

在滿足壓縮感知的數學模型的所有條件中，求取最稀疏解為求解最優化問題。

(2)

式(2)中：x為原始信號；‖x‖0為原始信號的L0范數。

求解最優解問題是一個非確定性多項式難題(non-deterministic polynomial hard, N-hard)問題，引入約束等距(restricted isometry property, RIP)性質，可表示為

0<δk<1

(3)

式(3)中：c為稀疏信號；δk為原始信號稀疏分解后的稀疏度。

在測量矩陣Φ滿足RIP性質的情況下，可以將求解L0范數問題轉化為求取L1范數問題，通過求取L1范數可以恢復出稀疏信號。

(4)

式(4)中：‖x‖1為原始信號的l1范數。

1.2 小波包分解

小波包分解可以在低頻和高頻完成對信號的分解，相比小波分解小波包分解有更高的頻率分辨率。小波包分解算法的遞推公式可表示為

(5)

小波包重構算法計算公式為

(6)

式(6)中：h1(l-2k)、g1(l-2k)分別為低通濾波器系數、高通濾波器系數。

通常情況下，小波包分解提取模式特征需要考慮兩個問題：一是小波基函數的選擇；二是分解層數的確定，層數過少難以提取信號的故障細節信息，層數過多則計算復雜、提取特征維度高，傳輸和存儲開銷大。

2 基于壓縮感知和小波信息熵的滾動軸承振動特征提取

所提出的滾動軸承振動信號特征提取方法，包括信號壓縮和特征選擇兩個基本環節：首先構造并應用PartHadamard測量矩陣完成振動信號的壓縮，然后對壓縮信號進行db5小波包分解，選擇歸一化小波信息熵作為滾動軸承振動信號特征。

2.1 構造測量矩陣

測量矩陣的選擇直接決定了采集的振動信號，采集的振動信號是否包含原始信號的大部分信息是衡量測量矩陣的重要準則，硬件設計能否易于實現是評判觀測測量矩陣的重要標準。根據壓縮感知理論，測量矩陣決定了獲取的壓縮信號，對基于壓縮信號的提取特征具有重要影響。構造滿足RIP性質的測量矩陣，才能保證獲取的壓縮信號包含原始信號大部分信息，滿足特征提取和后續故障診斷的需求。其次，還需考慮軟硬件成本和實現難易程度。

考慮PartHadamard矩陣滿足式(3)中所示的約束等距性質[16]，并且PartHadamard矩陣是由±1組成的正交矩陣，信號采集只涉及加減運算，易于硬件實現[17]，選取PartHadamard矩陣作為測量矩陣。構造PartHadamard矩陣需要借助哈達瑪(Hadamard)矩陣，選取Hadamard矩陣的前M′行向量構成測量矩陣。Hadamard矩陣構造過程如下。

最低階的Hadamard矩陣為二階Hadamard矩陣，其表達式為

(7)

四階Hadamard矩陣為

(8)

N階Hadamard矩陣為

(9)

2.2 db5小波包分解

由于采集的軸承故障信號包含很多沖擊成分，對滾動軸承振動壓縮信號小波包分解采用db5小波基函數[18]，通過滾動軸承故障信號的特征頻率分析確定分解層數[19]。具體方法是：首先計算滾動軸承各種故障狀態的故障頻率，然后求取滾動軸承振動信號的最佳分解層數J，其表達式為

(10)

式(10)中：fs為采樣頻率；f為故障頻率。

2.3 特征選擇

壓縮后的滾動軸承振動壓縮信號經J層小波包分解后，在頻域內被平均分解成2J個子頻帶，且每個子頻帶互不重疊，分解后得到2J個小波包系數，基于最后一層小波包分解每個子頻帶的小波包系數得到每個子頻帶包含的信號能量，根據每個子頻帶信號的能量比值得到每個子頻帶的小波信息熵，小波信息熵能夠反映滾動軸承運行狀態的不確定程度。將計算獲取的小波信息熵歸一化作為提取滾動軸承的特征參量，完成特征提取。軸承發生不同故障狀態產生的振動信號不同，在頻域內經小波包分解后的小波包系數不同，進而得到的歸一化小波信息熵也不同，得到的歸一化小波信息熵可以作為滾動軸承振動特征。歸一化小波信息熵的具體計算過程如下。

設第J層每個子頻帶的小波包系數分別為CJ(0)，CJ(1)，…，CJ(2J)。每個子頻帶的能量為

Ei=|DJ(k)|2

(11)

式(11)中：i為子頻帶編號，i=1,2,…,2J；k為小波系數的編號，k=1,2,…,2J，i和k一一對應；DJ為信號被分解后第k個子頻帶對應的小波系數數值。

第J層子頻帶的總能量為

(12)

第J層子頻帶的能量分布概率為

(13)

第J層每個子頻帶的小波信息熵為

Si=-PilgPi

(14)

歸一化小波信息熵為

(15)

式(15)中：minSi為第J層中子頻帶的最小小波信息熵；maxSi為第J層中子頻帶的最大小波信息熵。

3 實驗

選擇美國凱斯西儲大學(Case Western Reserve University, CWUR)電氣工程實驗室滾動軸承數據集進行數據壓縮、特征提取和故障診斷實驗，分析檢驗本文方法的可行性和有效性。

3.1 實驗數據

CWUR數據中的滾動軸承型號為SKF6205，軸承分別在內圈故障、外圈故障、滾動體故障和正常4種狀態下工作，以12 kHz的采樣頻率測量采集振動信號。選取滾動軸承空載條件下，故障尺寸直徑為0.177 8 mm的滾動軸承振動數據作為本次仿真的實驗數據集，4種狀態下均包含120個樣本，每個樣本有768個采樣點，實驗數據集中滾動軸承轉速等其他參數如表1所示。4種狀態下的振動數據波形示例如圖1所示。

表1 滾動軸承參數Table 1 Rolling bearing parameters

3.2 實驗過程與結果

實驗仿真包括構造不同PartHadamard測量矩陣壓縮振動信號、不同分解層數的特征提取和不同分類器的故障診斷，從壓縮性能和診斷精度兩個方面分析本文方法的整體性能。在實驗數據集上，在相同壓縮率下，分別采用不同構造方式的測量矩陣完成滾動軸承振動信號的壓縮，用以檢驗不同測量矩陣構造方法的壓縮性能。

然后完成壓縮信號的特征提取，采取不同的小波包分解層數來檢驗特征選擇的優劣。最后，基于每種狀態下特征提取后實驗數據的50%作訓練樣本，另外50%作測試樣本，分別采用PSO-SVM、KNN、DT和PNN4種分類算法進行故障診斷實驗，分析檢驗特征提取后在故障診斷中的精度。

首先進行不同測量矩陣構造方法的實驗仿真。測量矩陣的大小取決于信號的壓縮率，設壓縮率a=(Mlength-N)/N，其中M為壓縮信號每個樣本的數據長度，N為原始信號每個樣本的數據長度。PartHadamard測量矩陣由Hadamard矩陣的M′行向量構成，選取M′行向量的方式有多種：選擇前M′行向量、后M′行向量、選偶數行向量、選奇數行向量等。偶數行選取方式：如果選擇的偶數行大于要求的采樣點數，就繼續從偶數行中挑選偶數行；如果選擇的偶數行小于要求的采樣點數，就在剩余的奇數行中選擇偶數行補充[20]。選奇數行和選偶數行的方式類似。分別應用4種方法構成的PartHadamard矩陣采集滾動軸承振動信號。由于PartHadamard矩陣需要滿足：行數是2的整數次冪、被12整除后是2的整數次冪、被20整除后是2的整數次冪中的某一個條件，所以測量點數Mpoint不是連續的。

圖1 不同狀態下的滾動軸承信號Fig.1 Rolling bearing signals in different states

在選擇測量矩陣中，用小波包三層分解提取歸一化小波信息熵作為故障診斷的輸入向量。分別采用SVM、KNN(K=1)、DT和PNN分類方法進行故障診斷仿真，實驗結果依次記錄如表2所示。

表2 基于不同測量矩陣的診斷結果Table 2 Diagnosis results based on different measurement matrices

實驗結果顯示，在選擇Hadamard矩陣前M′行向量構造測量矩陣時，得到的壓縮信號在故障診斷中準確率顯著優于其他三種構造方法。在測量點數為128時，其他三種構造方法采集的信號，在故障診斷中準確率明顯降低。所以PartHadamard測量矩陣選擇Hadamard矩陣前M′行向量。

然后對壓縮信號特征提取，考慮到信號經小波包不同層數分解得到的特征參量和特征維度不同，導致故障診斷結果和計算時間不同。需要根據滾動軸承振動信號確定小波包最佳分解層數，分析確定壓縮信號小波包最佳分解層數需要考慮滾動軸承各狀態的故障頻率，實驗數據集中3種故障狀態的故障頻率計算方法如下。

內圈故障頻率計算公式為

(16)

外圈故障頻率計算公式為

(17)

滾動體故障頻率計算公式為

(18)

式中：r為被測滾動軸承的轉速；n為滾動軸承滾珠個數；d為軸承滾珠直徑；D為軸承節徑；α為軸承滾珠接觸角。

根據表1，由式(10)和式(16)～式(18)計算得出內圈故障、外圈故障和滾動體故障的小波包最佳分解層數依次分別是小于4、小于4和小于5，綜合分析最佳分解層數取3。

為進一步對比分析壓縮后振動信號基于不同小波包分解層數的特征提取性能。分別進行壓縮數據的2層小波包分解、3層小波包分解和4層小波包分解，提取滾動軸承壓縮信號的特征，繼而應用PSO-SVM、KNN(K=1)、DT和PNN 4種分類算法進行故障診斷。不同分解層數特征提取的故障診斷結果記錄如圖2所示。

由圖2可知，基于3層小波包分解的特征提取的故障診斷準確率普遍高于2層和4層小波包分解，在相同壓縮率下的故障診斷準確率最高。再考慮分解層數與特征參量的維度關系，4層小波包分解得到的是16維特征參量，3層小波包分解得到的是8維特征參量，3層小波包分解在特征提取和故障診斷中相比4層小波包分解降低了特征參量的存儲空間和節約了計算時間。本文方法兼顧了故障診斷精度和計算復雜度，特征提取性能優異。

為進一步檢驗本文方法的性能，對比文獻[21]，

圖2 不同分類方法的故障診斷結果Fig.2 Fault diagnosis results of feature extraction with different decomposition levels

在相同數據集上應用本文方法進行特征提取、相同比例訓練集和測試集、采用相同的PSO-SVM算法進行故障診斷實驗。文獻[21]通過計算原始信號和壓縮信號的13種敏感特征參量，選取壓縮前后特征參量變化不大的方差、峭度因子、波形因子作為故障診斷的輸入向量。診斷結果對比如表3所示。

可以看出，在相同壓縮率條件下，用本文方法相比文獻[21]，可以獲得更高的故障診斷精度。為分析本文方法對不同分類算法的普適性，基于上述提取特征和相同的訓練集與測試集劃分方法，分別應用KNN(K=1)、DT和PNN分類算法進行故障診斷仿真，結果記錄如表4所示，可以看出，本文方法提取特征在不同分類算法中都獲得了較高的故障診斷精度。

表3 不同特征提取方法的故障診斷精度Table 3 Fault diagnosis accuracy of different feature extraction methods

表4 不同分類算法故障診斷精度Table 4 Fault diagnosis accuracy of different algorithms

4 結論

針對滾動軸承在線監測和故障診斷應用，提出一種基于壓縮感知和小波信息熵的滾動軸承振動特征提取方法。基于不同測量矩陣構造方法、不同壓縮率和不同分類算法，在標準實驗數據上完成了特征提取和故障診斷實驗。結果表明所提方法能夠在較高的數據壓縮率條件下，保持較高的故障診斷精度。

(1)提出一種基于壓縮感知和小波信息熵的滾動軸承振動特征提取方法。其主要思想是針對在線監測和故障診斷實際需求，對壓縮后的監測信號直接求取其故障診斷特征。

(2)滾動軸承故障診斷仿真分析表明，在測量矩陣的行向量選取上，選擇前M′行向量構造順序Parthadmard矩陣作為測量矩陣，相比其他構造方法具有更高的診斷精度。

(3)仿真結果表明，傳統采樣方法中振動信號分析的小波包分解最佳層數方法，在壓縮后的信號處理中仍然適用。