基于貝葉斯法的定時(shí)截尾小樣本指數(shù)型裝備器材需求預(yù)測(cè)

王鐵寧， 吳龍濤， 楊 帆

(裝甲兵工程學(xué)院技術(shù)保障工程系， 北京 100072)

基于貝葉斯法的定時(shí)截尾小樣本指數(shù)型裝備器材需求預(yù)測(cè)

王鐵寧， 吳龍濤， 楊 帆

(裝甲兵工程學(xué)院技術(shù)保障工程系， 北京100072)

針對(duì)高新裝備器材故障數(shù)據(jù)少、需求規(guī)律不明確的問題，提出了一種定時(shí)截尾小樣本條件下指數(shù)型裝備器材的需求預(yù)測(cè)方法?；谪惾~斯法進(jìn)行了裝備器材壽命分布參數(shù)估計(jì)，討論了先驗(yàn)分布和損失函數(shù)的選擇問題。引入了K-S檢驗(yàn)法對(duì)壽命分布模型進(jìn)行擬合優(yōu)度檢驗(yàn)，并設(shè)計(jì)了壽命分布參數(shù)估計(jì)的卡方檢驗(yàn)方法。綜合考慮故障更換和定時(shí)更換，提出了部隊(duì)裝備器材年度需求預(yù)測(cè)方法，并采用蒙特卡洛法進(jìn)行了仿真驗(yàn)證。結(jié)果表明:基于貝葉斯法的參數(shù)估計(jì)結(jié)果能夠順利通過檢驗(yàn)，需求預(yù)測(cè)結(jié)果與仿真結(jié)果一致。

定時(shí)截尾小樣本; 裝備器材； 指數(shù)分布； 需求預(yù)測(cè)； 貝葉斯法

隨著部隊(duì)高新裝備列裝步伐的加快和換件修理方式的廣泛開展，裝備器材需求結(jié)構(gòu)逐漸發(fā)生變化，如對(duì)光電類及總成類裝備器材的需求量顯著增加。由于高新裝備列裝時(shí)間短，當(dāng)前高新裝備器材保障普遍存在消耗規(guī)律掌握不清、保障不及時(shí)的難題。

目前，我軍主要采用定額計(jì)算[1]方法管理庫(kù)存，各級(jí)倉(cāng)庫(kù)一般依據(jù)庫(kù)存標(biāo)準(zhǔn)，并結(jié)合現(xiàn)有庫(kù)存來制定需求計(jì)劃。但由于庫(kù)存標(biāo)準(zhǔn)本身就是粗略、不精確的估計(jì)[2]，因此定額計(jì)算方法難以精確測(cè)算裝備器材需求。同時(shí)，由于高新裝備服役時(shí)間短，可利用的歷史數(shù)據(jù)較少，故難以應(yīng)用傳統(tǒng)的基于大樣本的方法，如時(shí)間序列分析、回歸分析和神經(jīng)網(wǎng)絡(luò)方法等進(jìn)行器材需求預(yù)測(cè)?；疑A(yù)測(cè)模型[3]是小樣本預(yù)測(cè)經(jīng)常采用的方法，但其預(yù)測(cè)精度較低，且無法計(jì)算裝備器材滿足率。

基于可靠性的需求預(yù)測(cè)方法根據(jù)失效機(jī)理選擇裝備器材壽命分布類型[4]，并以此建立裝備器材需求計(jì)算模型。然而，由于裝備承制單位在交付裝備時(shí)往往無法提供其壽命分布的具體形式，且不同部隊(duì)的裝備使用和維護(hù)水平以及使用環(huán)境差異明顯，裝備承制單位提供的數(shù)據(jù)只能作為參考[5]。因此，筆者針對(duì)定時(shí)截尾小樣本，提出了基于貝葉斯法[6]的指數(shù)型裝備器材的需求預(yù)測(cè)方法，利用貝葉斯法估計(jì)裝備器材的壽命分布參數(shù)，然后綜合考慮故障更換和定時(shí)更換，給出了裝備器材的需求預(yù)測(cè)方法，并進(jìn)行了仿真驗(yàn)證。

