基于EEMD和IHHO-LSSVM滾動軸承剩余壽命預測方法的研究

胡豁然，李亞莎

(三峽大學電氣與新能源學院，湖北 宜昌 443002)

滾動軸承是電力設備中應用最廣的一種典型重要零部件，滾動軸承的運行狀況直接影響電氣設備的安全、可靠運行，壽命是衡量軸承的重要指標[1-2]。隨著機器學習的發展，基于數據驅動的軸承壽命預測逐漸成為國內外重點研究方向。但軸承原始信號中不可避免會有大量噪聲，影響了軸承壽命預測的準確性，因此需要對原始信號進行處理，文獻[3]針對原始信號中的大量噪聲，采用EEMD對原始數據進行分解、重構。為了進行退化特征的提取和壽命預測，文獻[4]提出了一種基于多頻率尺度模糊熵和ELM的滾動軸承剩余壽命預測；文獻[5-6]提出了一種基于改進HMM和Pearson相似度分析的滾動軸承自適應壽命預測方法；文獻[7-8]提出了一種基于小波包能量譜和改進FOA-GRNN的軸承壽命預測；文獻[9-10]提出了基于t-SNE和LSTM對滾動軸承剩余壽命的研究。上述方法為軸承退化特征的提取和壽命預測提供了有效借鑒，但仍有考慮不周的地方，例如對原始信號進行降噪的方法仍存在不足，壽命預測只能預測短期的退化趨勢且預測結果較為保守。

針對以上問題，提出一套新的滾動軸承壽命預測方法，先用EEMD對原始信號進行分解，然后根據相關性系數指標對IMF進行篩選并重構信號，并對特征指標集進行核主成分分析，選取累計貢獻率符合要求的核主成分進行分析，建立IHHO-LSSVM模型，采用改進的哈里斯鷹算法優化LSSVM參數，提高滾動軸承剩余壽命預測模型的性能和精度。

1 軸承性能退化評估模型構建

1.1 EEMD原理

WU[11]等提出的EEMD是EMD的一種改進算法，EMD在使用過程中存在因噪聲或脈沖產生的模態混疊現象，為了解決該問題，EEMD利用白噪聲的頻譜均勻性，通過每次給信號加入不同幅值的高斯白噪聲改變信號極值點特性，形成基本無擬合誤差的上、下包絡線，從而消除模態混疊現象。

EEMD具體分解過程如下。

a.在信號中加入一定幅值的白噪聲序列得到新的信號x(t)。

b.對新信號x(t)進行EEMD分解，獲得k個IMF分量cj(t)和一個余量r(t)，(j=1,2，…，k)。

(1)

c.每次加入相同振幅的不同白噪聲序列，依次重復步驟a、bN次。

d.對分解得到的各個IMF分量均值作為最終結果，最終得到的EEMD分解后的IMF為

(2)

e.選取合適的IMF分量進行重構數據。

1.2 核主成分分析

核主成分分析屬于非線性主元分析方法，采用非線性的方法提取主成分。

首先對時域、頻域、小波包能量譜所構成的樣本進行非線性變化，將其映射到一個高維特征空間F，協方差矩陣為

(3)

式中：xj為輸入空間的樣本(j= 1,2,…,m)；Φ(x)為在特征空間中的樣本點；C為協方差矩陣；m為輸入空間的樣本個數。

以上矩陣的特征向量和特征值滿足以下條件：

λV-CV=0

(4)

式中：V為協方差矩陣特征向量；λ為協方差矩陣特征值。

引入非線性函數Φ(xj)，可得：

λΦ(xj)V-Φ(xj)CV=0

(5)

式 (5)中得特征向量V可由Φ(x)線性表示，即：

(6)

將式(4)和式(6)代入式(5)，引入核函數Kij=K(xi,xj)=Φ(xi)Φ(xj),簡化后有：

mλα-Kα=0

(7)

式中：K為核矩陣；α為核矩陣K的特征向量。

對于任意樣本，在特征空間F中主元Φ(x)上的投影為

(8)

1.3 相關系數指標

分解結果中IMF分量容易產生虛假分量，因此引入相關系數辨別剔除虛假分量，設y1=(y11+y12+…+y1n)T，y2=(y21+y22+…+y2n)T為2個數據向量。則y1(t)、y2(t)的相關系數定義為

(9)

