劉　剛，黎　放

（1.海軍工程大學(xué)，武漢　430033；2.海軍工程大學(xué)艦船綜合保障工程研究室，武漢　430033）

?

基于貝葉斯模型的裝備剩余壽命預(yù)測(cè)研究

劉剛1，2，黎放1，2

（1.海軍工程大學(xué)，武漢430033；2.海軍工程大學(xué)艦船綜合保障工程研究室，武漢430033）

摘要：為了解決武器裝備日益復(fù)雜及維修工作日趨繁重的問題，運(yùn)用貝葉斯模型預(yù)測(cè)裝備修理后的剩余壽命，為合理安排其修理計(jì)劃提供依據(jù)。替代傳統(tǒng)的指數(shù)分布，用威布爾分布描述系統(tǒng)壽命特征，并運(yùn)用極大似然方法和貝葉斯方法估計(jì)威布爾分布的兩個(gè)未知參數(shù)，給出其置信區(qū)間。在此基礎(chǔ)上，對(duì)先驗(yàn)樣本和后驗(yàn)樣本兩種不同情況，分別運(yùn)用貝葉斯模型預(yù)測(cè)裝備修理后的剩余壽命，并給出實(shí)例。結(jié)果表明了該方法的有效性。

關(guān)鍵詞：威布爾分布，貝葉斯，剩余壽命，預(yù)測(cè)

0　引言

近年來，隨著科技水平的不斷發(fā)展，武器裝備現(xiàn)代化、模塊化程度不斷增加，采取傳統(tǒng)的維修模式已經(jīng)不能適應(yīng)新裝備的維修保障要求。對(duì)于裝備的使用人員或直接維修人員而言，傳統(tǒng)的預(yù)防性維修、視情維修策略越來越不能滿足裝備保障現(xiàn)代化保障需求。而采取基于狀態(tài)的維修（Condition-Based Maintenance，CBM）則能減少維修費(fèi)用，最大化利用裝備［1-2］。因此，對(duì)裝備的運(yùn)行狀況進(jìn)行監(jiān)控，提前發(fā)現(xiàn)裝備的故障征兆，采用合適的方式進(jìn)行修理，預(yù)測(cè)裝備的剩余壽命，這對(duì)于提高裝備可用性、降低維修消耗具有重要的理論和現(xiàn)實(shí)意義。

一般認(rèn)為復(fù)雜的機(jī)電產(chǎn)品壽命分布函數(shù)如同電子產(chǎn)品一樣，是服從指數(shù)分布的。但實(shí)踐證明，只有小部分復(fù)雜的機(jī)電產(chǎn)品壽命分布函數(shù)是服從指數(shù)分布的，大部分是服從威布爾分布的，通過威布爾分布的研究，指數(shù)分布只是其特殊形式。因此，對(duì)于武器裝備而言，運(yùn)用威布爾分布描述其壽命分布函數(shù)，更加合理、準(zhǔn)確。目前，威布爾分布被廣泛運(yùn)用于航天航空、機(jī)械制造、生物醫(yī)學(xué)等領(lǐng)域［3-5］。

因此，針對(duì)以上情況，本文研究運(yùn)用威布爾分布描述系統(tǒng)使用壽命，研究運(yùn)用極大似然方法和貝葉斯方法來求解威布爾分布的兩個(gè)未知參數(shù)，同時(shí)求得這兩個(gè)參數(shù)的貝葉斯置信區(qū)間；基于不同的分析樣本，運(yùn)用貝葉斯模型預(yù)測(cè)系統(tǒng)在修理后的剩余壽命，從而為使用和維修人員合理安排裝備的使用與維修計(jì)劃提供合理依據(jù)。

1　威布爾分布

威布爾分布的概率密度函數(shù)和分布函數(shù)如下所示：

其中，v是形狀參數(shù)，β是尺度參數(shù)。當(dāng)v＞1時(shí)，威布爾分布的值是遞減的；當(dāng)v＜1時(shí)，威布爾分布的值是遞增的；當(dāng)v=1，威布爾分布就成了指數(shù)分布；當(dāng)v≥3.5時(shí)，其曲線趨于正態(tài)分布的情況。因此，可以看出，威布爾分布是比指數(shù)分布、正態(tài)分布更一般的分布，可以更加準(zhǔn)確描述復(fù)雜系統(tǒng)的壽命特征。

