基于支持向量數(shù)據(jù)描述的剩余壽命預(yù)測方法

武千惠，黃必清

(清華大學 自動化系，北京 100084)

0 引言

設(shè)備故障不僅產(chǎn)生了大量的維修費用，降低了生產(chǎn)效率和經(jīng)濟效益，還會引發(fā)安全事故，造成人員傷亡。為了提高設(shè)備可靠性、保證人員安全、提高企業(yè)生產(chǎn)效率，確定維修保障資源、做出合理的維修決策十分重要。預(yù)測性維修是目前最為理想的維修策略[1]。預(yù)測性維修指在設(shè)備運行時對目標部件進行狀態(tài)監(jiān)測，判斷運行狀態(tài)并預(yù)測未來的發(fā)展趨勢，依據(jù)可能的故障時間制定維修計劃，安排維修活動。預(yù)測性維修的核心在于故障預(yù)測，通?；谑Ｓ鄩勖A(yù)測技術(shù)實現(xiàn)[2-3]。

目前，故障預(yù)測方法主要分為：基于模型的預(yù)測方法，基于經(jīng)驗的預(yù)測方法和數(shù)據(jù)驅(qū)動的預(yù)測方法[1,4]3大類。基于模型的方法通過應(yīng)用物理定律建立一組代數(shù)方程或微分方程作為預(yù)測模型，利用該模型表示退化現(xiàn)象(如變形、裂紋和腐蝕)。由于設(shè)備部件的退化過程通常都是非線性的，且與不確定的工況密切相關(guān)，因此該類方法的應(yīng)用場合十分受限。基于經(jīng)驗的預(yù)測方法利用經(jīng)驗反饋數(shù)據(jù)(如失效時間)調(diào)整一些分析模型(如威布爾分布、指數(shù)分布)的參數(shù)，通過這些參數(shù)模型估計系統(tǒng)的失效時間或剩余壽命，但往往存在預(yù)測精度不夠理想、過度依賴知識和經(jīng)驗等問題。數(shù)據(jù)驅(qū)動的故障預(yù)測方法不需要依靠物理參數(shù)(如工況、材料性質(zhì))，就可以通過傳感器狀態(tài)監(jiān)測數(shù)據(jù)建立預(yù)測模型。與基于模型和基于經(jīng)驗的方法相比，數(shù)據(jù)驅(qū)動的預(yù)測方法在預(yù)測準確度、計算復(fù)雜度和可推廣性等方面進行了折中，具有更強的可操作性。本文研究內(nèi)容即為數(shù)據(jù)驅(qū)動的故障預(yù)測方法。

研究人員在數(shù)據(jù)驅(qū)動的故障預(yù)測方面做了許多工作，主要包括人工神經(jīng)網(wǎng)絡(luò)(Artificial Neural Network, ANN)[5-7]、隱馬爾可夫模型(Hidden Markov Model, HMM)[1,8]、貝葉斯信念網(wǎng)絡(luò)(Bayesian Belief Network, BBN)[9]、支持向量機(Support Vector Machine, SVM)[10-11]和支持向量數(shù)據(jù)描述(Support Vector Data Description, SVDD)[12-14]等方法。例如，神經(jīng)網(wǎng)絡(luò)方面，Huang等[7]以軸承的退化期為研究對象，首先利用自組織映射(Self-Organizing Maps, SOM)神經(jīng)網(wǎng)絡(luò)生成的最小量化誤差(Minimum Quantization Error, MQE)訓練BP神經(jīng)網(wǎng)絡(luò)，然后通過失效權(quán)重適用(Weight Application to Failure Times, WAFT)技術(shù)構(gòu)建了剩余壽命預(yù)測模型。馬爾可夫模型方面，Tobon-mejia等[1]提出一種基于高斯混合隱馬爾可夫模型的預(yù)測方法，模型的隱狀態(tài)代表設(shè)備在失效過程中經(jīng)歷的退化狀態(tài)，剩余壽命可由每個退化狀態(tài)的逗留時間直接計算得到。貝葉斯置信網(wǎng)絡(luò)方面，Zhang等[9]提出一種利用高斯混合貝葉斯置信網(wǎng)絡(luò)識別軸承退化狀態(tài)，由退化狀態(tài)識別結(jié)果隨時間的變化預(yù)測剩余壽命的方法。由于設(shè)備的故障數(shù)據(jù)十分有限，基于神經(jīng)網(wǎng)絡(luò)、隱馬爾可夫模型和貝葉斯置信網(wǎng)絡(luò)的預(yù)測方法在實際應(yīng)用受到了很多限制，主要表現(xiàn)為容易出現(xiàn)過擬合、模型穩(wěn)定性較低、難以得到比較理想的預(yù)測結(jié)果等。

