基于高維隨機矩陣特征值之差的滾動軸承狀態異常檢測算法

朱文昌， 何雅娟， 王建波， 王 恒

(1. 南通大學 機械工程學院，江蘇 南通 226019；2. 南通大學 工程訓練中心，江蘇 南通 226019)

滾動軸承作為旋轉機械的關鍵部件，對其進行異常狀態檢測為設備預知維護提供了理論指導，避免嚴重故障的發生而造成的財產損失甚至人員傷亡。如何提高軸承異常狀態檢測結果的準確性成為研究的熱點和難點[1-2]。

Leite等[3]研究了基于熵的12個特征對軸承故障檢測性能的影響，使用不同的熵測度作為檢測熵特征變化的指標，提出一種非參數方法檢測軸承異常狀態，但在軸承早期異常狀態時，采集的信號成分復雜、故障信息微弱，易造成誤判。劉志亮等[4]從軸承的運行狀態角度考慮引入了安全域的概念，提出了一種基于核空間距離熵的安全域懲罰參數選擇算法，利用支持矢量數據描述(support vector data description, SVDD)對正常樣本構建安全域，結合核空間距離熵的安全域懲罰參數選擇算法尋找到最優懲罰參數提高了異常檢測準確率，但對高維數據處理能力不足。張西寧等[5]在格雷厄姆掃描法的基礎上，引入邊界軟化率增大了數據點外邊界的柔性，并結合射線法生成與輸入樣本反分布的數據集，使傳統模型成為擁有精細決策邊界的單類隨機森林，通過輸出待測數據的異常概率進行異常檢測，但選用多維縮放的降維方法進行特征約簡時可能會造成有用信息損失。Mao等[6]提出了一種使用半監督架構和深度特征表示的軸承早期故障在線檢測的新方法, 利用自動編碼器從目標軸承的正常狀態數據和輔助軸承的故障狀態數據中提取深度特征，引入一種安全的半監督支持向量機(safe semi-supervised support vector machine, S4VM)對軸承初期故障可有效識別，當數據樣本較大時仍可實現精確檢測，但對訓練集數據樣本要求較高，虛假特征易對S4VM造成干擾，影響最終檢測精度。

隨機矩陣理論主要研究隨機矩陣的樣本特征值、樣本特征向量及它們的統計函數在一般漸近體系下的漸近行為，目前已經在頻譜感知、電網配電、生物醫學等領域得到了成功應用[7-9]。Zeng等[10]基于隨機矩陣理論提出最大、最小特征值之比算法用于頻譜感知，克服了噪聲不確定性問題。劉威等[11]基于隨機矩陣單環定理對電網海量數據信息進行建模，實現擾動對電網影響程度分析和影響范圍評估。Koro?ak等[12]通過隨機矩陣理論的統計特性研究了胰島細胞的鈣離子分布，進而對疾病進行檢測。

為了適應設備健康監測在海量數據下得到較為精確的檢測結果，本文基于隨機矩陣理論，提出一種基于最大、最小特征值之差的異常檢測算法對軸承狀態進行檢測，應用結果表明該了方法的有效性和可行性。

1 基于高維隨機矩陣特征值之差的滾動軸承異常狀態檢測算法

1.1 高維隨機矩陣構造

根據滾動軸承健康監測的特點，采用大數據分析架構對其進行處理。設軸承健康監測數據維數為P、監測時間為T、監測節點數為Q，在采樣時刻ti，設備第j個節點(如測點)所監測的第n個運行狀態特征量定義為監測數據子空間

ljn(ti)(i=1,2,…,T;j=1,2,…,Q;
n=1,2,…,P)

(1)

對節點j而言，監測的所有Q個特征量可以構成一個列向量，即

lj(ti)=[lj1,lj2,lj3…,ljP]T

(2)

將不同采樣時刻的監測數據按照時間順序排列，構成一個時間序列矩陣，即

L=[lj(t1),lj(t2),…,lj(ti),…],
(j=1,2,3,…,Q)

(3)

該矩陣即對滾動軸承監測所得的數據源矩陣L∈cP×Q。

為了對滾動軸承特定時間節點數據進行實時分析，采用平移時間窗對軸承不同時刻的數據信息進行鎖定。設時間窗口H的規模為V×W，則在采樣時刻ti，對節點j構成矩陣為

Lj(ti)=[lj(ti-W+1),lj(ti-W+2),…,lj(ti)],
(i=1,2,…,P,j=1,2,…,Q)

(4)

在采樣時刻ti，窗口矩陣數據構造流程如圖1所示。

圖1 窗口矩陣數據構造流程圖Fig.1 Flow chart of window matrix data construction

(5)

式中，θj為隨機生成數，用于矩陣隨機化構造。

(6)

(7)

圖2 滾動軸承高維隨機矩陣構建流程圖Fig.2 Flowchart of constructing high-dimensional random matrix of rolling bearing

1.2 異常檢測特征值指標構造

由圖1、圖2所示方法構建軸承ti時刻的高維監測矩陣X，由式(8)計算其協方差矩陣

