基于MCMC 的MTBF 值區(qū)間估計方法研究

孟慶東，楊 珍

(佛山市南海區(qū)技師學院，廣東 528237)

0 引言

對機電產(chǎn)品進行可靠性分析最基本的工作是要確定其壽命分布參數(shù)。機電產(chǎn)品壽命服從兩參數(shù)Weibull 分布，其分布參數(shù)的貝葉斯估計需要對后驗分布進行二重積分。根據(jù)貝葉斯公式，在求解可靠性特征量MTBF 值區(qū)間估計時，需要進行內(nèi)積分限為函數(shù)的二重積分。在進行積分時，由于Bayes 后驗分布復雜且內(nèi)積分限為函數(shù)，數(shù)值積分方法失效［1-3］。

利用Bayes 理論進行統(tǒng)計推斷，需要積分運算。對于兩參數(shù)Weibull 分布，需要對Bayes 后驗分布進行二重積分。一般有兩種方法進行積分運算，即數(shù)值積分與Monte Carlo 積分。

運用數(shù)值積分計算時會遇到一些問題。當被積函數(shù)較復雜時，將導致計算精度降低。同時，積分區(qū)間的大小也會直接影響計算結果，使得Bayes 后驗計算出現(xiàn)不穩(wěn)定的趨勢。對于多重積分，因為每個外層積分都取決于內(nèi)層積分在一組點上的取值，所以會出現(xiàn)誤差的累積。

Monte Carlo 積分從目標分布中抽樣模擬，利用隨機抽樣所得點進行統(tǒng)計推斷，避免了數(shù)值積分中的求積節(jié)點計算與系統(tǒng)網(wǎng)格的劃分，使得積分計算的效率與可操作性顯著提高。但是，在Bayes 推斷中，運用Monte Carlo 方法是有困難的，具體表現(xiàn)在Bayes 后驗分布多為復雜、高維非標準的分布。當目標密度函數(shù)易于抽樣時，可以利用Monte Carlo 方法產(chǎn)生近似樣本，并利用樣本的函數(shù)值計算Bayes 后驗期望估計。當目標密度函數(shù)難以抽樣時，可以利用MCMC 方法思想直接構造出Bayes 后驗分布π(m，η|T)的馬爾科夫鏈，基于所得樣本可以進行各種統(tǒng)計推斷［4-6］。

1 MCMC 方法

MCMC 方法通過建立一個以后驗分布π(x)為平穩(wěn)分布的馬爾科夫鏈，對π(x)進行抽樣，然后根據(jù)所得樣本進行統(tǒng)計推斷［7］。如通過抽樣得到π(x)的樣本X(1)，…，X(n)，則:

其中D 為狀態(tài)空間，于是可得估計:

如果X(1)，…，X(n)是平穩(wěn)分布為π(x)的馬爾科夫過程的樣本時，上式也成立。

為了構造馬爾科夫鏈，設給定狀態(tài)Xt條件下，下一步狀態(tài)為Xt+1。則任意選擇一個不可約轉移概率q(* ，* )以及一個函數(shù)α(* ，* )，定義p(x，y)=q(x，y)α(x，y)，則p(x，y)為馬氏鏈一個轉移核。如果鏈在時刻t 處于狀態(tài)x，即X(t)= x，首先由q(* | x)產(chǎn)生一個潛在的轉移x →y，然后根據(jù)概率α(x，y)決定是否轉移。可以從Unif(0，1)中抽一個隨機數(shù)u，則:

通常，稱q(* ，* )為建議分布。

為使后驗分布π(x)成為平穩(wěn)分布，在有q(* ，* )后，應選擇一個α(* ，* )使相應的p(x，y)以π(x)為其馬氏鏈平穩(wěn)分布，可以選擇:

此時，可以證明產(chǎn)生的馬爾科夫鏈是可逆的，滿足平穩(wěn)方程，即: π(x)q(x，y)= π(y)q(y，x)

2 Weibull 分布下MTBF 值區(qū)間估計

2.1 先驗分布為伽馬－ 逆伽馬分布的MTBF 值區(qū)間