支持向量方面，陳斌等[12]針對轉(zhuǎn)子振動臺的不平衡故障預(yù)測問題，為不同的故障程度建立了不同的SVDD超球體，通過待測樣本到各超球體球心的相對距離確定其所屬的故障等級。因為隨著故障等級的增加，相對距離的方向具有指向性，所以可將其作為預(yù)測不平衡故障程度的準則。Shen等[14]考慮軸承性能退化及動態(tài)信息的隨機性和強模糊性，提出一種基于模糊SVDD和運行時間的單調(diào)性能退化評估指標，以有效刻畫軸承損傷情況隨工作時間的變化。Kim等[10]以軸承為研究對象，將其退化過程分為6個階段，利用“一對一”多分類SVM方法進行故障分類，通過訓練數(shù)據(jù)的分類情況求得軸承處于各個階段的概率，最后利用軸承的工作時長和概率值預(yù)測剩余壽命。已有的基于支持向量的預(yù)測方法雖然能夠?qū)﹄x散的目標部件的工作狀態(tài)進行判斷及預(yù)測，但存在難以預(yù)測連續(xù)的剩余壽命時長(如文獻[12-14])或剩余壽命預(yù)測精度不夠(如文獻[10-11])等問題，應(yīng)用于運維策略優(yōu)化等后續(xù)工作時存在困難。針對以上問題，本文提出一種基于SVDD的剩余壽命預(yù)測方法。

1 基本方法

1.1 支持向量數(shù)據(jù)描述

給定數(shù)據(jù)集{xi,i=1,2,…,N}，SVDD的優(yōu)化目標為尋找在特征空間包含所有樣本點的半徑最小的超球面[15]。令c為超球面的中心，R為超球面半徑，該優(yōu)化問題可以表示為如式(1)所示結(jié)構(gòu)風險最小化問題：

s.t.

(xi-c)T(xi-c)≤R2,i=1,2,…,N。

(1)

當訓練樣本中存在遠離其他普通樣本點的奇異點時，利用式(1)得到的超球面因包含大片空白區(qū)域而不能準確描述訓練樣本。為了提高算法魯棒性，引入松弛變量ξi和懲罰因子C，以在最小超球面半徑和最小非目標樣本分類誤差間尋求折中。C越大，上述奇異點對目標函數(shù)的影響越大，此時原問題轉(zhuǎn)換為：

s.t.

(xi-c)T(xi-c)≤R2+ξi,ξi≥0,i=1,2,…,N。

(2)

構(gòu)造Lagrange函數(shù)

αi≥0,γi≥0;

(3)

對L求偏導(dǎo)并令其等于0，有

(4)

將式(4)帶入式(3)，得到僅包含Lagrange乘子的原問題的對偶問題：

s.t.

(5)

為了使超球面對樣本點的包裹形式更加靈活，引入核函數(shù)將樣本點映射到高維空間。將式(5)中的內(nèi)積(xi·xj)用核函數(shù)K(xi·xj)替代：

s.t.

(6)

通過求解式(6)可以得到各αi的值，與0<αi

設(shè)f(z)為樣本z到超球面球心距離的平方，則

f(z)=‖z-c‖2=K(z·z)-

(7)

假設(shè)xsv為支持向量，由定義可知超球面半徑R即為支持向量到球心的距離，化簡后有

(8)

1.2 粒子群優(yōu)化算法原理