1 指數(shù)型壽命分布參數(shù)估計(jì)

1.1指數(shù)分布的一般形式

指數(shù)分布的重要特性之一是故障率λ(t)=λ為常數(shù)，即產(chǎn)品在下一時(shí)刻故障失效的概率與已使用時(shí)間無關(guān)，這也被稱作指數(shù)分布的“無后效性”[7]。設(shè)某裝備器材壽命T服從指數(shù)分布，其概率密度函數(shù)、累積分布函數(shù)和產(chǎn)品的可靠度函數(shù)分別為

f(t)=λexp(-λt),(t,λ>0)，

(1)

F(t)=1-exp(-λt),

(2)

R(t)=exp(-λt)。

(3)

E(T)=1/λ，為指數(shù)分布的期望，即產(chǎn)品的期望壽命(Mean Time To Failure，MTTF)。由指數(shù)分布的一般形式可知：對(duì)指數(shù)型裝備器材進(jìn)行可靠性估計(jì)就是估計(jì)故障率λ。極大似然估計(jì)法是一種在大樣本條件下簡(jiǎn)單有效的參數(shù)估計(jì)方法，但在截尾小樣本情形下偏差較大，且不易收斂[8]。貝葉斯法可充分利用專家經(jīng)驗(yàn)、裝備承制單位的參考數(shù)據(jù)等先驗(yàn)信息來提高參數(shù)估計(jì)的準(zhǔn)確性，同時(shí)隨著新信息的出現(xiàn)可不斷更新結(jié)果，是解決小樣本數(shù)據(jù)預(yù)測(cè)的有效方法[9]。

1.2故障率λ先驗(yàn)分布的選擇

應(yīng)用貝葉斯法進(jìn)行壽命分布參數(shù)估計(jì)時(shí)，首先要為參數(shù)選擇合適的先驗(yàn)分布。通?？蛇x擇無信息先驗(yàn)分布和有信息的共軛先驗(yàn)分布。

1)當(dāng)關(guān)于待估參數(shù)的先驗(yàn)信息很少時(shí)，一般選擇無信息先驗(yàn)分布。根據(jù)最常用的確定無信息先驗(yàn)分布的Jeffreys法則[9]，λ的先驗(yàn)分布為

π(λ)=λ-1。

(4)

2)共軛先驗(yàn)分布是指與后驗(yàn)分布屬于同一類分布的先驗(yàn)分布。伽馬分布為指數(shù)分布故障率λ常用的共軛先驗(yàn)分布，即λ～Γ(α,β)，其概率密度函數(shù)為

(5)

式中：α,β>0，為先驗(yàn)分布的超參數(shù)，其值通常根據(jù)先驗(yàn)信息使用矩估計(jì)法來確定；Γ(α)為伽馬函數(shù)。

1.3故障率λ的貝葉斯估計(jì)

在裝備維修過程中，裝備器材的更換原因大致分為2類：因故障失效被更換和到了指定的維修間隔期進(jìn)行等級(jí)維修時(shí)被更換。因此，可將裝備器材的更換時(shí)間樣本視為定時(shí)截尾樣本，截尾時(shí)間為定時(shí)維修間隔期。

設(shè)裝備的定時(shí)維修間隔期為t*，裝備器材更換時(shí)間樣本t的容量為n，故障失效時(shí)間數(shù)量為z，第i(1≤i≤z)個(gè)故障失效時(shí)間為ti，則樣本的聯(lián)合分布密度函數(shù)為

(6)

1)若選擇無信息先驗(yàn)，則樣本t和λ的聯(lián)合概率分布為

(7)

(8)

根據(jù)貝葉斯定理，在給定樣本t下，λ的后驗(yàn)分布密度函數(shù)為

(9)

即

(10)

2)若選擇共軛伽馬先驗(yàn)，λ的后驗(yàn)分布密度函數(shù)為

(11)

(12)

(13)

2 壽命分布參數(shù)估計(jì)檢驗(yàn)

2.1壽命分布模型的擬合優(yōu)度檢驗(yàn)

