基于皮爾遜相關系數的滾動軸承混合域特征選擇方法

肖 楊 李 亞 王海瑞 常夢容

（昆明理工大學信息工程與自動化學院）

時域、頻域、時頻域各類特征作為具有顯著類別差異信息的非平穩統計特征，能有效提高滾動軸承狀態監測和故障診斷的性能和效率，因此研究人員在此方面開展了廣泛的研究［1］。 熊鵬博和王曉東提出了一種基于多時域特征與支持向量機（Support Vector Machines，SVM）的單向閥故障診斷方法［2］。 馬欣欣和郭敏將采集到的信號進行集合經驗模態分解 （Ensemble Empirical Mode Decomposition，EEMD）， 得到若干個固有模態函數分量（IMF）［3］，然后在前三階IMF的基礎上提取時域、頻域和希爾伯特域特征，融合提取的特征組成多域特征向量組，最后送入支持向量機分類器中進行分類。 彭濤等對原始信號分別生成時域、頻域狀態特征［4］，并利用多分辨率小波分解生成時頻域狀態特征， 構建出144個表征原始振動信號特征的混合域特征集。 李大江提出一種基于局部均值分解（LMD）和共空間模式（CSP）的時-頻-空多域特征提取方法［5］。

上述研究中用于分析的數據集包含數百個特征（或屬性），其中大多數可能與故障診斷系統任務無關或冗余，因此，特征提取、特征降維和特征篩選顯得尤其重要。 戴豪民等采用加權最大相關最小冗余的特征選擇方法［6］，選取7個有效特征向量，輸入至SVM得到不錯的效果。 Tang X H等利用特征對特征的最大信息系數（MIC）得到的弱相關特征子集和特征對類別的最大信息系數（MIC）得到的強相關特征子集［7］，通過交集運算合并為最終的診斷特征集，在一定程度上減少了特征數量。 白麗麗等利用拉普拉斯（LP）對能表征狀態的特征進行選擇［8］，將選擇得到的數據輸入到鯨魚算法優化的SVM進行模式識別，證明了特征提取的有效性。 但這些降維方法對軸承故障的特征集線性相關性過高且包含大量的冗余信息，所得到的低維空間對原始信號的解釋具有一定的片面性。

為了解決單域特征難以表達原始信號的振動規律、高維特征容易發生過擬合并且容易引發維數災難的問題，同時處理好冗余性、相關性問題， 筆者提出了基于皮爾遜相關系數（Pearson Correlation Coefficient，PCC）的滾動軸承混合域特征選擇方法，通過多元信息特征向量組確定一組高質量特征來進行穩定的預測。 首先從原始信號中提取6個時域無量綱向量、10個時域有量綱向量、4個頻域特征向量、6個小波變換特征向量和10 個自適應噪聲的完整集成經驗模態分解（CEEMDAN）特征向量，結合提取出的特征參數，構造軸承故障混合域特征集。其次，運用PCC進行特征選擇， 對提取的混合域特征進行相關性分析，根據相關性，從特征集中剔除不相關和冗余的特征，提取出易于識別的低維主特征向量。 最后將低維特征集導入到隨機森林中作為模式識別的輸入。

1 混合域特征集的構成

1.1 時域特征集

時域信息是以時間為變量描繪出信號的波形，作為衡量信號特征的重要指標［9］。時域信號包括量綱特征參數和無量綱特征參數［10］。 筆者主要引入6個時域無量綱參數、10個時域有量綱參數，組成16維時域特征向量構成時域特征集，包括最大值、最小值、峰值、峰峰值、平均值、絕對平均值、方根幅值、方差、標準差、有效值、峭度、偏度、波形因子、峰值因子、脈沖因子和峪度因子。

1.2 頻域特征集

頻域信息是以頻率為變量描繪出頻率信號的幅度，作為衡量信號特征的重要指標。 筆者通過提取4個常用的頻域特征向量來構成頻域特征集，包括平均頻率、重心頻率、均方根頻率和頻率標準差。 在構造頻域特征集之前，采用傅里葉分析對原始信號進行處理。

1.3 時頻域特征集

1.3.1 小波時頻域特征集

小波分解主要是以短時傅里葉變換的理論為基礎通過小波函數對時間序列進行細致描述［11］，因此可在不同維度進行信號分析［12］。 若滾動軸承在某一時刻發生突變，單獨依靠原始信號并不能對故障點進行準確描述，需要從時間序列的不同維度對信號進行分析，包括整體性分析和局部分析，因此通過小波分解的方法能夠對時間序列進行全面的刻畫，從而準確定位故障振動沖擊時刻［13］。

1.3.2 CEEMDAN分解時頻域特征集

CEEMDAN方法通過自適應加入白噪聲，克服了EEMD方法的模態混疊問題， 獲得了較好的模態分離譜，同時提高了運算效率［14］，因此筆者提取CEEMDAN分解后的各分量的時頻域特征，包括排列熵和瞬時能量，確保達到一個良好的特征分析效果。