粒子群優(yōu)化(Particle Swarm Optimization, PSO)算法是Kennedy和Eherhart于1995年提出的一種智能優(yōu)化算法[16]。與遺傳算法等進化算法類似，PSO算法首先初始化為一群隨機粒子，每個粒子對應(yīng)一個可行解，然后通過不斷迭代找到優(yōu)化問題的最優(yōu)解。在每次迭代中，粒子通過跟蹤個體極值和全局極值更新自己在搜索空間中的位置和速度。

設(shè)在一個D維空間中有一個規(guī)模為m的粒子群，第i個粒子的位置為xi=(xi1,xi2,…,xiD)T，速度為vi=(vi1,vi2,…,viD)T。每個粒子的適應(yīng)度值由適應(yīng)度函數(shù)給出，由此可以確定粒子的個體極值pi=(pi1,pi2,…,piD)T和種群的全局極值pg=(pg1,pg2,…,pgD)T。在每次迭代中，粒子在個體極值和全局極值的影響下不斷優(yōu)化自己的位置，并以此影響其他粒子，使整個種群朝搜索空間中最優(yōu)解的方向前進。粒子的速度和位置分別按式(9)和式(10)更新：

(9)

(10)

式中：加速度因子c1和c2為非負常數(shù)，r1和r2為服從[0,1]均勻分布的隨機數(shù)，k為迭代次數(shù)，ω為慣性權(quán)重。慣性權(quán)重決定了粒子在第k次迭代時的速度對第k+1次迭代時的速度的影響程度，具有平衡算法的全局搜索能力和局部搜索能力的作用。Shi和Eberhart于1998年提出了線性減少慣性權(quán)重的自適應(yīng)策略[17]，即

(11)

線性減少慣性權(quán)重可以使PSO算法在收斂的同時，具有更強的跳出局部極值的能力。式中：ωmax和ωmin分別表示權(quán)重系數(shù)的最大值和最小值，kmax為最大迭代次數(shù)。

為了加快算法收斂速度，本文對粒子位置進行了限制。設(shè)粒子位置的最大值和最小值分別為xdmax和xdmin，則在位置更新后需要對其進行如下修正：

(12)

2 基于支持向量數(shù)據(jù)描述的剩余壽命預(yù)測方法

圖1所示為基于PSO-SVDD的剩余壽命預(yù)測流程，包括兩個階段：

(1)離線階段 首先通過小波包分解對歷史振動數(shù)據(jù)進行特征提取，然后利用PSO算法選擇合適的核函數(shù)參數(shù)，針對設(shè)備健康狀態(tài)的特征向量訓練SVDD模型，得到超球面的半徑和中心。

(2)在線階段 從實時傳感器信號提取到特征向量后，利用從離線階段得到的SVDD超球面信息求出待測設(shè)備的退化指數(shù)，根據(jù)退化的取值判斷目標部件的健康狀態(tài)，并在其進入非健康狀態(tài)后預(yù)測剩余壽命。

2.1 特征提取

小波包分解(Wavelet Packet Decomposition, WPD)是一種通過多次迭代的小波變換對輸入信號的細節(jié)部分逐步展開分析的方法，其迭代過程用一棵樹表示：根節(jié)點為原始信號，根節(jié)點的下一層是對原始數(shù)據(jù)集做了一次小波變換后的結(jié)果，之后的各層均由對上一次小波變換得到的低通和高通結(jié)果遞歸進行小波變換求得，分解過程完成后便可計算末層不同波帶的能量系數(shù)。通過迭代分解，WPD不但能夠?qū)⒋治鲂盘柕念l域進行逐級劃分，而且能夠以不同信號的特點為依據(jù)選擇與其信號頻譜匹配的頻帶，從而同時提高頻域和時域的分辨率。與離散小波變換相比，應(yīng)用小波包分解進行信號分析能夠得到更多細節(jié)信息，得到的頻域特征也更易于進行后續(xù)選擇和分類。

