王 璽,周 薇,胡昌華,張建勛,裴 洪,劉 軒

(1.海裝駐北京地區(qū)第一軍事代表室, 北京 100076； 2.戰(zhàn)略支援部隊某部， 北京 102209； 3.火箭軍工程大學導彈工程學院， 西安 710025)

隨著現代科學技術的飛速發(fā)展，復雜性已經成為工程設備的一個重要發(fā)展趨勢[1]。然而，由于外部環(huán)境和內部因素的相互作用，工程設備的性能會隨著運行的時間而逐漸產生退化，并在達到某種量級時，導致設備的失效。如果不能準確掌握設備的健康狀態(tài)，在運行中一旦發(fā)生突發(fā)失效，可能會對社會安全和人身財產造成難以預計的后果。因此，針對設備的退化特性，建立相應的退化模型，進而由監(jiān)測的退化數據實現對退化設備的剩余壽命(remaining useful life，RUL)預測是目前可靠性領域的一項研究熱門。

目前，基于退化建模的剩余壽命預測方法已經成功的應用在激光器[2]、陀螺儀[3-4]、綜合傳動裝置[5-6]等領域。例如，文獻[2]建立了同時考慮測量誤差和隨機效應的退化模型，利用經典Kalman濾波理論實現了激光器的剩余壽命預測。文獻[3]同樣地考慮了隨機效應對陀螺儀剩余壽命預測的影響，并且基于多組同類型陀螺儀的歷史退化數據實現了模型的參數估計。文獻[5]在退化建模中考慮了監(jiān)測數據的測量不確定性，進一步提高了綜合傳動裝置的剩余壽命預測精度。

根據現有文獻可知，工程設備在性能退化過程中普遍存在著個體差異不確定性和非線性[7-14]。其中，個體差異不確定性又稱為隨機效應，是指由于受到不同的外力影響，同類型退化設備在性能退化過程中存在的退化差異。因此，有必要在退化建模中充分地考慮個體差異不確定性對剩余壽命預測的影響。為此，文獻[10-12]展開了相關的研究工作，基于Wiener過程進行退化建模，并將模型中的漂移系數作為隨機參數，用以描述個體差異不確定性。通過利用經典Kalman濾波理論實現了退化設備的剩余壽命預測，并且證明了考慮個體差異不確定性和非線性可以有效地提高剩余壽命預測的精度。但是，大多數相關工作在參數估計中通常假定存在多組同類型設備的歷史退化數據，致使剩余壽命預測的精度受限于數據量[2,5-9,11,14]。然而，在實際中經常會遇到缺乏此類歷史退化數據的情況，如陀螺儀、綜合傳動裝置、航空發(fā)動機等，造價昂貴且數量少，工程上難以通過大量設備進行壽命試驗[6]。因此，以上方法具有一定的局限性，且目前相關的研究工作較少。

鑒于此，本文基于非線性Wiener過程建立了具有雙隱含狀態(tài)的隨機退化模型，并利用粒子濾波算法實現了對隱含狀態(tài)的實時估計，克服了經典Kalman濾波理論的限制，進而在首達時間的概念下推導出了剩余壽命的分布。針對缺乏歷史退化數據和先驗信息等問題，提出了一種基于粒子期望最大化(particle expectation maximization，PEM)算法的參數估計方法，實現了模型參數的自適應估計和剩余壽命分布的在線更新。最后，通過某型陀螺儀的實際退化數據驗證了本文方法的有效性和優(yōu)越性。

1 模型描述