(8)

式中，σ2為矩陣X的方差。

對Bn進行特征值分解，如式(9)所示

Bn=EΛET

(9)

式中：Λ為矩陣Bn特征值對應的對角矩陣；E為正交矩陣。

滾動軸承數據信號采集中，存在未知的噪聲信號，提取出Bn的特征值極值如式(10)所示

(10)

式中：λv-min與λv-max分別為振動信號協方差矩陣最小、最大值；λn為未知噪聲協方差矩陣特征值。

為消除噪聲干擾，提高檢測指標的魯棒性，采用最大、最小特征值之差構造異常檢測指標D

D=λmax-λmin=λv-max-λv-min

(11)

結合閾值對滾動軸承異常狀態進行檢測，判決方法如式(12)所示

(12)

式中：Y1為檢測出軸承為正常狀態；Y2為檢測出軸承為異常狀態；γ為異常檢測閾值。

1.3 異常檢測閾值確定

(13)

(14)

式中，Bn為X的協方差矩陣。

且滿足

(15)

(16)

將滾動軸承的實際工作狀態劃分為正常狀態S1，異常狀態S2。考慮到滾動軸承早期故障特征較為微弱，異常信號不明顯，在進行異常檢測時可能出現誤判情況，即軸承實際已進入異常狀態S2，但判斷為正常狀態Y1或軸承實際為正常狀態S1，判斷為異常狀態Y2，造成誤判，誤警率ηw定義為

ηw=P(Y1|S2)+P(Y2|S1)

(17)

為提高檢測的精度，利用誤警率ηw研究檢測閾值γ。由于軸承異常檢測時難以給出S2狀態下檢測指標D的概率密度分布，無法通過P(Y1|S2)得到一個較為準確的判定閾值。因此，本文基于P(Y2|S1)情況下推導出閾值γ，如式(18)所示。

(18)

將式(13)、式(15)和式(16)代入式(18)中，得到檢測閾值γ數學表達式如式(19)所示。

(19)

式中：FT-W(·)為Tracy-Widom第一分布累計函數；σ2為矩陣X的方差；N為矩陣X的列數。綜上所述，基于最大、最小特征值之差的滾動軸承異常檢測算法流程，如圖3所示。其具體步驟如下：

圖3 基于最大、最小特征值之差的滾動軸承異常檢測算法流程圖Fig.3 Flowchart of an algorithm for detecting abnormal state of rolling bearings based on the difference between the maximum and minimum eigenvalues

步驟1將采集到軸承各個時刻的信號按監測時間序列構造出滾動軸承數據源矩陣L；