為了研究的方便，令θ=1/βv，則式（1）和式（2）可表示為：

2　問題描述

假設(shè)系統(tǒng)使用壽命服從威布爾分布，內(nèi)嵌BIT監(jiān)控其運(yùn)行狀態(tài)，其中在運(yùn)行過程中，可能會(huì)發(fā)生n個(gè)故障，為了使系統(tǒng)運(yùn)行狀況最優(yōu)，其中m個(gè)故障可能需要提前修理。當(dāng)系統(tǒng)運(yùn)行到時(shí)間點(diǎn)x1時(shí)，r1個(gè)故障被修理；到時(shí)間點(diǎn)x2時(shí)，r2個(gè)故障被修理；直到時(shí)間點(diǎn)xm，rm=n-r1-r2-…rm-1-m個(gè)故障被修理為止。因此，修理的時(shí)間點(diǎn)集X=（x1，x2，…，xm），對(duì)應(yīng)時(shí)間點(diǎn)故障數(shù)量集為R=（r1，r2，…，rm）。在本文中，針對(duì)該系統(tǒng)，研究運(yùn)用極大似然方法和貝葉斯方法估計(jì)威布爾分布的參數(shù)值，并確定其置信區(qū)間，將其用于預(yù)測(cè)系統(tǒng)的剩余壽命。

3　參數(shù)估計(jì)

在運(yùn)用威布爾分布的過程中，由于θ和υ這兩個(gè)參數(shù)都是未知的，無法直接利用。一般情況下，可以利用系統(tǒng)運(yùn)行過程中積累的觀測(cè)數(shù)據(jù)來估計(jì)未知的參數(shù)，因此，文中基于這種思路，運(yùn)用極大似然方法來分析威布爾分布的兩個(gè)參數(shù)，同時(shí)利用貝葉斯模型分析威布爾分布的平方誤差損失和置信區(qū)間。

3.1極大似然估計(jì)

設(shè)系統(tǒng)運(yùn)行過程中的觀測(cè)數(shù)據(jù)為（x，r），則其似然函數(shù)為：

對(duì)于極大似然函數(shù)，分別對(duì)θ和v求導(dǎo)，并使求導(dǎo)后的公式為0，得到：

將式（9）代入式（8），可以得到的極大似然估計(jì)可以通過下式求出：

通過式（11）～式（13）可以得到矩陣I（θ，v），在一般情況下，（，）的極大似然估計(jì)近似于均值為（θ，v），協(xié)方差為I-1（θ，v）的正態(tài)分布，即：

其中I0（θ，v）為以下矩陣：

3.2貝葉斯估計(jì)

極大似然估計(jì)所利用的信息一般都是已有的信息，俗稱先驗(yàn)信息，這些信息反映的是系統(tǒng)歷史運(yùn)行狀況，雖然一定程度上能反映系統(tǒng)未來的趨勢(shì)，但是不夠全面；而貝葉斯估計(jì)除了能利用先驗(yàn)信息以外，還能模擬產(chǎn)生系統(tǒng)未來的運(yùn)行信息，俗稱后驗(yàn)信息，將先驗(yàn)信息同后驗(yàn)信息結(jié)合起來，共同估計(jì)θ和v的參數(shù)值，得到的參數(shù)值能更加全面反映系統(tǒng)的壽命特征［6-9］。

假設(shè)θ和v有如下先驗(yàn)形式：

其中a，b，c，d都是正實(shí)數(shù)。

通過式（5）和式（14），可以得出參數(shù)θ和υ的先驗(yàn)分布如下：

其中：

同時(shí)分別求出θ和v后驗(yàn)邊緣分布：

基于平方誤差損失，參數(shù)θ和v的貝葉斯估計(jì)值為：

可以看出，估計(jì)值中含有復(fù)雜的積分函數(shù)Φ（a，c，x），因此，考慮運(yùn)用林德利近似法來計(jì)算貝葉斯估計(jì)值，計(jì)算公式如下：

