基于貝葉斯的航空備件需求研究與應(yīng)用

陳鳳騰，左洪福

（1.徐州工程學(xué)院，江蘇 徐州 221008；2.南京航空航天大學(xué)民航學(xué)院，南京 210016）

基于貝葉斯的航空備件需求研究與應(yīng)用

陳鳳騰1，左洪福2

（1.徐州工程學(xué)院，江蘇 徐州 221008；2.南京航空航天大學(xué)民航學(xué)院，南京 210016）

鑒于目前航空備件需求計(jì)算應(yīng)用泊松公式存在計(jì)算結(jié)果誤差較大的問(wèn)題，依據(jù)歷史故障數(shù)據(jù)并應(yīng)用廣義更新過(guò)程建立反映歷史故障數(shù)據(jù)以及維修能力的可靠性分布模型，同時(shí)針對(duì)目前歷史故障數(shù)據(jù)主要是小樣本事件，應(yīng)用貝葉斯方法并采用切片采樣法對(duì)可靠性分布模型進(jìn)行參數(shù)估計(jì)。最后通過(guò)實(shí)例進(jìn)行分析，結(jié)果表明實(shí)際故障次數(shù)和實(shí)際故障情況較為吻合，增大了模型的應(yīng)用范圍。

貝葉斯；航空備件；故障率；切片采樣法

備件需求預(yù)測(cè)直接影響航空公司安全運(yùn)營(yíng)和成本管理，是維修活動(dòng)的基礎(chǔ)，是關(guān)系到航空公司信譽(yù)和效益的一個(gè)關(guān)鍵環(huán)節(jié)。航空備件需求和預(yù)測(cè)問(wèn)題涉及了可靠性、維修性、航空保障以及系統(tǒng)安全等領(lǐng)域，成為國(guó)內(nèi)外學(xué)者關(guān)注的課題之一[1-2]。

目前飛機(jī)制造商（如波音、空客等）向用戶提供備件推薦數(shù)量的計(jì)算方法主要是齊次泊松公式[3]。該方法主要假設(shè)各種零部件的可靠性模型服從理想的指數(shù)分布，這在計(jì)算初始備件需求數(shù)量時(shí)具有方便、易于計(jì)算等特點(diǎn)，但未能充分反映實(shí)際設(shè)備故障特性，而且應(yīng)用該方法要假設(shè)部件的修理和更換服從更新過(guò)程，這和實(shí)際情況脫節(jié)，其預(yù)測(cè)結(jié)果和實(shí)際需求偏差較大。

面對(duì)這些問(wèn)題，文獻(xiàn)[4]引進(jìn)維修度并應(yīng)用廣義更新過(guò)程進(jìn)行可靠性建模，不僅可以反映設(shè)備實(shí)際更換的特點(diǎn)，而且還可以反映設(shè)備維修能力以及修理后恢復(fù)的能力。

對(duì)于歷史故障、壽命數(shù)據(jù)為大樣本時(shí)，應(yīng)用上述方法并通過(guò)似然函數(shù)進(jìn)行參數(shù)估計(jì)，可以獲得較高的精度。但由于目前使用年限以及管理能力、水平的限制，獲得的歷史故障數(shù)據(jù)往往是小樣本，直接應(yīng)用似然函數(shù)計(jì)算方法存在其計(jì)算精度較低的問(wèn)題，為此，本文根據(jù)航空裝備6種故障率特性，首先建立相應(yīng)的故障分布模型，然后根據(jù)歷史故障數(shù)據(jù)和先驗(yàn)知識(shí)，初步確定模型的分布和范圍。由于使用年限和管理水平的限制，設(shè)備故障多為小樣本事件，故本文選取適于處理小樣本事件的貝葉斯方法，同時(shí)對(duì)于高維復(fù)雜的無(wú)法直接計(jì)算的數(shù)學(xué)公式應(yīng)用切片采樣法進(jìn)行數(shù)值仿真，并將其應(yīng)用到任意分布函數(shù)。為了反映同類(lèi)設(shè)備歷史故障和航空公司對(duì)該類(lèi)設(shè)備的實(shí)際維修能力，本文借助廣義更新過(guò)程進(jìn)行可靠性建模。

1 廣義更新過(guò)程