2 特征選擇

特征選擇的目的不僅僅是為數據降維，還要消除冗余和無關的特性［15］。 通過度量特征間的相關性，可以消除冗余特征。 兩個特征之間的相關性越強，它們之間的冗余性和可替代性就越強［16］。 此外，通過測量特征與類別之間的相關性，可以消除不相關的特征。 特征選擇一般包括3個步驟：

a. 搜索。 在特征空間中搜索特征子集，每個子集被稱為一個狀態，由選定的特性組成。

b. 評價。 輸入一個狀態（子集），通過評價函數或預置的閾值輸出評價值，使評價值達到最優值。

c. 分類。 使用最終的特征集完成分類算法。

皮爾遜相關系數是由卡爾·皮爾遜提出的，定義為秩變量之間的相關系數［17］。 對于容量為n的樣本，將n個原始數據轉換為等級數據，相關系數為：

其中，rxy表示兩個變量x、y之間的線性相關程度，rxy的值在-1和+1之間；x=［x1，x2，…，xn］，y=［y1，x2，…，yn］；、分別為xi、yi的平均值。

若rxy＞0，表示兩個變量正相關，即一個變量的值越大，另一個變量的值也越大；若rxy＜0，則表示兩個變量負相關；當rxy=0時，表示x和y不相關。 相關系數絕對值越大，相關性越強；相關系數絕對值越接近0，相關性越弱［18］。 一般情況下，變量的相關強度由以下取值范圍來判斷：相關系數絕對值在0.8～1.0，非常強相關；在0.6～0.8，強相關；在0.4～0.6，中等相關；在0.2～0.4，弱相關；在0.0～0.2，非常弱相關或不相關。

3 實驗

3.1 軸承故障數據說明

實驗采用的是美國凱斯西儲大學軸承數據中心采集的軸承故障數據集。 數據集包含正常、內圈故障、滾動體故障和外圈故障（6點鐘方向）4種不同狀態的數據， 除正常數據外每種狀態有3種故障深度類型， 直徑分別為0.177 8、0.355 6、0.533 4 mm，即共10類故障類別。軸承電機載荷為0，軸承轉速為1 797 r/min。 每類數據劃分為115個分類樣本，10類總共1 150個樣本。 訓練集大小為700，即每類70個，測試集為450個，每類45個。 分類情況見表1， 其中IF、RF、OF分別為內圈、 滾動體、外圈故障（6點鐘方向）。

表1 軸承數據集描述

3.2 實驗特征提取

本實驗中， 每類原始時間序列的總長度為117 760，將其切分為115份，每份長度為1 024，分別提取每段時間序列的時域、頻域以及時頻域共36個特征。 原始特征向量記為A1，A2， …，An，n=1024，提取出新的特征向量表示為B1，B2，…，Bm，m=36（m＜n），fi為其對應映射函數，可表示為：

Bi=fi（A1，A2，…，An），i∈［1，m］

1～10為有量綱向量時域特征， 分別為最大值、最小值、峰值、峰峰值、平均值、絕對平均值、方根幅值、方差、標準差和有效值；11～16為無量綱向量時域特征，分別為峭度、偏度、波形因子、峰值因子、 脈沖因子和裕度因子；17～20為頻域特征，分別為平均頻率、重心頻率、頻率均方根和頻率標準差。

時頻特征主要提取小波變換和CEEMDAN相關特征。 其中， 小波變換將原始振動信號進行3層分解，劃分成8個子頻帶。 由于前4個子頻帶包含了原始信號的大部分能量， 因此提取前4個子頻帶小波尺度熵， 再提取信號的小波能量譜熵和小波奇異熵， 構成其中一個時頻特征子集。

分別對10類信號原始信號進行CEEMDAN分解，得到8個固有模態分量，如圖1所示。 第1個分量信號的振動頻率比其他分量的大， 而第2～6分量相比其他分量振動頻率更符合高次諧波的特征，所有分量能量大多集中在第2～6個分量中，并且在時間周期范圍內具有正弦波的特性，屬于有效分量。 而剩余分量含有的能量信息較少，振動沖擊特征不突出，與原始信號無太大關系，可視為無效分量。 因此，選取第2～6分量，并提取出各個模態的排列熵和瞬時能量，作為整個故障診斷數據的一個子集。

圖1 CEEMDAN分解圖

3.3 實驗分析

本次實驗將提取到的混合域特征進行相關性分析，分別計算每個特征向量與其他特征的皮爾遜相關系數值，并求其絕對值。 由于離群值對皮爾遜相關性分析較為敏感，若特征中存在離群點，則計算結果將小于實際計算結果，從而對特征分析產生錯誤的判斷， 圖2為隨機抽取部分特征散點分布圖。

本實驗用中位數對離群點進行替換處理，以確保皮爾遜相關系數的有效性。 圖2a顯示未處理前皮爾遜相關系數值為0.434， 由圖明顯可得出兩個特征向量存在相關關系， 但由于離群值的干擾，相關系數減小，使得計算結果與實際分布產生較大誤差，由中位數替代之后，重新計算皮爾遜相關系數為0.632， 計算結果符合特征規律。 表2為圖2a～d未處理和已處理的相關系數值變化情況。