通過式(9)可知，計算各IMF與原始數據之間的相關系數，絕對值越大表示在該IMF分量中所包含的原始信息量越大。峭度指標是一個無量綱參數，它對于振動沖擊信號非常敏感，通過峭度指標可以反映軸承振動信號中的沖擊成分，因此適用于對機械系統進行故障診斷。當滾動軸承零部件在正常運行的過程中，對采集后的振動信號進行峭度值的計算。

經過峭度指標和有效系數的辨別篩選，選取合適的IMF分量重構信號。在重構信號中選取合適的特征值能夠反映軸承的退化性能評估指標。

2 軸承壽命預測模型的建立

2.1 改進SVM預測模型

在支持向量機的基礎上引入最小二乘線性理論就是最小二乘支持向量機，它是SVM的改進算法，具有和標準支持向量機互不相同的約束條件，它是用等式約束代替SVM中的不等式約束，將不等式約束的二次規劃問題轉化成線性矩陣求解問題，與標準的支持向量機相比，改進的LSSVM求解速度得到改善，精度得到提高。LSSVM的函數估計問題如下。

(10)

此時約束條件：

s.t.yi-ωΦ(xj)-b=0

(11)

式中：ω為權向量；γ為懲罰參數；b為偏差向量；ei為松弛變量，表示第i個數據的預測輸出和實際輸出的誤差值。

式(11)和式(12)對應的拉格朗日函數為

b+ei-yi]

(12)

式中：αi為Lagrange乘子，依據KKT條件進行求解，即對式(12)中的ω、b、e、αi求偏導數，并令其偏導數為零。

消掉變量ω、e，得到方程組如下。

(13)

yn]T。

σ為核函數寬密度。

最終得到LSSVM預測回歸模型為

(14)

建立LSSVM輸入x=(xn-m,xn-m+1,…,xn-1)與輸出y=(xn)一種非線性的映射關系：Rm→R,m為輸入函數，然后得到LSSVM的訓練樣本和測試樣本對為

式中：l為訓練樣本的對數；s為訓練數據的起點；n為預測點數。輸入為Xtrain，輸出為Ytrain，從而訓練出預測模型，求出LSSVM中預測函數(14)的具體數學表達式，將測試樣本Xtest帶入到式(14)中可以得到對應的n個預測值。

2.2 改進哈里斯鷹優化算法(IHHO)

2.2.1 HHO優化算法

HHO算法[12]是Heidari在2019年受哈里斯鷹捕食獵物啟發所提出的一種新型智能優化算法，主要由搜索階段、搜索與開發轉化階段和開發階段組成。

哈里斯鷹棲息在搜索空間[lb,ub]觀察某位置圍捕獵物；在迭代時根據其他個體與獵物的位置以概率p選擇更新。

(15)

式中：X(t+1)和X(t)表示第t+1次和t次迭代時個體的位置；Xt為獵物時的位置；Xd為個體的位置；ri(i=1,…,4)和p均為[0,1]內的一個隨機數，Xm(t)為t次迭代時的N只個體表達平均位置。

HHO算法根據獵物的逃逸能量在搜索和不同的開發行為之間轉換，用能量E線性動態遞減調控算法的全局搜索和局部開發性能，其定義式為

(16)

式中：E0為獵物的初始能量，為[-1,1]之間的隨機數；T為最大迭代數。當|E|≥1，HHO進入全局搜索階段；當|E|<1時，HHO算法則進入局部開發過程。

在HHO算法中，通過因子β∈[0,1]來描述獵物是否逃脫成功：當β<0.5時表示獵物逃脫成功，反之則失效：同時哈里斯鷹根據獵物能量E與0.5的相對大小來采取不同的圍捕策略，包含軟圍攻、硬圍攻、漸進式快速俯沖圍捕和漸進式快速俯沖硬捕4種情形。

Case1：軟圍攻，當β≥0.5且|E|≥0.5，表示獵物的能量E充足，并且以隨機跳躍的方式逃脫圍捕，但最終失敗并被捕獲，哈里斯鷹的位置更新為

X(t+1)=ΔXt(t)-E|JXt(t)-Xt(t)|

(17)

式中：ΔXt為第t次迭代時鷹與獵物的位置偏差；J=2(1-r5) 為獵物逃跑時隨機跳躍能力，r5∈(0,1)的一個隨機數。

