李 瑞， 汪立新， 劉 剛

（第二炮兵工程大學，陜西西安 710025）

0 引 言

慣性平臺是導彈控制系統(tǒng)的關(guān)鍵設(shè)備，其性能好壞直接關(guān)系到導彈武器系統(tǒng)的可靠性和制導精度，因此，需要對慣性平臺進行實時準確的健康狀態(tài)預測和管理［1］。工程實踐表明，預測與健康管理（Prognostics and Health Management，PHM）是提高設(shè)備可靠性和安全性、減小失效風險、降低維護費用的一項重要新興技術(shù)［2］。PHM包括兩層含義：一是預測，通過對設(shè)備的故障診斷和健康狀態(tài)監(jiān)測，預測設(shè)備的剩余壽命；二是健康管理，根據(jù)診斷／預測信息、可用資源、使用要求等信息，做出最優(yōu)的維護決策。其中，基于狀態(tài)監(jiān)測的剩余壽命預測是后期確定最優(yōu)維護時機、優(yōu)化監(jiān)測間隔、制定備件訂購策略以及提供延壽依據(jù)的關(guān)鍵［3］。因此，通過對慣性平臺的健康狀態(tài)監(jiān)測，及時準確地預測剩余壽命，對提高整個導彈武器系統(tǒng)的可靠性、安全性和打擊精度具有十分重要的意義。

鑒于慣性平臺價格昂貴的特點，難以得到足夠多的失效數(shù)據(jù)，使得傳統(tǒng)的以歷史失效數(shù)據(jù)為基礎(chǔ)的方法不再適用。因此，依據(jù)設(shè)備的退化機理（如機械零件磨損、疲勞裂紋擴展、絕緣材料老化等）建立退化規(guī)律模型，進而估計其剩余壽命成為一條經(jīng)濟可行的途徑［4－5］。由于陀螺漂移是慣性平臺的主要失效形式，占整個慣性平臺失效的70%左右，因此通常采用陀螺漂移系數(shù)表征慣性平臺的健康狀態(tài)［1］。漂移系數(shù)的測試值越大，平臺系統(tǒng)的工作性能越壞，健康狀態(tài)越差，當超過規(guī)定的技術(shù)指標后，慣性平臺失效，不能繼續(xù)使用。由于漂移系數(shù)的增長過程具有增加或減小的非單調(diào)特性，而現(xiàn)有的退化建模方法，如Wiener過程、Gamma過程、Markov鏈、隱含馬氏過程等，只有Wiener過程具有非單調(diào)的特性［6］。因此，現(xiàn)有針對慣性平臺的退化建模和剩余壽命預測都采用了基于Wiener過程的方法［7－9］。

Wiener過程的非單調(diào)和無窮可分特性，使其能夠更好地描述系統(tǒng)的動態(tài)特性，并在多種設(shè)備的退化建模和剩余壽命預測中有著廣泛的應(yīng)用［7－18，19］，如慣性平臺［7－9］、激光發(fā)生器［11］、鋰電池［12］、梁橋［13］、滾動軸承［14－16］、連續(xù)攪拌器［17］?；赪iener過程的離線剩余壽命預測方法已經(jīng)相對成熟，研究熱點主要集中在剩余壽命的在線預測。文獻［14-16］首次采用貝葉斯方法實現(xiàn)了基于Wiener過程方法的實時剩余壽命預測，該方法需要較多的歷史退化數(shù)據(jù)作為先驗信息，不適用于慣性平臺這種小子樣場合。文獻［17］將Kalman算法與期望最大化算法結(jié)合起來實時更新Wiener過程中的漂移參數(shù)，并考慮了Wiener過程首達時間的概念，然而該方法同樣需要較多的先驗信息。文獻［18］在此基礎(chǔ)上考慮了漂移參數(shù)估計不確定性對剩余壽命預測的影響，并采用強跟蹤濾波算法和期望最大化算法實現(xiàn)了所有參數(shù)的在線更新，從而減小了對先驗信息的依賴。然而，現(xiàn)有的這些方法沒有考慮測量誤差對剩余壽命預測的影響，而測量誤差普遍存在于退化過程的狀態(tài)檢測中，需要在建模過程中加以考慮［9－12］。從文中分析可知，考慮測量誤差能夠降低預測的不確定性，因此，文中將測量誤差引入到對慣性平臺的退化建模中，得出了考慮測量誤差時剩余壽命的解析解，建立了性能退化的狀態(tài)空間預測模型，并采用Kalman算法與期望最大化算法結(jié)合，實時更新模型參數(shù)。最后通過慣性平臺的退化測量數(shù)據(jù)驗證了文中提出的方法。