步驟2利用平移時間窗鎖定數據源矩陣各個時間段軸承數據，將鎖定的數據矩陣通過分段、隨機化、矩陣擴增及維度重構等方法構造出行列之比c(0

步驟3提取高維隨機矩陣最大、最小特征值構造出特征值之差指標D，利用隨機矩陣理論中的M-P及Tracy-Widom第一分布計算檢測閾值γ，結合檢測指標對軸承運行狀態進行判斷。

2 應用研究

2.1 數據來源

試驗數據來自辛辛那提大學智能維護系統(intelligent maintenance system，IMS)中心的滾動軸承全壽命試驗。試驗采用4個Rexnord ZA-2115雙列滾動軸承安裝在同一主軸上，軸承徑向施加2 722 kg恒定載荷，電機轉速2 000 r/min。試驗結束時軸承1外圈發生故障并失效，此過程共經歷164 h，每10 min對軸承進行一次數據采集，采樣頻率20 kHz，得到采集數據文件984組。本文對軸承1數據進行分析研究。

對采集到的軸承1數據通過式(1)～式(3)可得到監測數據源矩陣L∈c20 480×984，設平移時間窗口H的規模為20 000×2，對軸承不同時刻的數據信息進行鎖定；經過4次模擬矩陣的構建，并通過矩陣分段、隨機化、擴增及維度重構，最終得到ti時刻高維隨機特征矩陣X∈c400×500(行列之比c=0.8)，按照2.2節、2.3節所述方法求得軸承1不同時刻的檢測指標D與相應的檢測閾值γ，利用式(12)對軸承狀態異常檢測及判決。

2.2 誤警率對異常狀態檢測的影響

由式(19)可知，檢測閾值γ與誤警率ηw有關，不同誤警率會影響軸承異常檢測結果。取不同誤警率ηw對軸承1異常點進行檢測，結果如圖4～圖6所示。

圖4 當ηw=0.10時軸承1早期異常檢測結果Fig.4 Early abnormal detection results of bearing 1 when ηw=0.10

圖5 當ηw=0.05時軸承1早期異常檢測結果Fig.5 Early abnormal detection results of bearing 1 when ηw=0.05

圖6 當ηw=0.01時軸承1早期異常檢測結果Fig.6 Early abnormal detection results of bearing 1 when ηw=0.01

由圖4～圖6可知，當ηw=0.10時特征值指標與閾值在文件號532出現交點，可判斷早期異常的發生。文件號532前特征值指標位于閾值下方，進入異常狀態后特征值指標整體位于閾值上方。但在交點附近，出現誤判現象，即通過閾值判定軸承已進入異常狀態Y2，但仍存在少數特征值指標D小于閾值γ。當ηw=0.01時，此時閾值曲線可較好地將軸承1正常與異常狀態進行劃分，閾值與特征值指標的“混疊”現象基本消除，提高了滾動軸承早期異常檢測的準確率，但檢測的早期異常點時間滯后。

為進一步研究不同誤警率ηw取值時閾值γ與特征值指標D的關系，分別取ηw為0.01，0.05，0.10，對動態閾值曲線與特征值指標曲線之間的絕對誤差值之和進行分析，如圖7所示。隨著誤警率ηw減小，檢測閾值與特征值指標之間的絕對誤差值在逐漸減小，閾值曲線從整體上更逼近于檢測指標，提高了檢測的精度。因此，可根據采集到的軸承數據復雜程度及實際生產對異常檢測準確率進行綜合考慮，選擇合適的誤警率對檢測閾值進行調節，以達到最佳的檢測效果。

圖7 不同ηw下檢測閾值與檢測指標絕對誤差值之和Fig.7 The sum of the detection threshold and the absolute error value of the detection index under different ηw

2.3 滾動軸承異常狀態檢測結果分析

由2.2節可知，當ηw=0.01時，軸承1的檢測效果最佳，檢測結果見圖6。軸承1檢測指標與檢測閾值在文件號543時相交，在交點之后軸承1特征值指標上升趨勢明顯，閾值曲線位于特征值指標曲線下方，判定文件號543為軸承1早期異常點，與文獻[16]檢測出的早期異常點一致，與峭度指標(如圖8所示)檢測結果(文件號649)相比，本文所提算法可提前17.67 h檢測出軸承異常狀態。

圖8 基于峭度指標法檢測軸承異常狀態檢測結果Fig.8 Detection results of bearing abnormal state based on kurtosis index method

文獻[17]提出了采用隨機矩陣最大、最小特征值之比的檢測指標，并給出了檢測閾值t′，如式(20)所示

(20)

式中：M和N分別為矩陣的行和列；c為矩陣行、列之比；FT-W(·)為Tracy-Widom第一分布累計分布函數。

取誤警率ηw=0.01，c=0.8，基于最大、最小特征值之比的特征指標與閾值對軸承1進行異常狀態檢測，并與基于特征值之差的檢測結果進行比較，如圖10所示。采用特征值之比在文件號598時檢測出軸承異常狀態，與本文所提出算法相比滯后約9.2 h，說明最大、最小特征值之差算法對軸承早期異常更為敏感。

圖9 基于不同特征值指標的軸承異常狀態檢測結果比較Fig.9 Comparison of bearing abnormal state detection results based on different eigenvalue indexes

最大、最小特征值之比指標在軸承早期異常點后波動劇烈，單調性較差，這是由于軸承進入退化階段，提取的最小與最大特征值之間差異逐步增大，λmax/λmin將此差異放大了，造成相鄰時刻指標變化加劇，降低了指標的單調性。而采用特征值之差指標對軸承1運行狀態進行表征時，數據曲線波動小、穩定性高，可有效消除噪聲帶來的干擾，降低了軸承在退化階段數據差異對指標單調性的影響，進一步提高檢測指標魯棒性。

此外，由式(20)可知，最大、最小特征值算法給出的異常閾值主要由矩陣規模c和誤警率ηw確定，在兩者確定的條件下，異常狀態檢測閾值為常數。但是，滾動軸承健康監測是動態過程，不同的時間段采集的數據不同，在保證準確率的情況下，檢測閾值應隨著軸承監測信號變化而改變。本文利用平移時間窗對不同時刻的數據信息進行鎖定，結合當前時刻數據方差σ2與規模c給出了動態閾值的計算公式，當軸承處于正常狀態時，采集到的數據分布穩定，方差σ2變化較小，閾值接近于常數。當滾動軸承進入異常狀態時，隨著退化程度的加劇，不同采樣時刻數據間差異性越來越大，方差σ2變化較大，閾值隨之增加，其變化趨勢為隨特征指標變化的動態曲線，可實現對軸承異常狀態的有效檢測。

3 結 論

(1) 利用平移時間窗對滾動軸承不同時刻數據進行鎖定，并通過分段、隨機化、擴增和維度重構等方法進行高維隨機特征矩陣的構造。

(2) 基于隨機矩陣理論提出一種最大、最小特征值之差的滾動軸承異常狀態檢測算法，利用特征值之差構建出的檢測指標可降低未知噪聲帶來的干擾，提高指標的表征能力，結合誤警率所推導出動態檢測閾值，可實現對早期異常狀態的檢測。

(3) 利用最大、最小特征值之差算法對美國IMS軸承1進行異常檢測，結果表明，與最大、最小特征值之比算法比較，不僅可及早檢測出滾動軸承早期異常的發生，且檢測效果更優，更符合設備實際退化過程。