通過貝葉斯估計(jì)確定了裝備器材壽命分布的形式后，就需要通過壽命分布的擬合優(yōu)度檢驗(yàn)來驗(yàn)證分布模型是否適用于樣本。χ2檢驗(yàn)是一種常用的擬合優(yōu)度檢驗(yàn)方法，但只適用于大樣本(n≥50)完全數(shù)據(jù)的情形，針對(duì)小樣本定時(shí)截尾數(shù)據(jù)，筆者引入了K-S檢驗(yàn)法[5]。

設(shè)裝備器材壽命T的真實(shí)分布函數(shù)為Ψ(t)，則檢驗(yàn)假設(shè)為

H0:Ψ(t)=F(t)；H1:Ψ(t)≠F(t)。

K-S檢驗(yàn)的統(tǒng)計(jì)量為

(14)

式中：Ψz(t)為壽命T關(guān)于給定失效樣本的經(jīng)驗(yàn)分布函數(shù)。由于Ψz(t)是小樣本，其計(jì)算公式[11]為

(15)

2.2參數(shù)估計(jì)的假設(shè)檢驗(yàn)

若默認(rèn)裝備器材壽命T服從指數(shù)分布，則只需要對(duì)λ的估計(jì)結(jié)果進(jìn)行假設(shè)檢驗(yàn)。且檢驗(yàn)假設(shè)為

3 指數(shù)型裝備器材需求預(yù)測(cè)

圖1為部隊(duì)裝備器材需求發(fā)生的一般過程。

圖1 部隊(duì)裝備器材需求發(fā)生一般過程

假設(shè)某部隊(duì)有A型現(xiàn)役裝備m臺(tái)，每臺(tái)裝備有1件器材a。根據(jù)滾動(dòng)式循環(huán)動(dòng)用原則，一般每臺(tái)裝備的年度動(dòng)用時(shí)間不同，設(shè)Ti(1≤i≤m)為第i臺(tái)裝備的年度動(dòng)用時(shí)間。為了便于分析計(jì)算，進(jìn)行如下假設(shè)：

1)裝備器材故障失效時(shí)采取換件維修，換件時(shí)間忽略不計(jì)；

2)裝備的動(dòng)用時(shí)間即為裝備器材的動(dòng)用時(shí)間；

3)各裝備中器材a的壽命相互獨(dú)立且服從同一指數(shù)分布。

在計(jì)算裝備器材年度需求時(shí)，需要分別考慮故障更換需求和定時(shí)更換需求。已知壽命服從指數(shù)分布的器材a在年度動(dòng)用時(shí)間TS內(nèi)的故障次數(shù)N(Ti)服從泊松分布P(λTi)，即

(16)

式中：k∈N。根據(jù)泊松分布的可加性可得m件器材a的年度故障次數(shù)為N(TS)=N(T1)+N(T2)+…+N(Tm)，服從泊松分布N(TS)～P(λT1+λT2+…+λTm)，即

(17)

若規(guī)定器材a的年度裝備器材滿足率為ρ，則因故障更換產(chǎn)生的裝備器材需求量S1滿足

(18)

根據(jù)裝備維修計(jì)劃，一般每年都會(huì)安排摩托小時(shí)已到達(dá)定時(shí)維修間隔期的裝備進(jìn)行等級(jí)維修。在等級(jí)維修時(shí)，器材a可能采取不更換、必須更換或視情更換3種維修策略[2]，對(duì)于視情更換，通常指定一個(gè)更換概率p。假設(shè)在m臺(tái)裝備中有m0臺(tái)需要進(jìn)行等級(jí)維修，如器材a不更換，其定時(shí)更換需求量S2=0；如必須更換，S2=m0；如視情更換，則S2服從二項(xiàng)分布S2～B(m0,p)，同式(18)，滿足率為ρ時(shí)的S2為

(19)

裝備器材a的年度需求量為S=S1+S2。

4 算例分析