廣義更新過(guò)程（GRP）首先被Kijima和Sumita提出并應(yīng)用到復(fù)雜可修系統(tǒng)的分析[4-5]，根據(jù)維修方式的不同將該系統(tǒng)分5種狀態(tài)，同時(shí)引入表征維修能力的維修度（q）概念和其對(duì)應(yīng)：修復(fù)如新（q=0）、修復(fù)如舊（q=1）、修復(fù)介于新舊之間（01）。通常航空備件可分為可修件和不可修件兩種，對(duì)于可修件，可以直接進(jìn)行上述分類(lèi)；至于不可修件，考慮到出現(xiàn)故障后直接用新件更換，故也將其視為可修件，其特性符合“修復(fù)如新”，因而，可將廣義更新過(guò)程應(yīng)用于所有的航空備件需求量預(yù)測(cè)分析和應(yīng)用。

對(duì)于航空備件需求模型，對(duì)應(yīng)不同維修策略，根據(jù)該類(lèi)設(shè)備的故障特性建立起相應(yīng)的故障分布模型，限于篇幅，本文主要對(duì)可靠性分布模型為Weibull進(jìn)行計(jì)算說(shuō)明。

2 極大似然法

當(dāng)歷史故障數(shù)據(jù)為大樣本時(shí)，往往用極大似然法對(duì)分布模型進(jìn)行參數(shù)估計(jì)，由Weibull分布函數(shù)和公式可求得第i（1≤i≤n）次的概率密度函數(shù)如式（1）和式（2）所示

在應(yīng)用廣義更新過(guò)程進(jìn)行參數(shù)估計(jì)時(shí)存在以下兩種情況：故障截尾和時(shí)間截尾。故障截尾表示既定故障次數(shù)出現(xiàn)之后即可結(jié)束故障事件統(tǒng)計(jì)，而時(shí)間截尾表示既定時(shí)間T到達(dá)之后即可結(jié)束統(tǒng)計(jì)。對(duì)于對(duì)這兩種截尾進(jìn)行參數(shù)估計(jì)主要通過(guò)如下兩種似然函數(shù)計(jì)算[6-7]。

1）故障截尾的似然函數(shù)為

為了求q的參數(shù)估計(jì)值，現(xiàn)對(duì)ln L求參數(shù)q的偏導(dǎo)數(shù)，并令其為0，聯(lián)立上面3個(gè)公式，應(yīng)用擬牛頓迭代算法可得到參數(shù)估計(jì)值。

2）時(shí)間截尾的似然函數(shù)為

式中：δi表示截尾與否的變量，對(duì)q的參數(shù)估計(jì)方法同上。

當(dāng)歷史數(shù)據(jù)較多時(shí)，應(yīng)用極大似然函數(shù)效果較好，對(duì)于小樣本數(shù)據(jù)來(lái)說(shuō)，其計(jì)算精度較低，而采用Bayes法基于先驗(yàn)知識(shí)，處理小樣本精度較高。

3 Bayes法

根據(jù)圖1首次壽命周期概率密度分布如式（5）所示

其中：f（1θT1）是 θ的后驗(yàn)分布在周期1給定故障經(jīng)歷；L（T1θ）是似然函數(shù)。

式（7）求得的后驗(yàn)分布在第2次循環(huán)周期內(nèi)則變成先驗(yàn)分布，具體如式（6）所示

第n次循環(huán)的壽命周期內(nèi)應(yīng)用Bayes定理如式（7）所示

或展開(kāi)如式（8）所示

其中：n是故障次數(shù)；k是計(jì)算因子。具體如式（9）所示

4 切片樣本法