1 基于Wiener過程的退化建模

文中主要考慮線性退化過程的剩余壽命預測問題。令X（t）表示某設(shè)備的實際退化過程，則基于Wiener過程的退化模型可表示為

式中：λ——漂移系數(shù)；

σB——擴散系數(shù)；

x0——初始狀態(tài)；

B（t）——標準布朗運動，表征衰退過程的動態(tài)特性。

從式（1）可得

可以看出，退化過程的速度與漂移系數(shù)λ密切相關(guān)，退化過程的不確定性隨著時間的增長而增長，可以描述設(shè)備失效產(chǎn)生機理和動態(tài)運行環(huán)境的變化。由于退化過程主要受漂移系數(shù)λ的影響，在離線預測中，通常假設(shè)λ為隨機變量，用來表示不同設(shè)備之間的個體差異［14－16］。對于在線預測，漂移系數(shù)λ也是影響預測結(jié)果的主要因素，因而有必要考慮漂移系數(shù)λ參數(shù)估計的不確定性對評估結(jié)果的影響。類似于文獻［8，18］，文中也把漂移系數(shù)λ看成隨機變量，以表示參數(shù)估計過程中的不確定性。

在狀態(tài)測試過程中，由于受到儀器或環(huán)境噪聲的影響，真實的退化過程往往不能精確測量，觀測得到的測量值不可避免地存在測量誤差?？紤]測量是離散的，則在時間ti時的測量值可以表示為：

式中：yi——測量值；

εi——測量誤差，假設(shè)是獨立同分布，服從N（0，）；

εi——測量過程中的不確定性，可以由儀器標定數(shù)據(jù)和歷史數(shù)據(jù)統(tǒng)計得出。

為了減小實時預測對先驗信息的依賴，有必要通過實際的狀態(tài)監(jiān)測值實時更新模型參數(shù)。為此，構(gòu)建如下狀態(tài)空間模型：

其中，t0＝0，ζi～N（0，ti－ti－1），η服從正態(tài)分布，表示漂移參數(shù)隨著時間的變化特征。

然而，由于實際數(shù)據(jù)經(jīng)常出現(xiàn)xi＜xi－1的情況，容易造成漂移系數(shù)λ出現(xiàn)負值，從而導致預測性能變差。針對這些不足，文中將狀態(tài)空間模型修正如下：

其中，εi～N（0，σ2），η～N（0，Q），λi表示給定狀態(tài)真實值X0：i下的漂移系數(shù)，X0：i＝｛x0，x1，x2，…，xi｝。

顯然，以上模型為線性狀態(tài)空間模型，可以通過Kalman濾波在線估計參數(shù)，具體如下：

1）初始化參數(shù)：

2）狀態(tài)預測：

3）方差更新：

可以得出，在時刻ti，λi～N，Pi｜i）。

2 剩余壽命預測

考慮測量誤差情況下剩余壽命分布的求解?；赪iener過程進行壽命分析時，通常將壽命T定義為狀態(tài)監(jiān)測值首次到達失效閾值w的時間，即T＝inf｛t：X（t）≥w｜x0＜w｝。剩余壽命Lt則定義為Lt＝｛lt：T－t｜T＞t｝表示在使用一段時間t后的剩余壽命。對于進行定期狀態(tài)檢測的設(shè)備，剩余壽命的預測還需要利用當前時刻的狀態(tài)監(jiān)測值。設(shè)某個測試點的監(jiān)測數(shù)據(jù)為xi，則

顯然，隨著測試點的增多，系統(tǒng)動態(tài)特性的不確定性也會越來越小。在時刻ti，設(shè)xi＜w，剩余壽命可以表示為衰退過程首次達到w的剩余時間，即

式中：w——失效閾值。

為了得到考慮測量誤差情況下的剩余壽命分布，需要給出以下引理。