選取均勻分布分布為建議分布，則:

其中，α(x，y)為Metropolis-Hastings 比率，my和ηy是在一步迭代中，從建議分布中抽取的建議值，mx和ηx是一步迭代的初始值。迭代過程如下:

①θ = (m，η)= (θ(t－1)0，θ(t－1)1);

②從建議分布Unif(－0.5，0.5)中產(chǎn)生候選點;θ'

③令θ″ = θ + θ'，即讓建議值隨機游動一個距離;

④根據(jù)上式計算接受概率α(θ，θ″);

⑤以概率α(θ，θ″)接受θ″ = θ(t)，否則令θ =θ(t)。

以某機電新產(chǎn)品運行時間為現(xiàn)場樣本，用上述方法進行計算，得到結果如表1 所示。

表1 MCMC 方法計算結果

圖1 和圖2 分別是參數(shù)m、η 以及MTBF 值的全部迭代軌跡和后驗分布的密度直方圖。

圖1 參數(shù)m 和η 及MTBF 值迭代軌跡圖

圖2 參數(shù)m 和η 及MTBF 值的后驗密度直方圖

由于計算可靠性特征量MTBF 值區(qū)間估計時，數(shù)值積分方法失效，依據(jù)工程經(jīng)驗用參數(shù)點估計的區(qū)間估計端點來近似MTBF 值的置信區(qū)間。但是這種方法存在一定問題，一是求得的置信區(qū)間精度不高，有擴大的趨勢;二是近似方法并沒有可靠的理論支撐。通過表1 的計算結果，可以看出MCMC 方法解決了上述問題，使得可靠性評估在質(zhì)量與效率上都得到提高。

2.2 先驗分布為二元正態(tài)分布的MTBF 值區(qū)間估計

先驗分布為二元正態(tài)分布時，計算結果如表2 所示。

表2 MCMC 方法計算結果

圖3 和圖4 分別是參數(shù)m、η 以及MTBF 值的全部迭代軌跡和后驗分布的密度直方圖。

圖3 參數(shù)m 和η 及MTBF 值的迭代軌跡圖

圖4 參數(shù)m 和η 及MTBF 值的后驗密度直方圖

3 結論

通過MCMC 構建馬爾科夫鏈，使Bayes 后驗計算運用的數(shù)值積分方法轉化成從簡單的分布中抽樣并推斷。基于抽樣所得分布參數(shù)樣本，MTBF 值、可靠度、失效率等各種可靠性特征量的求解都可以有效實施，不再受數(shù)值積分的限制，提高了模型計算的穩(wěn)定性、可操作性與適應性。MCMC 方法解決了數(shù)值積分問題帶來的不便，有利于貝葉斯理論在可靠性評估中的推廣。

［1］林靜. 基于MCMC 的貝葉斯生存分析理論及其在可靠性評估中的應用［D］. 南京:南京理工大學，2008.

［2］李斌全，戴怡. 基于馬氏鏈蒙特卡洛方法的數(shù)控系統(tǒng)可靠性評估［J］. 組合機床與自動化加工技術，2011(10):69－71.

［3］潘海濤，溫小霓. 基于MCMC 方法的GARCH 模型參數(shù)估計［J］. 統(tǒng)計與信息論壇，2009，24(4):12－16.

［4］曹小群，宋君強，張衛(wèi)民，等. 基于MCMC 方法的Lorenz混沌系統(tǒng)的參數(shù)估計［J］. 國防科技大學學報，2010，32(2):68－72，145.

［5］李凌，徐偉. 基于MCMC 方法的繼電器加速壽命試驗分析［J］. 低壓電器，2010(2):13－16，35.

［6］王進玲，曾聲奎，馬紀明，等. 混合變量系統(tǒng)基于MCMC的自適應重要抽樣法［J］. 宇航學報，2012，33(1):94－101.

［7］黎光明，張敏強. 先驗信息對MCMC 方法估計概化理論方差分量變異量的影響［J］. 統(tǒng)計與決策，2012(7):27－30.