直接計(jì)算式（8）比較困難，可以采用切片取樣法[8]對(duì)其進(jìn)行求解，切片取樣法主要是根據(jù)多變量函數(shù)的超矩陣進(jìn)行切片取樣，該方法主要根據(jù)Markov-Chain Monte-Carlo（MCMC）樣本技術(shù)，假定從某分布中對(duì)變量x取樣，在Rn的某個(gè)子集中取值，其密度與函數(shù)f（x）成比例，在函數(shù)f（x）的n+1維區(qū)域內(nèi)進(jìn)行均勻取樣，通過(guò)輔助實(shí)值變量y，定義在區(qū)域U={（x，y）：0

其中：x的邊際密度函數(shù)如式（11）所示

為了對(duì)變量x取樣，可以在區(qū)域（x，y）內(nèi)進(jìn)行聯(lián)合取樣，然后忽略變量y。定義一個(gè)收斂于Markov-Chain的均勻分布并從U中均勻產(chǎn)生獨(dú)立點(diǎn)。通過(guò)Gibbs樣本，從給定x值對(duì)y的條件分布進(jìn)行取樣，即讓其在區(qū)間（0，f（x））內(nèi)是均勻的；然后從給定y值再對(duì)x的條件分布進(jìn)行取樣，也即在區(qū)域S={x：y

其中：π（θ，{t1，…，tn}）的先驗(yàn)密度函數(shù)π0（α，β，q）各個(gè)參數(shù)相互獨(dú)立，具體表示如式（13）所示

其中：參數(shù)α和q視作均勻分布，故其密度函數(shù)為1；參數(shù)β服從正態(tài)分布，由于0≤β≤10，故可取μβ=0.63和 σβ=0.5。

應(yīng)用多變量分布使用切片取樣法，對(duì)于多參變量，定義區(qū)間 I=（L，R）通過(guò)超矩陣 H={x：Li

1）在（0，f（x0））上均勻取一個(gè)實(shí)值y，然后定義一個(gè)水平切片：S={x：y

2）尋求 x0附近超矩陣，H=（L1，R1）× …×（Ln，Rn），包含全部或部分的切片；

3）在超矩陣從切片中繪制新點(diǎn)x1。

第2步的超區(qū)間值，Neal通過(guò)4種方法尋找這個(gè)區(qū)間。對(duì)于Weibull分布，考慮所獲得先驗(yàn)知識(shí)對(duì)參數(shù)進(jìn)行限制，分別取 tmin≤α≤tmax、0≤β≤10 以及-2≤q≤5。

5 故障次數(shù)計(jì)算

由切片取樣法產(chǎn)生參數(shù)α、β和q后，代入到式（1）和式（2），即建立反映壽命、維修能力的故障壽命分布概率函數(shù)，然后再采用文獻(xiàn)[4]的計(jì)算方法，也即應(yīng)用式（14）和式（15）計(jì)算故障時(shí)間，進(jìn)而求得既定時(shí)間內(nèi)的故障次數(shù)[4]。

計(jì)算過(guò)程如下：

1）如果Σti≥T，則停止產(chǎn)生隨機(jī)數(shù)；

2）如果Σti＜ T，按照式（14）和式（15）繼續(xù)產(chǎn)生下一個(gè)隨機(jī)數(shù)。

6 實(shí)例分析

本文取5 000 h內(nèi)應(yīng)用在5架飛機(jī)上的某個(gè)零件故障數(shù)據(jù)記錄，共有25個(gè)數(shù)據(jù)，具體如表1所示，其中失效時(shí)間是設(shè)備發(fā)生故障的時(shí)間，這不僅包括首次故障時(shí)間，也包括發(fā)生故障進(jìn)行修理后重新使用直至下一次發(fā)生故障的時(shí)間。另外，表1中的δi數(shù)值1表示故障發(fā)生沒(méi)有發(fā)生截尾，0表示截尾。對(duì)于首先應(yīng)用數(shù)理統(tǒng)計(jì)對(duì)上述歷史故障數(shù)據(jù)進(jìn)行分析，由分析知該類(lèi)航空設(shè)備的壽命服從二參數(shù)Weibull分布（α，β），然后進(jìn)一步考察該類(lèi)設(shè)備受到維修能力（即維修度）對(duì)產(chǎn)品壽命的影響。

對(duì)歷史故障取樣本容量N=25，在觀察周期內(nèi)記錄個(gè)體失效的時(shí)間，記為t（i），對(duì)上述數(shù)據(jù)進(jìn)行初步數(shù)理統(tǒng)計(jì)分析后，初步確立分布參數(shù)的范圍為：－3≤q≤4，0≤β≤10，min（ti）≤α≤max（ti），應(yīng)用式（12）先驗(yàn)分布作為三參數(shù)模擬的f（x），應(yīng)用多變量切片取樣法，分別對(duì)參數(shù)α給兩個(gè)初值3 492和3 501，進(jìn)行多參數(shù)切片樣本法模擬參數(shù)，模擬過(guò)程如圖2所示，從圖2中可以清楚地看出參數(shù)α在迭代次數(shù)為540次時(shí)，兩條鏈?zhǔn)諗康? 499.27 h左右；類(lèi)似地，參數(shù)β賦予兩個(gè)初值3和4，其模擬過(guò)程如圖3所示，在迭代次數(shù)為610次時(shí)，兩條鏈?zhǔn)諗康?.75左右。維修度q賦予兩個(gè)初值0.1和3.2，模擬的參數(shù)如圖4所示，參數(shù)q在迭代次數(shù)為710次時(shí)，兩條鏈?zhǔn)諗康剑?.287 1左右。從這3個(gè)參數(shù)的迭代圖中可以清楚地看到每個(gè)參數(shù)賦予兩個(gè)初值，經(jīng)過(guò)若干次迭代后具有收斂的特性，故取α為3 499.27，β為2.75，q為0.287 1。然后將這3個(gè)參數(shù)代入式（1）和式（2），得到了反映故障、維修能力的壽命分布概率密度函數(shù)；通過(guò)該概率密度函數(shù)，應(yīng)用式（14）和式（15）分別仿真5架飛機(jī)的該零件在 1 000、2 000、3 000、4 000、5 000 和 6 000 h 內(nèi)的故障次數(shù)，具體如圖5所示，圖5中的實(shí)際故障次數(shù)即表1中的故障，其中貝葉斯計(jì)算方法表示用該方法仿真出的數(shù)據(jù)。

表1 零件故障數(shù)據(jù)的監(jiān)控記錄Tab.1 Monitoring records of fault data for part

通過(guò)圖5的比較，二者的數(shù)據(jù)基本吻合，可以用于小樣本故障數(shù)據(jù)的故障壽命分布概率密度函數(shù)建模，同時(shí)還能反映維修能力的高低。

7 結(jié)語(yǔ)

針對(duì)航空裝備歷史故障的小樣本事件，應(yīng)用廣義更新過(guò)程建立反映該設(shè)備歷史故障以及該設(shè)備發(fā)生故障后的維修能力的模型，對(duì)于分布函數(shù)比較復(fù)雜的可以通過(guò)切片樣本法完成，不僅可以處理常用的分布函數(shù)，還可處理任意分布函數(shù)，另外，本方法最大的特點(diǎn)是建立在已有經(jīng)驗(yàn)和先驗(yàn)知識(shí)的基礎(chǔ)上，可以充分反映實(shí)際運(yùn)行情況，在航空裝備可靠性建模、故障預(yù)測(cè)和備件需求量的計(jì)算方面使精度得到很大的提高。另外，隨著計(jì)算機(jī)技術(shù)的發(fā)展，使得應(yīng)用本方法不僅提高了精度，而且既簡(jiǎn)單又方便，極大地增加了模型的應(yīng)用范圍。

[1]DINESH KUMAR U，KNEZEVIC J.Availability based spare optimization using renewal process[J].Reliability Engineering and System Safety，1998，59：217-223.

[2] 伏宏勇，趙 宇.航空備件計(jì)劃模型與方法研究[D].北京：北京航空航天大學(xué)，2003.

[3] BOEING COMPANY.Spares Provisioning ProductsGuide[M].Seattle：Boeing Company，2003.

[4] 陳鳳騰，左洪福.基于廣義更新過(guò)程的航空備件需求研究和應(yīng)用[J].應(yīng)用科學(xué)學(xué)報(bào)，2007，10：47-49.

[5]YANEZ M，JOGLAR F，MODARRES M.Generalized renewal process for analysis of repairable systems with limited failure experience[J].Reliability Engineering and System safety，2002，77：167-180.

[6] 鄧永錄，梁之舜.隨機(jī)點(diǎn)過(guò)程及其應(yīng)用[M].北京：科學(xué)出版社，1992.

[7]EDWARD P C KAO.An Introduction to Stochastic Processes[M].Beijing：China Machine Press，2003.

[8] ANDREW G J.Generalisation and bayesian solution of the general renewal process for modeling the reliability effects of imperfect inspection and maintenance based on imprecise data[D].USA：University of Maryland，2005.

[9]RADFORD M NEAL.Slice Sampling[J].The Annals of Statistics.2003，31（3）：705-767.

Research and Application of Aero-Spare Based on Bayesian Method

CHEN Feng-teng1，ZUO Hong-fu2
（1.Xuzhou Institute of Technology， Xuzhou 221008，China；2.Civil Aviation College， Nanjing University of Aeronautics and Astronautics， Nanjing 210016，China）

There are much more deviations in the calculation of the aero-spares demand by the Poisson formula in recently decades，so we establish the reliability distribution model to reflect the historical failure data and the maintenance capabilities with the application of the generalized renewal process.At the same time for the historical fault data for the current principal is a small sample of events，Bayesian method and slice sampling method were applied to estimate the parameter of the reliability distribution.Finally，with the example analysis， the results show that the number of actual faults and actual failures of the more consistent，increasing the application of this model.

Bayesian； aero-spares； failure rate； slice sampling

TP3

A

1674－5590（2011）02－0013－05

2010-11-22；

2011-02-27 基金項(xiàng)目：國(guó)家自然科學(xué)基金項(xiàng)目（60672164）；徐州市科技計(jì)劃項(xiàng)目（XJ09070）

陳鳳騰（1971—），男，江蘇贛榆人，講師，博士，研究方向?yàn)榭煽啃?、維修性與故障診斷.

（責(zé)任編輯：楊媛媛）