表2 圖2a～d相關系數在離群點處理前后的變化

圖2 隨機抽取部分特征散點分布

設置皮爾遜相關系數閾值為0.450，再分別統計每類特征與其他特征的相關性的強弱， 大于0.450說明具有強相關性， 小于0.450說明特征之間弱相關，并根據弱相關性大小進行排序。 表3為計算的平均皮爾遜相關系數值，按由弱到強進行排序。

表3 計算的平均皮爾遜相關系數值

（續表3）

皮爾遜相關系數越大相關性越強，對于故障診斷越不利，因此需要篩選出相關性較弱的特征向量。 第1個特征向量與其他特征向量的相關性最弱，相關系數平均值僅為0.036，接近于零，對于分類具有較好的表現。而相關性最高的幾個特征向量的值已經超過閾值，表明其本身包含的特征信息與其他特征向量重復概率較大， 可以剔除。

為了更直觀地比較特征間的相關關系， 同時選取具有最小相關關系的特征向量（CEEMDAN第一分量的排列熵值） 和具有最大相關關系的特征向量（峰峰值）進行相關性實驗，如圖3、4所示，每個特征隨機抽取4個特征與之進行分析，繪制散點圖。為了排除隨機實驗的影響，本次實驗共進行了5次，統計每次實驗的平均相關系數值，結果見表4。

圖3 CEEMDAN第一分量排列熵與隨機特征相關關系散點圖

圖4 峰峰值與隨機特征相關關系散點圖

表4 5次實驗平均相關系數

由圖3、4可知，基于CEEMDAN第一分量的排列熵與其他特征均無明顯的相關性特點，點的分布規律較為均勻， 大多數表現為弱相關關系，是較為理想的分類特征。 而峰峰值與其他特征的分布特點是由左下角分布到右上角，呈現較強的正相關關系，點的分布擬合接近一條直線，不具有分類利用價值。

由表4數據可計算出峰峰值與其余全部特征的平均相關系數SAVG為0.462，已經大于0.450的閾值，且在進行隨機實驗時，相關系數浮動較大。 而CEEMDAN第一分量排列熵的SAVG僅為0.036 1，遠低于峰峰值的，并且在進行隨機抽取時，其實驗結果均在平均值上下較小范圍內浮動，具有較強魯棒性。 同時，將篩選后的特征向量集進行故障診斷，建立7個診斷模型進行準確率對比，模型分別為基于原始時域的故障診斷方法、基于原始頻域的故障診斷方法、基于原始時頻域的故障診斷方法、基于最大信息系數（MIC）的故障診斷方法、基于PCA降維的故障診斷方法以及基于原始混合域特征的故障診斷方法，分類方法采用隨機森林（rf）進行分類，MIC-rf模型和PCA-rf模型均在本實驗特征向量集下進行特征降維或篩選，具體數據見表5。

表5 7種分類模型性能比較

首先從準確率方面進行分析，基于原始頻域的故障診斷方法準確率約為81.00%，說明頻域特征表現力不強，對故障不能進行很好的識別。 而基于MIC的故障診斷方法和基于PCA的故障診斷方法由于特征選擇錯誤，導致重要性較高的特征被剔除，準確率僅約為85.00%。 識別準確率最高的為筆者所提出的方法和基于原始混合域的故障診斷方法，準確率可達約97.00%，因此可看出混合域特征集經過PCC篩選后， 重要性較高的特征向量得以保留，相關性較強的特征被剔除。

再從運行時間角度進行分析， 通過表5可看出基于MIC的故障診斷方法由于其近似算法計算時間較長， 導致其診斷時間遠遠超過其他模型。剩余模型的診斷時間相比MIC模型較好， 但與筆者提出的方法也有一定差距。 同時，從表中可看出雖然基于原始混合域的故障診斷模型準確率較高，但筆者提出的方法僅需約2 s便達到較高準確率，與原始混合域的故障診斷模型相比診斷時間縮短近2倍。 因此，從診斷準確率和運行時間進行綜合考慮，筆者提出的基于皮爾遜相關系數的研究方法具有更大的優勢。

4 結論

4.1 針對單個特征故障診斷精度不高、特征提取和特征集構建困難的問題，分別從原始信號的時域、 頻域和時頻域提取各個維度的綜合特征參數，充分利用了不同維度的有效信息。

4.2 為解決高維數據冗余性過高、相關性較強的缺點，提出利用皮爾遜相關系數對混合域特征集進行特征篩選，降低特征集相關性，為后續故障識別提供較為合理干凈的數據。

4.3 所提的基于皮爾遜相關系數的混合域軸承故障診斷方法，分類準確率可達97.32%，相比其他方法有較為明顯的優勢， 在進行特征選擇后，準確率未出現明顯下降，具有較高的工程應用價值。