考慮采用具有一般性的非線性Wiener過程來描述復雜退化設備的隨機退化過程{X(t),t≥0}。具體地，{X(t),t≥0由標準布朗運動(standard Brownian motion，SBM)驅動，表示為：

(1)

基于首達時間(first hitting time，FHT)的概念[9-15]，認為退化過程{X(t),t≥0} 首次穿過設定的失效閾值ω時，設備的壽命終止。因此，將設備的壽命T定義為：

T=inf{t∶X(t)≥ω|X(0)<ω}

(2)

式(2)中：inf表示求集合的下確界，即首達時間；ω為預先設定的失效閾值，一般由工業(yè)標準確定。

在以上的框架下，假定在離散時刻點0=t0

Lk=inf{lk>0∶X(lk+tk)≥ω}〗

(3)

(4)

2 基于粒子濾波的剩余壽命預測

2.1 基于PF的雙隱含狀態(tài)估計

(5)

(6)

(7)

(8)

然而，在執(zhí)行PF算法時會遇到一個嚴重的粒子退化現象，即隨著迭代次數的增加，只有少量粒子的重要性權重比較大，大部分粒子的重要性權重很小。因此，粒子重要性權重的方差隨著運行時間不斷增大，狀態(tài)空間中的有效粒子數較少，使得大部分計算時間浪費在更新重要性權重小的粒子上，造成估計性能的下降[19-20]。為了衡量粒子的退化程度，引入有效抽樣尺度，定義為：

(9)

式(9)中：Neff越小，意味著粒子退化現象越嚴重。通常，抑制粒子退化的重要手段是利用重采樣算法，通過對粒子的重新采樣，剔除重要性權重小的粒子，“繁殖”重要性權重大的粒子，從而降低粒子退化現象。一般，設定重采樣閾值為Nth，當Neff

圖1 重采樣的示意圖

(10)

下面，將基于雙隱含狀態(tài)的后驗分布，實現退化設備的剩余壽命預測。

2.2 考慮個體差異不確定性的剩余壽命預測

根據式(3)中定義的剩余壽命，在考慮個體差異不確定性的情況下，剩余壽命的PDF可以計算如下：

fLk|Θ,X1∶k(lk|Θ,X1∶k)=Ezk|Θ,X1∶k[fLk|zk,Θ,X1∶k(lk|zk,Θ,X1∶k)]

(11)

式(11)中：Ezk|Θ,X1∶k[·]是在給定模型參數向量Θ和退化數據X1∶k的條件下，關于zk的條件數學期望算子。

為了計算式(11)，首先給出以下結果[16]。

(12)

基于PF算法，聯合式(11)和式(12)，在tk時刻預測退化設備的剩余壽命的PDF可以近似計算為：

(13)

基于以上結果，在獲得新的退化數據后，可以通過式(13)預測退化設備的剩余壽命。但需要注意的是，上述結果是在模型參數向量Θ給定的條件下推導的，因此必須先對模型參數向量Θ進行參數估計。然而，現有大多數相關工作存在著一個共性問題，即假定存在多組同類型設備的歷史退化數據，基于極大似然方法離線地估計出模型參數，并且這些參數一旦確定，就不再隨著新的退化數據進行在線更新[2,5-9,11,14]。此外，在實際中經常會遇到缺乏此類歷史退化數據的情況，例如新列裝的武器裝備、新型工程設備等。鑒于此，本文提出了一種基于PEM算法的參數估計方法，利用實時監(jiān)測的退化數據來實現模型參數的自適應估計和剩余壽命分布的在線更新，使得預測的剩余壽命更加準確，有利于掌握退化設備的健康狀態(tài)。

3 基于PEM算法的自適應參數估計

3.1 總體框架

根據極大似然估計(maximum likelihood estimation，MLE)理論，由于隨機參數λk和θk為雙隱含狀態(tài)，因此需要基于可觀測的退化數據X1∶k來估計出模型參數向量Θ。這里，自然可以采用一種廣泛應用于缺失或不完全數據條件下的EM算法來實現模型的參數估計。EM算法是通過最大化聯合似然函數p(zk,X1∶k|Θ)去逼近模型參數的MLE，具體過程通過迭代以下兩步來實現。

E步：

(14)

M步：

(15)

首先，基于貝葉斯定理和建立的退化模型(4)，E步中的對數聯合似然函數可以通過式(16)形式進行計算，即：

lnp(zk,X1∶k|Θ)=lnp(X1∶k|zk,Θ)+lnp(zk|Θ)=

(16)

然后，通過E步對式(16)求條件數學期望，有：

(17)

式(17)中：

3.2 基于粒子平滑的雙隱含狀態(tài)平滑

(19)

p(zm|Θ,X1∶k)=p(zm|Θ,X1∶m)×

(20)

(21)

因此，基于式(20)和式(21)，可以得到：

(22)

到此，基于PS算法解決了非線性平滑問題。下面，將給出本文所提出的PEM算法的具體實現過程。

3.3 基于PEM算法的實現過程

基于PS算法，可以直接得到式(17)中的Λ1和Λ3，結果如下：

(23)

(24)

(25)

(26)

式(26)中：

(27)

基于以上的結果，可以將本文所提出的剩余壽命自適應預測方法歸納如下：

(3) M步。根據式(26)和式(27)，計算

如果無法給出解析解，則通過搜索算法進行求解。

(4) 迭代。如果EM算法已經收斂，則終止此步驟，否則令q=q+1，并返回步驟(2)。

(5) 剩余壽命預測。

② 基于式(5)，在tk時刻預測退化設備的剩余壽命的PDF可以由式(13)計算。

由此可見，本文所提出的剩余壽命自適應預測方法，可以利用每一個退化數據對模型參數進行自適應地估計，進而實現剩余壽命分布的在線更新，克服了非線性退化設備缺乏歷史退化數據和先驗信息的問題，這對工程實際具有重要意義。

4 實例驗證

固定在慣性平臺上的陀螺儀是航空航天及導彈武器系統(tǒng)的關鍵部件，其工作性能直接影響了導航的精度。當慣性平臺工作時，陀螺儀的轉子高速旋轉勢必會造成旋轉軸機械磨損，隨著磨損的積累，表現為漂移系數增大，最終導致設備失效。而且對于安裝在慣性平臺上的陀螺儀，數量少且測試成本高，一旦發(fā)生失效，所造成的后果將難以估計。因此，準確估計陀螺儀的剩余壽命，對提高整個系統(tǒng)的安全性和可靠性至關重要。

下面通過某型陀螺儀在實際使用中監(jiān)測的退化數據來驗證本文方法的有效性。如圖2所示，在實驗中該型陀螺儀根據技術指標，設定失效閾值ω=0.3(°/h)，每隔2 h監(jiān)測一次，直到372 h為止，共監(jiān)測了187次。這里，選擇以下冪函數形式的隨機退化模型進行建模，即：

X(t)=X(0)+λ·tθ+σBB(t)

(28)

到目前為止，此類模型已經在陀螺儀、航天發(fā)動機、鋰電池等復雜設備的退化建模中得到了廣泛的應用[8-13]。

圖2 陀螺儀漂移的退化軌跡

從圖2中可以看出，隨著監(jiān)測時間的增加，陀螺儀的漂移系數整體呈上升趨勢。基于該陀螺儀的實際退化數據，應用所提出的PEM算法可以實現模型參數的自適應估計，結果如圖3所示。從圖3中可以看出，所提方法可以在每一個退化數據可用時估計模型的參數，而且估計的模型參數會隨著退化數據的積累快速收斂。

為了進一步驗證本文方法的有效性和優(yōu)越性，選擇文獻[11]和文獻[12]中的方法進行對比研究。

1) M1：假定存在多組同類型設備的歷史退化數據，基于極大似然的方法離線地估計出模型參數，并且這些參數一旦確定，就不再隨著新的退化數據進行更新[11]。