考慮到有2個(gè)參數(shù)的情況，令u（λ）=（λ1，λ2），其平均值近似為：

其中：

此處，ij為先驗(yàn)信息矩陣的倒置矩陣的第i行第j列的值，L為先驗(yàn)數(shù)據(jù)的極大似然函數(shù)。式（20）中未知參數(shù)（λ1，λ2）的極大似然估計(jì)值為（1，2）。

令u（λ）=（θ，v），則有：

其中Iij為先驗(yàn)信息矩陣的第i行第j列的值。

根據(jù)林德利近似法，參數(shù)θ和v的貝葉斯估計(jì)值為：

3.3貝葉斯置信區(qū)間

假設(shè)函數(shù)g（θ，v）的置信水平為100（1-γ）%，其貝葉斯置信區(qū)間為（Lg，Ug），則（Lg，Ug）需滿足以下條件：

實(shí)際應(yīng)用中一般用γ/2代替γ。

通過式（16）和前面的計(jì)算，假設(shè)參數(shù)v的置信水平為100（1-γ）%，則其置信區(qū)間（Lv，Uv）可以通過下式求解：

其中，

同理，θ的置信區(qū)間（Lθ，Uθ）可以通過下式求解：

其中，

4　基于多樣本的貝葉斯預(yù)測(cè)

貝葉斯估計(jì)可以同時(shí)利用先驗(yàn)信息和后驗(yàn)信息，因此，在樣本的選取中，基于先驗(yàn)信息的樣本可以簡(jiǎn)稱為先驗(yàn)樣本，基于后驗(yàn)信息的樣本可以簡(jiǎn)稱為后驗(yàn)樣本［10-14］。在本節(jié)，將分別討論基于先驗(yàn)樣本和后驗(yàn)樣本來開展系統(tǒng)運(yùn)行時(shí)間的貝葉斯預(yù)測(cè)。

4.1基于先驗(yàn)樣本的預(yù)測(cè)

令Ys表示系統(tǒng)更換備件后的運(yùn)行時(shí)間，其中Ys=Xm+s，s=1，2，…，rm。當(dāng)X=x時(shí)，Ys的概率密度函數(shù)為：

其中，R（·）=1-F（·|θ，v）為服從威布爾分布的生存函數(shù)，κ為歸一化常數(shù)，滿足

通過式（23）和式（15），可以得到基于先驗(yàn)樣本的分布函數(shù)為：

因此，基于先驗(yàn)樣本預(yù)測(cè)的生存函數(shù)為：

因此，通過對(duì)先驗(yàn)樣本分析，及求解下列公式，可以得到置信水平為100（1-γ）%時(shí)，Ys的置信區(qū)間為（L1，U1）。

P（Ys＞U1|x）=γ/2，P（Ys＞L1|x）=1-γ/2

4.2基于后驗(yàn)樣本的預(yù)測(cè)

在很多工程應(yīng)用領(lǐng)域，分析人員一般希望在現(xiàn)有數(shù)據(jù)的基礎(chǔ)上，對(duì)系統(tǒng)未來的運(yùn)行狀況進(jìn)行預(yù)測(cè)，即運(yùn)用來自相同的壽命時(shí)間分布的一組規(guī)模為N的后驗(yàn)樣本預(yù)測(cè)修理后系統(tǒng)的運(yùn)行時(shí)間。令Ti（i=1，2，…，N，T1＜…＜TN）為基于規(guī)模為N的后驗(yàn)樣本的系統(tǒng)運(yùn)行時(shí)間，如式（3）所示。s個(gè)順序統(tǒng)計(jì)量的概率密度函數(shù)如下所示：

通過式（25）和式（15），可以得到基于后驗(yàn)樣本的分布函數(shù)為：

相應(yīng)地，其生存函數(shù)為：

因此，通過后驗(yàn)樣本分析，求解下列公式，可以得到置信水平為100（1-γ）%時(shí)，Ts的置信區(qū)間為

5　算例