Case2：硬圍攻，當β≥0.5且 |E|<0.5，表示獵物的能量E較低而被鷹直接捕獲，位置更新為

X(t+1)=Xt(t)-|EΔX(t)|

(18)

Case3：漸進式快速俯沖圍捕，當β<0.5且|E|≥0.5，表示獵物能量E充足可保證逃脫成功，但鷹以最優方向俯沖并軟圍捕獵物，其位置更新為

(19)

Y=Xt(t)-E|JXt(t)-Xt(t)|

(20)

Z=Y+S+LF(D)

(21)

式中:D為問題維度；S為D維隨機向量；LF為Levy函數。

Case4：漸進式快速俯沖硬捕，當β<0.5且|E|<0.5，表示獵物能量E較低，而鷹俯沖硬圍捕獵物，其位置更新式為

(22)

2.2.2 改進搜索階段

為了促進搜索階段的全局搜索，式(8)哈里斯鷹移動取決于獵物和其他個體的位置，而對于IHHO，第1種策略中(p<0.5)，其他個體鷹捕食獵物的移動被認為是模仿哈里斯鷹曲折運動飛行，采用Levy飛行模擬。而在第2種策略中(p>0.5)，哈里斯鷹個體移動在內部是隨機的。因此為了改善這一階段，將哈里斯鷹與其他3只隨機鷹初始值的平均位置來升級更新位置。

(23)

式中：Xj、Xk、Xl為哈里斯鷹隨機初始選擇的位置；Xm(t)為獵物的位置；LF為Levy函數。當概率因子p<0.5，哈里斯鷹的位置基于Levy飛行；當概率因子p≥0.5，哈里斯鷹遵循式(23)中第2行公式，隨機更改內部的位置。

2.2.3 改進能量調控機制

在HHO算法中，獵物能量E的大小反映哈里斯鷹對模型最優解的搜尋捕獲能力，并不能有效描述哈里斯鷹經多輪協同圍捕獵物的情況[13]，本文提出一種能量周期性調控機制，實現了“全局+局部”的尋優搜索。以余弦函數表示獵物能量E的周期遞變性，定義式為

(24)

式中：k=0,1,2,…，為獵物能量E的周期遞減數。

3 IHHO-LSSVM模型流程圖

結合IHHO 算法和改進LSSVM預測模型，其IHHO-LSSVM的流程如圖1所示。

圖1 基于IHHO-LSSVM的流程

選取KPCA降維后的軸承退化特征作為輸入向量，軸承剩余壽命值作為輸出向量，以預測壽命均方誤差(MSE)最小化為IHHO尋優目標，以核參數δ和正規化參數C為優化參數。

(25)

4 基于軸承全壽命預測的驗證

4.1 試驗裝置

采用來自美國凱西西儲大學電氣工程實驗室的滾動軸承全壽命周期加速軸承性能退化數據。試驗平臺如圖2所示。

圖2 滾動軸承全壽命試驗平臺

軸承試驗臺上的主軸上安裝了4個Rexnord-ZA-2115型號的雙列滾動軸承。試驗過程中，每個軸承將承受26 670 N的徑向載荷，每個軸承的X和Y方向各安裝了1個加速度傳感器，采樣頻率為20 kHz，每隔10 min采集1次滾動軸承振動信號，一共采集984個樣本，每個樣本長度20 480個點。共20 152 320個數據，軸承運行了1周，直到軸承退化失效，如圖3所示。

圖3 軸承原始信號

4.2 基于EEMD的軸承信號降噪

采集的信號中一般含有噪聲，不處理將會對軸承的性能退化指標產生較大影響。文中采用EEMD方法分解、篩選、重構，得到重構信號(每隔10個數據抽取1個，每個樣本長度2048個點)。具體操作：對原始振動信號進行EEMD計算，其中EEMD算法設置為Rstd=0.2，添加高斯白噪聲次數為1000，分解得到如圖4所示的21個IMF向量，計算21個IMF向量的峭度值和相關系數。

從圖4可知，IMF1-IMF12振動沖擊信號較多，而后續IMF缺少信息，利用相關系數對12個IMF分量進行篩選，前12個相關系數和峭度值的計算值，如表1所示。

圖4 各個IMF分量

表1 前12個相關系數和峭度值的計算