引理1 設(shè)截止當前時刻的歷史狀態(tài)真實值為X0：i，在時刻ti，漂移系數(shù)λi～N，Pi｜i），相應(yīng)的剩余壽命分布可表示為

上述定理的證明可以參考文獻［8］的定理5和文獻［18］的定理1，這里不再贅述。假設(shè)λi為一確定值，即Pi｜i＝0，則上述定理轉(zhuǎn)化為一般Wiener過程首達時間的逆高斯分布形式。由于引理1將漂移參數(shù)的不確定性考慮進來，預測精度得以進一步提高［8，18］。為推導文中模型下的剩余壽命分布，先給出以下引理，具體證明可以參考文獻［8］的定理4和文獻［12］的引理2。

引理2 設(shè)Z～N（μ，σ2），A，B，C為常數(shù)，則有

式中：E（）——表示求期望。

然而，由于真實的退化過程往往不能精確測量，只能得到觀測值，由式（2）可得xi＝y(tǒng)i－εi，則xi～N（yi，σ2），即只能得到真實狀態(tài)量的概率分布。據(jù)此，可以得到如下結(jié)論。

定理1 設(shè)截止當前時刻的歷史狀態(tài)觀測值為Y0：i，Y0：i＝｛y0，y1，y2，…，yi｝。在當前時刻ti，漂移系數(shù)λi～N，Pi｜i），相應(yīng)的剩余壽命分布可表示為

證明 由上文分析可知：

那么令

根據(jù)引理2對式（7）計算關(guān)于xi的期望就可得式（9）。證畢。

可以看出，當不存在測量誤差，即σ2＝0時，式（9）轉(zhuǎn)變?yōu)槭剑?），則文中結(jié)果是現(xiàn)有基于Wiener過程的在線剩余壽命預測方法的一般化推廣。

3 模型參數(shù)在線估計

考慮到部分設(shè)備具有價格昂貴、壽命長的特點，往往難以獲得足夠多的先驗信息，那么先驗參數(shù)Θ＝（a0，P0，Q，，σ2）會存在較大的不確定性。再者，先驗參數(shù)也包含了不同樣本的差異性，因而需要在線更新先驗參數(shù)，即出現(xiàn)新的狀態(tài)觀測值后，在更新漂移參數(shù)的同時，也實現(xiàn)先驗參數(shù)的在線更新。目前主要通過極大似然估計法在線更新先驗參數(shù)，在當前時刻ti，Y0：i對數(shù)似然函數(shù)為

式中：p（Y0：i｜Θ）——在給定Θ情況下Y0：i的概率分布函數(shù)。

由式（1）～式（3）可知：

又因為λ0～N（a0，P0），式（10）可以轉(zhuǎn)化為

可以通過最大化對數(shù)似然函數(shù)估計未知參數(shù)Θ的值，即

由于漂移參數(shù)Υ0：i＝｛λ0，λ1，λ2，…，λi｝是隱含的隨機變量，難以直接最大化對數(shù)似然函數(shù)來估計Θ。然而，期望最大化（Expectation Maximization，EM）算法是計算包含隱性變量概率模型中參數(shù)極大似然估計的常用迭代算法［21］。初始化參數(shù)Θ＝Θ（0）后，EM算法主要包括以下兩步：

1）E步，計算包含隱含變量的對數(shù)似然函數(shù)期望值，即計算式（11）關(guān)于漂移參數(shù)Υ0：i的期望值。從式（11）可知，計算式（11）期望值之前，首先需要計算λi，λi－1λi，的期望值，這些期望值可以通過Kalman平滑的方法計算［17－18，20］，即

其向后平滑步驟如下：

其中協(xié)方差Mj｜i的初始化值為Mi｜i＝（1－Kiti）Pi－1｜i－1。令Cj－1，j｜i＝E（λj－1λj｜Y0：i，Θ），Cj｜i＝E（｜Y0：i，Θ），則由式（13）可得對數(shù)似然函數(shù)關(guān)于漂移參數(shù)Υ0：i的期望

2）M步，極大化估計后的對數(shù)似然函數(shù)。先驗參數(shù)Θ中，與設(shè)備性能退化個體差異密切相關(guān)的是（a0，P0，Q，），而σ2表示測量過程當中的不確定性，可由儀器標定數(shù)據(jù)和歷史數(shù)據(jù)統(tǒng)計得出，因此在參數(shù)迭代過程中，只需要更新參數(shù)（a0，P0，Q，）。最大化式（14）得