因為小波包分解具有足夠的高頻分辨率，而設(shè)備部件振動信號中的故障信息大多包含于該頻段[18]，所以本文采用WPD進行特征提取。假設(shè)分解層數(shù)為l，則最后一層的節(jié)點數(shù)為L=2l。令fit表示最后一層第i個節(jié)點的能量，則t時刻的特征向量可以表示為

ft=(f1t,f2t,…fLt)T。

(13)

為了提高系統(tǒng)的魯棒性，應(yīng)對WPD結(jié)果進行Z-score標準化處理，再將標準化后的結(jié)果作為特征進行后續(xù)處理與分析。

2.2 PSO-SVDD模型

設(shè)備部件在其使用階段會經(jīng)歷由健康到失效的轉(zhuǎn)變過程。健康狀態(tài)的特征總是相似的，而隨著退化的產(chǎn)生及發(fā)展，與健康狀態(tài)的特征相比，新的特征在特征空間中的差異越來越顯著，分布也越來越分散。

第1章已經(jīng)說明，SVDD模型的建立準則為找到在特征空間包含所有樣本點的半徑最小的超球面。因此，如果用設(shè)備部件健康狀態(tài)的歷史數(shù)據(jù)訓練得到SVDD超球體，則待測樣本點與超球體球心間的距離可作為反映目標部件健康狀態(tài)及退化發(fā)展程度的指標。

以SVDD作為健康狀態(tài)評估和剩余壽命預(yù)測的基本模型，首先需要確定核函數(shù)類型及參數(shù)。與多項式核函數(shù)相比，本文采用高斯核函數(shù)，不但能夠?qū)崿F(xiàn)非線性映射，而且參數(shù)數(shù)量較少，能夠有效降低模型的復(fù)雜度。高斯核函數(shù)的表達式為

(14)

由于高斯核函數(shù)的參數(shù)σ對SVDD描述樣本數(shù)據(jù)的準確性影響很大，且對于不同的樣本集，很難預(yù)先確定合適的σ取值。為此，本文提出一種利用PSO選擇合適的核函數(shù)參數(shù)并由此訓練SVDD超球面的方法。

定義退化指數(shù)

(15)

DI表征退化程度的大小，DI越大，目標部件退化越嚴重。

假設(shè)目標部件于t1時刻進入非健康狀態(tài)，于t2時刻完全故障停機。因為設(shè)備部件進入非健康狀態(tài)后的故障特征隨時間呈指數(shù)增長[19-20]，所以由更優(yōu)核函數(shù)參數(shù)計算所得的DI值隨時間的變化應(yīng)該更加接近指數(shù)函數(shù)?；谝陨戏治?，本文利用DI(t1)和DI(t2)構(gòu)造指數(shù)函數(shù)g(t)=aeb(t-t1)，進而定義適應(yīng)度函數(shù)

(16)

核函數(shù)參數(shù)的PSO算法流程如圖2所示。

得到PSO算法優(yōu)化的核函數(shù)參數(shù)后，利用SMO算法即可通過歷史數(shù)據(jù)的特征訓練SVDD超球面[15]。對于每個訓練集，計算包含設(shè)備部件全生命周期的DI取值，由DI的變化情況求得該訓練集進入非健康狀態(tài)的時間t1和完全故障停機的時間t2，進而通過DI(t1)和DI(t2)求出表征該訓練集退化過程的指數(shù)函數(shù)，為在線階段預(yù)測剩余壽命提供依據(jù)。

2.3 退化檢測和剩余壽命預(yù)測

在在線階段，設(shè)t時刻由傳感器監(jiān)測數(shù)據(jù)經(jīng)特征提取后得到的特征向量為z，則可由式(17)求得此時該設(shè)備部件的退化指數(shù)DI(t)。DI的取值越大，設(shè)備部件的退化越嚴重，剩余壽命越短。

一種退化檢測是直接將DI的零點作為設(shè)備退化的起始時刻，當DI≤0時，則說明目標部件處于健康狀態(tài)；反之，若DI>0，則說明目標部件進入非健康狀態(tài)。由于實際測量所得振動信號存在噪聲，可選擇μ+3σ作為退化閾值，當檢測窗口中DI的均值達到該閾值時，定義當前時刻為退化起始時刻，并從此時刻開始預(yù)測剩余壽命。