為了驗(yàn)證本文提出的方法的有效性，對(duì)某型設(shè)備開展具體試驗(yàn)分析，下面是該型設(shè)備在不同條件作用下記錄下來的數(shù)據(jù)，該數(shù)據(jù)服從威布爾分布［15］。從觀測(cè)樣本n=19中隨機(jī)挑選規(guī)模m=8的樣本，具體數(shù)據(jù)如表1所示。

利用式（9）和式（10）進(jìn)行計(jì)算，可以分別求出θ和v的極大似然估計(jì)值=0.114 8，=0.974 3。假設(shè)（θ，v）的先驗(yàn)分布是式（14）所描述的兩變量分布，令參數(shù)a=0.1，b=0.5，c=1.5，d=15。利用式（18）和式（19）進(jìn)行計(jì)算，可以分別求出θ和v的貝葉斯估計(jì)值B=0.116 7，B=0.973 6。利用式（21）和式（22）進(jìn)行計(jì)算，運(yùn)用林德利近似法，分別求出θ和v的貝葉斯估計(jì)值L=0.112 1，L=0.984 9。同時(shí)，可以求得在置信水平為95%時(shí)，θ和v的置信區(qū)間分別是（0.0370，0.245 8）和（0.515 5，1.527 7）。

表1　試驗(yàn)中挑選的樣本數(shù)據(jù)

置信水平為95%時(shí)，在先驗(yàn)樣本的數(shù)據(jù)條件下，系統(tǒng)的貝葉斯估計(jì)剩余壽命分布和壽命時(shí)間間隔。通過先驗(yàn)樣本，得到剩余壽命分布函數(shù)，從該分布中生成一組規(guī)模N=5的后驗(yàn)樣本。通過對(duì)后驗(yàn)樣本進(jìn)行分析，s（1≤s≤5）個(gè)故障被修理后，系統(tǒng)剩余壽命的貝葉斯估計(jì)間隔預(yù)測(cè)與實(shí)際值的對(duì)比如下頁表2所示。

可以清晰地看到，隨著s的增加，估計(jì)的間隔長(zhǎng)度也隨之增大，即系統(tǒng)的剩余壽命也越長(zhǎng)；通過對(duì)比分析，根據(jù)后驗(yàn)樣本預(yù)測(cè)的系統(tǒng)的剩余壽命更接近實(shí)際值，證明了后驗(yàn)樣本是先驗(yàn)樣本的有力補(bǔ)充，對(duì)兩者綜合分析更能準(zhǔn)確預(yù)測(cè)系統(tǒng)剩余壽命。

表2　預(yù)測(cè)結(jié)果與實(shí)際值對(duì)比

6　結(jié)論

相比指數(shù)分布，威布爾分布能更加真實(shí)地反映復(fù)雜系統(tǒng)的壽命特征，運(yùn)用威布爾分布描述系統(tǒng)壽命時(shí)間，使得對(duì)系統(tǒng)剩余壽命的預(yù)測(cè)能夠更加貼近實(shí)際，便于使用和維修人員安排合理的修理計(jì)劃。本文運(yùn)用極大似然方法和貝葉斯模型估計(jì)威布爾分布的兩個(gè)未知參數(shù)，并給出這兩個(gè)參數(shù)的貝葉斯置信區(qū)間，為威布爾分布在復(fù)雜系統(tǒng)剩余壽命預(yù)測(cè)中的應(yīng)用打下基礎(chǔ)。然后，在基于先驗(yàn)樣本和后驗(yàn)樣本的不同情況下，運(yùn)用貝葉斯模型對(duì)該系統(tǒng)剩余壽命進(jìn)行了預(yù)測(cè)研究。最后開展了實(shí)例分析，結(jié)果表明，該方法與實(shí)際試驗(yàn)結(jié)果更加接近，與傳統(tǒng)基于先驗(yàn)樣本的方法相比，該方法是有效的。

參考文獻(xiàn)：

［1］楊志波，董明.動(dòng)態(tài)貝葉斯網(wǎng)絡(luò)在設(shè)備剩余壽命預(yù)測(cè)中的應(yīng)用研究［J］.計(jì)算機(jī)集成制造系統(tǒng)，2007，13（9）：1811-1815.