將式（16）和式（17）計算得到的Θ（k）值重新代入到E步中，繼續(xù)重復E步和M步，直到求得似然函數(shù)的最大值。

4 慣性平臺剩余壽命預測實驗

4.1 慣性平臺陀螺漂移退化數(shù)據(jù)

慣性平臺在使用中，由于陀螺轉(zhuǎn)子的高速旋轉(zhuǎn)造成轉(zhuǎn)軸的磨損，以及隨著貯存時間的延長平臺系統(tǒng)自身性能的退化［1，9］，引起陀螺漂移系數(shù)的增長，進而影響整個導航系統(tǒng)的可靠性和精度。文中使用文獻［18］提供的某型號慣性平臺漂移系數(shù)退化數(shù)據(jù)進行算法比較研究。該平臺共進行了73次周期性工作測試，采樣間隔為2.5 h，具體的漂移監(jiān)測數(shù)據(jù)如圖1所示。

圖1 慣性平臺陀螺漂移退化數(shù)據(jù)

根據(jù)技術(shù)指標規(guī)定，失效閾值為0.37°／（h·g）。對于該慣性平臺，在時刻t73的監(jiān)測數(shù)據(jù)為0.356 6°／（h·g），在該時刻的剩余壽命估計為6 h，則實際剩余壽命估計為186 h。

4.2 實驗比較分析

從文獻［18］中的慣性平臺壽命預測試驗可以看出，文獻［14-16］提出的方法極大地依賴于先驗信息，然而慣性平臺性能退化的先驗信息相對較少，因而文中主要與文獻［17-18］提出的方法進行比較，記文中方法為M0，文獻［17］的方法為M1，文獻［18］的方法為M2。為了比較這些方法的優(yōu)劣，引入AIC準則和關(guān)于剩余壽命預測的總均方誤差（TMSE）作為評價方法好壞的比較。AIC準則用來評價模型擬合的優(yōu)良性，TMSE用來評價對剩余壽命預測的準確性，二者的計算公式如下：

式中：k——參數(shù)數(shù)量；

依據(jù)慣性平臺歷史退化信息，取先驗參數(shù)＝（2e－3，1e－4，1e－4，2e－4，3e－5），可得采用方法M0，M1和M2預測的剩余壽命概率密度函數(shù)，分別如圖2和圖3所示。

圖2 M0，M1的剩余壽命概率密度估計結(jié)果比較

圖3 M0，M2的剩余壽命概率密度估計結(jié)果比較

從圖2和圖3可以看出，3種方法的剩余壽命分布都能夠覆蓋實際剩余壽命，而文中方法計算的剩余壽命概率密度函數(shù)更窄，也更集中于實際剩余壽命。

為了進一步比較3種方法的效果，計算3種方法的95%置信區(qū)間，如圖4所示。

可以看出，得出的剩余壽命置信區(qū)間最窄，說明文中方法的不確定小，精度高。原因是方法M1只更新了漂移參數(shù)，沒有實時更新其它參數(shù)，而方法M2由于沒有考慮測量誤差的影響，使得退化過程的不確定性（σ2B）增大，從而增加了預測的不確定性。

圖4 M0，M1，M2的剩余壽命預測置信區(qū)間結(jié)果比較

文中方法與文獻方法的比較見表1。

表1 文中方法與文獻方法的比較

從表1可以看出，考慮測量誤差后，表示退化過程不確定性的σ2B值有所降低，預測不確定性也隨之降低。比較TMSE可知，文中算法的準確性都顯著優(yōu)于方法M1和M2。綜上所述，文中算法能夠降低剩余壽命預測的不確定性，提高預測精度，進而提高維護決策的置信度。

5 結(jié) 語

1）建立了考慮測量誤差情況下基于Wiener過程的線性退化模型，修正了用于漂移系數(shù)遞歸估計的狀態(tài)空間模型，并推導出了剩余壽命分布的解析表達式。

2）為了實現(xiàn)所提模型的參數(shù)實時更新，通過Kalman濾波算法實時估計Wiener過程中的漂移系數(shù)，并通過期望最大化算法實現(xiàn)了其它相關(guān)參數(shù)實時估計。