以某型裝甲車輛的“火控計(jì)算機(jī)通信板”(簡(jiǎn)稱“通信板”)為例，利用本文提出的方法進(jìn)行可靠性統(tǒng)計(jì)分析和需求預(yù)測(cè)。已知單臺(tái)裝甲車輛安裝一塊通信板，通過裝備承制單位獲知通信板的工作壽命“大約為1 200 h”。某部隊(duì)在裝備動(dòng)用過程中得到一組更換時(shí)間樣本為81 h，230 h，274 h，387 h，600 h，688 h，798 h，974 h，1 000*h，1 000*h(帶*的為截尾數(shù)據(jù))。

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

首先，計(jì)算通信板故障率λ的后驗(yàn)分布。若選擇無信息Jeffreys先驗(yàn)，由式(10)可得λ的后驗(yàn)分布為Γ(8,6 032)。若選擇伽馬先驗(yàn)，以裝備承制單位提供的數(shù)據(jù)作為先驗(yàn)信息，設(shè)通信板壽命均值為1 200 h，標(biāo)準(zhǔn)差也設(shè)為1 200 h，當(dāng)λ～Γ(α,β)時(shí)，MTTF服從IΓ(α,β)分布。應(yīng)用矩估計(jì)法可得

(20)

由式(20)可得先驗(yàn)分布的超參數(shù)α=3，β=2 400。再由式(11)可得λ的后驗(yàn)分布為Γ(11,8 432)。圖2 為λ先驗(yàn)分布和后驗(yàn)分布的概率密度函數(shù)對(duì)比分析?？梢钥闯觯?種后驗(yàn)分布得到的故障率都高于裝備承制單位給出的故障率；同時(shí)，由Jeffreys先驗(yàn)得到的后驗(yàn)故障率高于由伽馬先驗(yàn)得到的后驗(yàn)故障率，這是因?yàn)镴effreys先驗(yàn)為無信息先驗(yàn)，其后驗(yàn)分布取決于樣本。一方面，由于裝備承制單位無法預(yù)知部隊(duì)的裝備使用環(huán)境及裝備使用、維修人員的水平，因而高估了通信板的實(shí)際使用壽命；另一方面，該樣本僅為1次小樣本抽樣，存在一定的偶然性和不確定性。由伽馬先驗(yàn)得到的后驗(yàn)分布綜合考慮了裝備承制單位先驗(yàn)信息和樣本的新信息，因此更加合理。

圖2 λ先驗(yàn)分布和后驗(yàn)分布概率密度對(duì)比

表1 K-S模型擬合優(yōu)度檢驗(yàn)結(jié)果

4.2需求預(yù)測(cè)

假定該單位有該型裝甲車輛100臺(tái)，年度動(dòng)用計(jì)劃為(40,15)、(30,25)、(20,60)、(10,150)，并計(jì)劃安排10臺(tái)車輛進(jìn)行等級(jí)維修，等級(jí)維修時(shí)通信板采取視情維修策略，更換概率為0.2，需預(yù)測(cè)滿足率為0.95時(shí)通信板的需求量。

為了驗(yàn)證預(yù)測(cè)結(jié)果的準(zhǔn)確性，采用蒙特卡洛法對(duì)裝甲車輛通信板需求預(yù)測(cè)過程進(jìn)行仿真。圖3為基于MATLAB平臺(tái)的裝備器材需求預(yù)測(cè)仿真流程。

圖3 裝備器材需求預(yù)測(cè)仿真流程

圖4 10 000次仿真的年度裝甲車輛通信板需求預(yù)測(cè)樣本

圖5 仿真樣本的經(jīng)驗(yàn)分布函數(shù)

5 結(jié)論

筆者以指數(shù)壽命型高新裝備器材為對(duì)象，結(jié)合實(shí)際情況，對(duì)小樣本條件下裝備器材的需求預(yù)測(cè)問題進(jìn)行了研究。使用貝葉斯法綜合了廠家先驗(yàn)信息和使用過程中產(chǎn)生的樣本信息，對(duì)裝備器材的壽命分布參數(shù)進(jìn)行了估計(jì)，并給出了分布模型和分布參數(shù)的檢驗(yàn)方法。在此基礎(chǔ)上，綜合考慮故障更換和定時(shí)更換，提出了裝備器材年度需求預(yù)測(cè)方法，并采用蒙特卡洛法進(jìn)行了仿真驗(yàn)證，為解決當(dāng)前高新裝備器材消耗規(guī)律掌握不清、保障不及時(shí)的問題提供了一種可行的方法。