［2］孫紹輝，王華偉，李偉.潛在故障期內(nèi)航空發(fā)動(dòng)機(jī)的剩余壽命預(yù)測(cè)［J］.航空計(jì)算技術(shù)，2012，42（1）：8-11.

［3］WANG W B.A two-stage prognosis model in condition based maintenance［J］.European Journal of Operational Research，2007，12（182）：1177-1187.

［4］NAGI G.Mark lawley residual life predictions from vibration-based degradation signals：a neural network approach a neural network approach［J］.IEEE Transactions on Industrial Electronics，2004，51（3）：694-700.

［5］YUEN H K，TSE S K. Parameters estimation for weibull distributed lifetimes under progressive censoring with random removals［J］. Statist Comput Simul，1996，3（55）：57-71.

［6］FERNáNDEZ A J. On estimating exponential parameters with general type II progressive censoring［J］. Statist Plann Inference，2004，1（121）：135-147.

［7］SOLIMAN A A. Estimation of parameters of life from progressively censored data using burr-XII model［J］. IEEE Trans Reliab，2005，5（54）：34-42.

［8］ASGHARZADEH A. Point and interval estimation for a generalized logistic distribution under progressive type II censoring［J］. Commun Statist Theory Meth，2006，3（35）：1685-1702.

［9］KUS C，KAYA M F. Estimation for the parameters of the pareto dsitribution under progressive censoring［J］. Commun Statist Theory Meth. 2007，2（36）：1359-1365.

［10］MOUSA M A，SAGHEER S A. Bayesian prediction for progressively type-II censored data from the rayleigh model ［J］. Commun Statist Theory Meth，2005，1（34）：2353-2361.

［11］WU S J，CHEN D H，CHEN S T. Bayesian inference for rayleigh distribution under progressive censored sample［J］. Appl Stochast Models Bus Ind，2006，2（22）：269-279.

［12］ELRAHMAN A M，SULTAN K S. Reliability estimation based on general progressive censored data from the Weibull model：Comparison between Bayesian and classical approaches［J］. METRON - Int Statist LXV，2007，6 （16）：239-257.

［13］WU S J. Estimation of the two-parameter bathtub-shaped lifetime distribution with progressive censoring［J］. Appl Statist. 2008，3（35）：1139-1150.

［14］KUNDU D. Bayesian inference and life testing plan for the weibull distribution in presence of progressive censoring［J］. Technometrics，2008，9（50）：144-154.

［15］尹斌，鄭先斌.武器裝備全壽命管理中的標(biāo)準(zhǔn)化體系研究［J］.四川兵工學(xué)報(bào)，2013，24（10）：63-66.

［16］NELSON W. Applied life data analysis［M］. NewYork：Wiley Press，1982，3（176）：514-532.

Bayesian Model for Predicting Residual Life of Equipment

LIU Gang1，2，LI Fang1，2
（1. Naval University. of Engineering，Wuhan 430033，China；2. Lab of Naval Vessel Integrated Logistics Support，Naval University. of Engineering，Wuhan 430033，China）

Abstract：To solve the problem of weapon equipment is the growing complexity and the maintenance work becoming heavy daily，the Bayesian model is applied to predicting residual life of equipment after repair，the foundation was proposed for reasonably arranging repair plan. The traditional exponential distribution is replaced，then Weibull distribution is used for describing life characteristic of system，furthermore，and the methods of Maximum likelihood and Bayesian are applied to estimating the two parameters of Weibull distribution，and the credible intervals are given. On the basis，according to the two situations that prior sample and posterior sample，the residual life after repair is predicted by Bayesian model respectively. The illustrative example is given in paper，and the proposed method is proved to be effective.

Key words：weibull distribution，bayesian，residual life，predict

中圖分類號(hào)：TP391.73

文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1002-0640（2016）05-0019-06

收稿日期：2015-04-05修回日期：2015-05-07

作者簡(jiǎn)介：劉剛（1982-），男，湖北武漢人，碩士，講師。研究方向：裝備綜合保障工程、可靠性維修性測(cè)試性。