圖3所示為剩余壽命的預(yù)測方法，其中指數(shù)曲線由訓練集的DI(t1)和DI(t2)求得。

由于被測部件在進入非健康狀態(tài)后，其特征向量的變化幅度很大，為了消除噪聲和不穩(wěn)定因素對預(yù)測精度的影響，最終輸出的剩余壽命預(yù)測值由該時刻和此前(win_size-1)個時刻的剩余壽命的初始預(yù)測值經(jīng)移動平均求得，以達到抑制預(yù)測值的大幅波動、更好地反映剩余壽命真實值及其變化趨勢的目的。其中，win_size為移動窗口大小。

3 實驗案例分析

本章的實驗數(shù)據(jù)來自FEMTO-ST研究院搭建的PRONOSTIA平臺。PRONOSTIA是一個為軸承故障檢測、故障診斷及故障預(yù)測研究服務(wù)的實驗平臺，能夠提供球軸承從投入使用至失效的全部測量數(shù)據(jù)，平臺情況如圖4所示[21]。PRONOSTIA平臺主要由旋轉(zhuǎn)部分、負載部分和測量部分3部分構(gòu)成。旋轉(zhuǎn)部分通過電動機、變速系統(tǒng)和聯(lián)軸器等帶動軸承支撐軸旋轉(zhuǎn)，進而使內(nèi)圈與軸承支撐軸固定的軸承產(chǎn)生旋轉(zhuǎn)運動。負載部分向測試軸承施加徑向壓力，加速軸承壽命周期。測量部分通過鉑熱電阻測量軸承工作溫度，通過微加速度計測量軸承在水平方向和垂直方向的振動信號。

本文以數(shù)據(jù)集[22]中軸承B1_1和B1_3的水平方向振動信號作為訓練數(shù)據(jù)和測試數(shù)據(jù)，對提出的剩余壽命預(yù)測方法進行驗證。其中，采樣頻率為25.6 kHz，采樣間隔為10 s，轉(zhuǎn)速為1 800 rpm，負載為4 000 N[21]。實驗首先通過WPD從訓練數(shù)據(jù)和測試數(shù)據(jù)中提取特征。將小波包基設(shè)為db4，分解層數(shù)設(shè)為3，得到WPD分解樹末層8個節(jié)點的能量系數(shù)如圖5和圖6所示。

圖中不同曲線表示不同節(jié)點的能量系數(shù)隨時間的變化情況。可以看出，實驗初期各節(jié)點的能量系數(shù)很小，且能夠在較長時間內(nèi)保持平穩(wěn)，說明軸承處于健康狀態(tài)的時間在其整個生命周期內(nèi)所占比例很大。隨著時間的推移，能量系數(shù)開始增長并發(fā)生顯著變化，說明軸承進入非健康狀態(tài)，且故障程度逐步加深。因此，WPD樹末層節(jié)點的能量系數(shù)能夠反映軸承故障的演變情況，是一種有效的特征提取方法。

在利用PSO算法優(yōu)化SVDD核函數(shù)參數(shù)σ時,取搜索空間維度D=1，其他參數(shù),即種群大小m，最大迭代次數(shù)kmax，加速度因子c1和c2，權(quán)重系數(shù)最大值ωmax和最小值ωmin，粒子位置區(qū)間xmax和xmin，停止條件閾值threshold的具體設(shè)定如表1所示，得到如圖7所示的優(yōu)化結(jié)果。

表1 PSO參數(shù)取值

圖7說明限制粒子活動范圍后，PSO算法在迭代不到50次后即停止搜索，得到全局最優(yōu)解p=42.96。實驗還發(fā)現(xiàn)，當不限制粒子位置時，PSO算法在迭代70余次后滿足終止條件，停止搜索，可見限制粒子活動范圍能夠有效減小迭代次數(shù)，加快算法收斂速度。