[1] 劉旭陽，吳龍濤，周萬里.基于ARIMA模型的裝備器材需求預(yù)測(cè)方法[J].裝甲兵工程學(xué)院學(xué)報(bào)，2016，30(6)：21-25.

[2] 劉慎洋，高崎，李志偉，等.基于使用壽命的裝備備件消耗預(yù)測(cè)模型[J].裝甲兵工程學(xué)院學(xué)報(bào)，2015，29(4)：27-30.

[3] 郭曉君，劉思峰，吳利豐. 污染物減排預(yù)測(cè)的灰色Markov組合模型與算法[J]. 計(jì)算機(jī)應(yīng)用研究，2013，20(12)：3670-3673.

[4] 李源，楊建軍.導(dǎo)彈武器裝備備件需求量分析與計(jì)算方法[J].兵工自動(dòng)化，2007，26(8)：22-23.

[5] 蔡景，左洪福，王華偉.基于成本的民用航空發(fā)動(dòng)機(jī)維修方案優(yōu)化研究[J].機(jī)械科學(xué)與技術(shù)研究，2007，26(2)：167-171.

[6] 曾志國(guó).貝葉斯可靠性[M].北京：國(guó)防工業(yè)出版社，2014：67-70.

[7] 張春曉，周美茵，王志平.基于貝葉斯更新的飛機(jī)結(jié)構(gòu)腐蝕可靠度模型[J].航空學(xué)報(bào)，2014，35(7)：1931-1940.

[8] 馮自立，陳曉陽，顧家銘.Weibull分布下基于矩法的小樣本量定時(shí)截尾數(shù)據(jù)的估計(jì)[J].四川兵工學(xué)報(bào)，2010，31(3)：93-97.

[9] SARHAN A.M. Empirical bayes estimates in exponential reliability model [J]. Applied mathematics and computation 2003，135(1)：319-332.

[10] 尹慧琳，楊筱菡，陸恒.Wiener過程性能退化產(chǎn)品可靠性評(píng)估新Bayes方法[J].同濟(jì)大學(xué)學(xué)報(bào)(自然科學(xué)版)，2015,43(8)：1234-1237.

[11] EBELING C.E.An introduction to reliability and maintainability engineering [M].New York：Waveland Pr Inc,2009:52-53.

[12] 劉天華，張志華，李慶民，等. 威布爾型多不可修部件備件需求確定方法[J].系統(tǒng)工程理論與實(shí)踐，2012，32(9)：2010-2015.

(責(zé)任編輯: 王生鳳)

DemandForecastingforEquipmentMaterialswithExponentialLifeDistributionBasedonBayesianEstimationunderTypeICensoredSmallSample

WANG Tie-ning, WU Long-tao, YANG Fan

(Department of Technical Support Engineering, Academy of Armored Force Engineering, Beijing100072, China)

In order to address the problem that the demand discipline of new and high-tech equipment materials cannot be mastered well because of lack of failure data, a method of demand forecasting is proposed for equipment materials with exponential distribution under small failure samples. In a view of small samples, the parameter of equipment material life distribution is estimated by Bayesian estimation, the option of prior distribution and loss function is discussed. Then, statistical tests for goodness-of-fit of the life distribution and parameter estimation are performed by K-S test and Chi-Square test respectively. In consideration of fault replacement and timely replacement, the annual demand forecasting method of army equipment is put forward and devised with a Monte-Carlo simulation test. As the example shows, the estimation result performs well in the test, and the forecasting value is identical to that from the simulation.

type I censored small sample; equipment materials; exponential distribution; demand forecasting; Bayesian estimation

1672-1497(2017)04-0029-06

2017-03-28

軍隊(duì)科研計(jì)劃項(xiàng)目

王鐵寧(1962-)，男，教授，博士。

E92;N945.24

：ADOI：10.3969/j.issn.1672-1497.2017.04.006