利用相關系數和峭度值的計算，峭度值偏離越大，說明軸承的故障信息越多。因此選取了1、3、4、5、6、7、8、9總計8個作為重構信號分量。重構信號中沖擊的成分更明顯，同時高頻和部分低頻中噪聲被去除。

4.3 驗證結果與分析

針對軸承時域，頻域中每個特征所呈現出的退化能力各不相同，選取能夠較好地反映出軸承在全壽命周期過程中退化趨勢的特征，剔除在全壽命周期過程中幾乎沒有任何變化的特征。因此，在時域退化特征中，篩選了均方根值、標準差、峰峰值為時域特征，均值頻率、標準差為頻域特征?？偣策x取了22個退化特征，利用核主成分分析方法(KPCA方法)將特征退化指標加權融合，提取KPCA的第一主成分作為性能退化評估指標，既能夠得到軸承性能退化趨勢，又使維數得到簡化。將KPCA第一主成分進行標準化處理，結果如圖5所示。

圖5 KPCA第一主成分

對比KPCA第一主成分退化趨勢和時域、頻域、小波包能量譜可以發現，KPCA第一主成分在500點(早期故障)后就已經呈現優良的上升趨勢，而單一的時域、頻域、小波包能量譜則不具有這種優勢。因此，利用KPCA加權融合提取主成分呈現軸承退化性能趨勢比單一的特征要全面，既有良好的上升趨勢，又對軸承故障早期具有敏感性。將KPCA各個成分的貢獻率從大到小依次排列，選取累計貢獻率達到85%的前n個主成分，如表2所示。

表2 部分KPCA的貢獻率

表2中，前3個主成分已經達到85%，故選取軸承退化期的第1、2、3個主成分和對應的剩余壽命組成訓練樣本，建立IHHO-LSSVM壽命預測模型，利用500點后的數據訓練，等間隔取7個作為預測樣本模型輸入，將LSSVM的核參數δ和正規化參數C通過改進的哈里斯鷹算法來進行優化，設置種群空間維度dim=2，最大迭代次數maxgen=100,設置哈里斯鷹種群數目sizepop=30，將未改進的HHO算法與IHHO算法的適應度值(適應度值選擇均方誤差)曲線迭代進行對比，如圖6所示。

由圖6可知，IHHO算法在收斂精度上提高了許多，為了對比和評估預測軸承壽命預測效果，將未采用EEMD-KPCA處理的IHHO-LSSVM模型和采用EEMD-KPCA處理的LSSVM模型進行壽命預測對比。如圖7所示，各個模型的預測誤差對比如表3所示。

圖6 2種方法尋優LSSVM參數的適應度值收斂趨勢

(a)EEMD-LSSVM模型

(b)IHHO-LSSVM模型

(c)EEMD-KPCA和IHHO-LSSVM模型圖7 各個模型壽命預測結果

表3 各個模型的預測誤差對比

由圖7可知，縱坐標表示每個樣本對應的剩余壽命時間序列，橫坐標表示軸承的使用時間，3種模型預測的軸承剩余壽命與真實壽命變化趨勢一致。對照圖7、表3可知，利用EEMD降噪去除噪聲、分解、重構信號，用KPCA融合前3個主成分作為軸承性能退化指標，最后采用IHHO-LSSVM作為壽命預測模型的3類誤差均小于前2種，說明EEMD分解重構具有消除噪聲的影響；KPCA消除了時域、頻域等各特征指標間的相關性導致的冗余數據；IHHO-LSSVM模型預測準確率對于其他支持向量機模型有所提高。

5 結論

a.采用EEMD-KPCA方法提取軸承性能退化指標，用 EEMD降噪重構信號，針對時域、頻域及小波包能量譜各特征指標間相關性而導致的數據冗余現象，提出采用KPCA提取主成分。該指標能夠兼顧軸承早期故障敏感性和全面表征軸承退化情況，為建立預測模型提高了良好的基礎。

b.提出了一種改變搜索階段和能量調控機制的新型哈里斯鷹群智能算法，與原算法相比，更容易找到LSSVM預測模型的核參數δ和正規化參數C的最優解。

c.建立一種新的IHHO-LSSVM滾動軸承壽命預測模型，與EEMD-LSSVM和HHO-LSSVM預測模型相比，滾動軸承的壽命預測精度得到了提高，為滾動軸承壽命預測提供了一種新的方法，確保電力設備全壽命周期健康管理，保證生產正常連續運行。