2) M2：在退化建模中忽略了非線性函數中參數的隨機性，利用經典的Kalman濾波理論實現了漂移系數的實時估計，并且在新的退化數據可用時，可以實現模型參數的自適應估計和剩余壽命分布的在線更新[12]。

圖3 模型參數的自適應估計曲線

在實驗中，為了驗證本文方法可以有效地應用于缺乏先驗信息的退化設備，進行如下設定：本文方法和M2方法選用隨機的模型初始參數，而M1方法選用合適的模型初始參數。圖4給出了三種方法在監(jiān)測時間300～360 h時，每隔5個監(jiān)測時間點預測的剩余壽命的PDF對比圖。從圖4中可以看出，三種方法預測的剩余壽命的PDF曲線均覆蓋了實際的剩余壽命，并且隨著監(jiān)測時間的增加，PDF曲線愈來愈尖銳，這意味著預測的剩余壽命結果愈來愈準確，不確定性愈來愈低。此外，由于本文方法考慮了非線性函數中參數的隨機性，并且在每一個監(jiān)測時間點都可以實現模型參數的自適應估計和剩余壽命分布的在線更新，因此本文方法得到的PDF曲線更加尖銳，說明預測的不確定性更低，相比之下具有更好的預測性能。

圖4 三種方法預測的剩余壽命的PDF曲線

為了進一步量化剩余壽命預測的精度，引入可靠性領域中常用的均方誤差(mean square error，MSE)指標，該指標同時考慮了剩余壽命預測的精度和不確定性，因此在剩余壽命預測中被廣泛應用。具體地，在tk時刻的定義式如下：

(29)

圖5給出了三種方法在所有監(jiān)測時間點的MSE結果。從圖5中可以看出，由于本文方法和M2方法選用隨機的模型初始參數，因此在監(jiān)測初期退化數據較少時，得到的結果并不理想。但是，隨著監(jiān)測時間的增加，退化數據的積累，本文方法的MSE明顯低于M1方法和M2方法，說明本文方法在預測剩余壽命的性能上具有明顯的優(yōu)勢。

圖5 三種方法在所有監(jiān)測時間點的MSE曲線

以上的實驗結果表明，本文方法可以顯著地提高剩余壽命預測的精度，降低預測的不確定性，在實際應用中具有更好的預測性能，并且更加適應于缺乏歷史退化數據和先驗信息的非線性退化設備。

5 結論

1) 本文在狀態(tài)空間模型的框架下建立了具有雙隱含狀態(tài)的非線性隨機退化模型，并且在經典Kalman濾波理論不適用的情況下，基于粒子濾波算法實現了對雙隱含狀態(tài)的實時估計。

2) 基于粒子濾波算法，在首達時間的概念下推導出了剩余壽命的分布，并且在獲得新的退化數據時，可以實現剩余壽命分布的在線更新。

3) 提出了一種基于PEM算法的參數自適應估計方法，克服了非線性退化設備缺乏歷史退化數據和先驗信息的問題，實現了模型參數的自適應估計。