3）基于慣性平臺陀螺漂移退化數(shù)據(jù)的剩余壽命預測實驗表明，文中算法具有更好的建模能力和更高的剩余壽命預測準確性，具有一定的工程應(yīng)用價值。

［1］ 胡昌華，司小勝.基于信度規(guī)則庫的慣性平臺健康狀態(tài)預測參數(shù)在線估計［J］.航空學報，2010，31（7）：1454-1465.

［2］ Pecht M G.Prognostics and health management of electronics［M］.New Jersey：John Wiley，2008.

［3］ 張旭剛，張華，江志剛.基于剩余使用壽命評估的再制造方案決策模型及應(yīng)用［J］.機械工程學報，2013，49（7）：51-57.

［4］ 張小麗，陳雪峰，李兵，等.機械重大裝備壽命預測綜述［J］.機械工程學報，2011，47（11）：100-116.

［5］ 敖銀輝，汪寶生.鉆頭磨損檢測與剩余壽命評估［J］.機械工程學報，2011，47（1）：177-181.

［6］ Si X S，Wang W B，Hu C H，et al.Remaining use-ful life estimation：A review on the statistical data driven approaches［J］.European Journal of Operational Research，2011，213（1）：1-14.

［7］ Si X S，Wang W B，Hu C H，et al.Remaining useful life estimation based on a nonlinear diffusion degradation process［J］.IEEE Transactions on Reliability，2012，61（1）：50-67.

［8］ Si X S，Wang W B，Chen M Y，et al.A degradation path-dependent approach for remaining useful life estimation with an exact and closed form solution［J］.European Journal of Operational Research，2013，226（1）：53-66.

［9］ 司小勝，胡昌華，周東華.帶測量誤差的非線性退化過程建模與剩余壽命估計［J］.自動化學報，2012，38（1）：1471-1484.

［10］ Whitmore G A.Estimating degradation by a wiener diffusion process subject to measurement error［J］.Lifetime Data Analysis，1995，1（3）：307-319.

［11］ Peng C Y，Tseng S T.Mis-specification analysis of linear degradation models［J］.IEEE Transactions on Reliability，2009，58（3）：444-455.

［12］ Feng L，Wang H L，Si X S，et al.A state spacebased prognostic model for hidden and age-dependent nonlinear degradation process［J］.IEEE Transactions on Automation Science and Engineering，2013，10（4）：1072-1086.

［13］ Wang X.Wiener processes with random effects for degradation data［J］.Journal of Multivariate Analysis，2010，101：340-351.

［14］ Gebraeel N Z，Lawley M A，LI R，et al.Residual-life distributions from component degradation signals：A Bayesian approach［J］.IIE Transactions，2005，37：543-557.

［15］ Kaiser K A，Gebraeel N Z.Predictive maintenance management using sensor-based degradation models［J］.IEEE Transactions on Systems，Man，and Cybernetics Part A：Systems and Humans，2009，39（4）：840-849.

［16］ Elwany A H，Gebraeel N Z.Sensor-driven prognostic models for equipment replacement and spare parts inventory［J］.IIE Transactions，2008，40（6）：629-639.

［17］ Wang W B，Carr M，Xu W，et al.A model for residual life prediction based on Brownian motion with an adaptive drift［J］.Microelectronics Reliability，2011，51：285-293.

［18］ Si X S，Wang W B，Hu C H，et al.A wiener process based degradation model with a recursive filter algorithm for remaining useful life estimation［J］.Mechanical Systems and Signal Processing，2013，35（2）：219-237.

［19］ Sun B，Zeng S K，Kang R，et al.Benefits and challenges of system prognostics［J］.IEEE Transactions on Reliability，2012，61（2）：323-334.

［20］ 馮磊，王宏利，韓曉霞，等.考慮監(jiān)測噪聲的陀螺儀剩余壽命在線預測［J］.長春工業(yè)大學學報：自然科學版，2013，34（2）：155-159.

［21］ Dempster A P，Laird N M，Rubin D B.Maximum likelihood from incomplete data via the EM algorithm［J］.Journal of the Royal Statistical Society，Series B.，1977，39（1）：1-38.