圖3和圖8所示分別為訓練軸承和測試軸承的退化指數(shù)DI隨時間的變化情況。實驗初期DI的取值與波動幅度很小，說明軸承處于健康狀態(tài)，且健康狀態(tài)的時長在其整個生命周期中占有很大比例；進入退化期后，軸承的退化指數(shù)取值開始逐漸增大且產(chǎn)生較大波動，總體呈現(xiàn)出類指數(shù)函數(shù)的增長趨勢，說明被測軸承的故障特征在加速發(fā)展，與軸承投入使用后的客觀規(guī)律一致。因此，本文定義DI的取值能夠準確反映軸承健康狀態(tài)的變化情況，DI的取值越大，軸承的故障程度越深。

RMSE=

(17)

表2 剩余壽命預(yù)測結(jié)果對比

σ39.0042.9646.00RMSE18.635 512.449 414.913 8

圖9～圖11和表2表明，本文所提方法能夠有效評估軸承的健康狀態(tài)并預(yù)測其剩余壽命，特別是經(jīng)PSO算法選擇核函數(shù)參數(shù)后，軸承剩余壽命的預(yù)測準確度得到了顯著提升，比較結(jié)果說明PSO算法在優(yōu)化SVDD核函數(shù)的選擇時具有一定的效果和優(yōu)勢。

需要說明的是，因為設(shè)備部件的故障數(shù)據(jù)十分有限，而文獻[5-7]中提出的ANN模型通常需要更多的訓練數(shù)據(jù)才能得到相對比較理想的模型，所以在訓練數(shù)據(jù)不足時易于產(chǎn)生較大預(yù)測誤差而無法實際應(yīng)用。本文提出的剩余壽命預(yù)測算法算法參數(shù)較少，易于實現(xiàn)，且SVDD模型在訓練數(shù)據(jù)偏少的情況下與神經(jīng)網(wǎng)絡(luò)模型相比更加穩(wěn)定，因此具有更高的實際應(yīng)用價值。另外，文獻[8-12]中健康狀態(tài)的識別結(jié)果是離散的，剩余壽命由離散的健康狀態(tài)識別結(jié)果估計得到，因此具有一定的階躍特性，而當訓練數(shù)據(jù)與測試數(shù)據(jù)的來源部件工作總時長相差較大時，由于目標部件在每個健康狀態(tài)的持續(xù)時間難以確定，會出現(xiàn)較大的預(yù)測誤差。由圖3和圖8可知，本文定義的DI不但能夠表征設(shè)備健康狀態(tài)的變化，而且在進入退化期后具有類似指數(shù)函數(shù)的變化趨勢，因此能夠通過DI的變化趨勢得到更加連續(xù)、更加準確的剩余壽命預(yù)測結(jié)果。

4 結(jié)束語

基于SVDD模型,本文定義了能夠表征設(shè)備健康狀態(tài)的DI，并由此提出一種剩余壽命預(yù)測方法。首先利用WPD技術(shù)，從歷史傳感器狀態(tài)監(jiān)測數(shù)據(jù)中提取包含設(shè)備部件故障信息的特征向量；其次，引入PSO算法，求得能夠使訓練集DI取值的變化趨勢更加接近指數(shù)規(guī)律的核函數(shù)參數(shù)；然后利用由PSO算法求得的核函數(shù)參數(shù)和設(shè)備健康狀態(tài)的歷史數(shù)據(jù)訓練SVDD模型，通過待測樣本點和SVDD超球面信息計算描述待測設(shè)備退化程度的DI；最后，通過DI隨時間的變化檢測退化起始時刻，利用退化期故障特征與時間的近似指數(shù)關(guān)系預(yù)測剩余壽命。

本文通過軸承加速壽命實驗數(shù)據(jù)對所提方法進行了驗證，實驗結(jié)果表明，該方法能夠有效預(yù)測軸承剩余壽命，而且利用PSO算法能夠避免人為選擇SVDD核函數(shù)參數(shù)時的盲目性，使得DI取值隨時間的變化能夠更好地描述目標部件退化過程，從而有效提升預(yù)測準確率。