三參數(shù)威布爾分布參數(shù)估計(jì)及在可靠性分析中的應(yīng)用

第一作者鄭銳男，碩士，工程師，1985年4月生

三參數(shù)威布爾分布參數(shù)估計(jì)及在可靠性分析中的應(yīng)用

鄭銳

(中國(guó)船舶重工集團(tuán)公司第七一六研究所，江蘇連云港222006)

摘要：威布爾分布是可靠性中應(yīng)用最廣泛的分布之一。三參數(shù)威布爾分布尤其適用于在開始使用時(shí)有一時(shí)間段內(nèi)不發(fā)生故障的情況。由于該分布的位置參數(shù)不等于0，在參數(shù)估計(jì)時(shí)不能采用簡(jiǎn)單的參數(shù)估計(jì)方法實(shí)現(xiàn)，限制了該分布形式在可靠性分析中的應(yīng)用。根據(jù)三參數(shù)威布爾分布的特點(diǎn)提出了一種綜合圖解法和遺傳算法的參數(shù)估計(jì)方法，應(yīng)用該方法可以獲得更精確的參數(shù)估計(jì)值。隨后應(yīng)用于某系列數(shù)控車床計(jì)算機(jī)數(shù)控系統(tǒng)的故障分析中，驗(yàn)證了參數(shù)估計(jì)方法的可行性。

關(guān)鍵詞：三參數(shù)威布爾分布；參數(shù)估計(jì)；圖解法；遺傳算法；數(shù)控機(jī)床

基金項(xiàng)目：國(guó)家自然科學(xué)

收稿日期：2013-09-03修改稿收到日期：2014-03-27

中圖分類號(hào):TG659文獻(xiàn)標(biāo)志碼：A

Parameter estimation of three-parameter Weibull distribution and its application in reliability analysis

ZHENGRui(The 716th Research Institute of China Shipbuildnig Industry Corporation, Lianyungang 222006,China)

Abstract:Weibull distribution is one of the most widely used distributions in reliability analysis. However, the application of three-parameter Weibull distribution is often limited due to the fact that its positional parameter is not equal to 0 and its parameters cannot be estimated with simple methods. Actually, the three-parameter Weibull distribution is particularly suitable for the case that there is no fault in a period of time at the beginning of use. Here, a parameter estimation method synthesizing graphical method and genetic algorithm for three-parameter Weibull distribution was presented according to the characteristics of the distribution. Then the three-parameter Weibull distribution with the proposed parameter estimation method was applied in failure analysis of a series of CNC lathes and successful results were achieved.

Key words:three-parameter Weibull distribution; parameter estimation,graphical method; genetic algorithm; CNC machine tools

可靠性分析的首要問(wèn)題是尋找能夠確切反映系統(tǒng)故障機(jī)理并與故障數(shù)據(jù)的分析結(jié)果相符合的故障分布規(guī)律。常用的方法是將故障數(shù)據(jù)擬合成某種分布形式，在確定出各分布參數(shù)之后，對(duì)故障數(shù)據(jù)進(jìn)行可靠性評(píng)估和預(yù)測(cè)[1]。

威布爾分布是描述機(jī)械系統(tǒng)及其零部件壽命數(shù)據(jù)分布規(guī)律最常用的一種分布形式[2]。兩參數(shù)威布爾分布憑借其參數(shù)估計(jì)簡(jiǎn)單和適應(yīng)性較強(qiáng)的特點(diǎn)在故障分析中占據(jù)重要的地位。兩參數(shù)威布爾分布描述的設(shè)備在t>0的任意時(shí)間內(nèi)都可能發(fā)生故障。然而，大量的工程實(shí)踐證明，許多設(shè)備在投入使用后一段時(shí)間段內(nèi)不會(huì)發(fā)生任何故障，這些設(shè)備的故障數(shù)據(jù)經(jīng)過(guò)威布爾變換后在威布爾概率紙上呈現(xiàn)的不再是一條直線，此時(shí)如果仍然采用兩參數(shù)威布爾分布模型擬合故障數(shù)據(jù)就可能會(huì)給可靠性分析帶來(lái)較大的誤差。三參數(shù)威布爾分布尤其適用于這種情況，能夠提供更精確的分析結(jié)果。

三參數(shù)威布爾分布的參數(shù)估計(jì)相比兩參數(shù)威布爾分布的參數(shù)估計(jì)要復(fù)雜得多，其最常用的參數(shù)估計(jì)方法是極大似然法[3]。三參數(shù)威布爾分布的似然方程組是非線性的，結(jié)構(gòu)十分復(fù)雜，常采用Newtom-Raphson迭代等方法求解。用這種方法求解時(shí)，初選值的選擇不當(dāng)會(huì)導(dǎo)致迭代過(guò)程發(fā)散，無(wú)法得到最終結(jié)果。另外，三參數(shù)威布爾分布的參數(shù)估計(jì)方法還有雙線性回歸法、概率權(quán)重矩法、相關(guān)系數(shù)優(yōu)化法、灰色估計(jì)法、Bayes統(tǒng)計(jì)分析法等[4-8]，這些方法針對(duì)不同樣本容量的適應(yīng)能力各不相同，但都難免復(fù)雜的計(jì)算。本文提出一種綜合圖解法和遺傳算法的三參數(shù)威布爾分布參數(shù)估計(jì)方法，其中圖解法為遺傳算法提供搜索范圍，而遺傳算法可通過(guò)迭代獲得更精確的參數(shù)估計(jì)值。

1參數(shù)估計(jì)方法

1.1圖解法

三參數(shù)威布爾分布的可靠度函數(shù)表達(dá)式為

(1)

式中：β為形狀參數(shù)，η為尺度參數(shù)，γ為位置參數(shù)(t<γ時(shí)無(wú)故障)。

對(duì)式(1)兩次取自然對(duì)數(shù)后可以轉(zhuǎn)換為

ln{-ln[R(t)]}=βln(t-γ)-βln(η)

(2)

對(duì)式(2)作威布爾變換，即令y=ln{-ln[R(t)]}，x=lnt，并不能得到線性函數(shù)，但可以利用威布爾概率紙進(jìn)行參數(shù)估計(jì)。根據(jù)式(2)的性質(zhì)，當(dāng)t-γ=0，即x=ln(γ)時(shí)，y趨向于-；而當(dāng)t趨向于時(shí)，y=βx-βln(η)。即x=lnγ和y=βx-βln(η)是式(2)經(jīng)過(guò)威布爾變換后的得到的函數(shù)的兩條漸近線。由此可得三參數(shù)威布爾分布參數(shù)估計(jì)的步驟如下:

(2)在故障數(shù)據(jù)的WPP圖上擬合兩條漸近線：其中一條在所有散點(diǎn)左側(cè)，垂直于橫坐標(biāo)，其表達(dá)式為x=x0；另一條是當(dāng)x趨向于∞時(shí)的漸近線，其表達(dá)式為y=kx-b。根據(jù)三參數(shù)威布爾分布的圖形性質(zhì)可得以下方程組。

(3)

1.2遺傳算法

圖1　遺傳算法流程圖 Fig.1Flow Chart of GA

遺傳算法(GA)具有較強(qiáng)的全局優(yōu)化能力，是一種自適應(yīng)的、智能的搜索技術(shù)，其最成功的應(yīng)用領(lǐng)域是復(fù)雜的非線性優(yōu)化問(wèn)題。遺傳算法是一種群體型操作，該操作以群體中的所有個(gè)體為對(duì)象，它的基本流程如圖1所示，包含6個(gè)基本要素。

(1)參數(shù)編碼。編碼是解空間向GA空間的映射，是連接問(wèn)題與算法的橋梁。編碼的方法有二進(jìn)制編碼，實(shí)數(shù)編碼，區(qū)間值編碼等。編碼方法的選擇將直接影響選擇、雜交、變異等算子的設(shè)計(jì)，也影響算法的復(fù)雜性。

(2)生成初始群體。遺傳迭代前，需要首先建立一個(gè)由若干初始解(也稱為個(gè)體)組成的群體，即初始群體。初始群體中的每個(gè)個(gè)體都是通過(guò)隨機(jī)方法得到的參數(shù)編碼。

生成初始群體前通常要確定每個(gè)參數(shù)的范圍，稱作搜索空間。搜索空間必須足夠大，能夠涵蓋每一個(gè)適合的個(gè)體，但是太大的搜索空間又會(huì)影響收斂的速度，因此搜索空間的定義非常重要。對(duì)于三參數(shù)威布爾分布的各個(gè)參數(shù)的搜索空間，這里利用圖解法產(chǎn)生的參數(shù)估計(jì)值獲得。首先給定寬松系數(shù)θ(0<θ<1)，定義分布中每個(gè)參數(shù)的搜索空間為

(4)

式中:θ值越大，搜索范圍越大；反之，越小。

(3)適應(yīng)度評(píng)估。遺傳算法在遺傳迭代時(shí)基本不利用外部信息，僅以適應(yīng)度函數(shù)為依據(jù)，利用種群中的每個(gè)個(gè)體的適應(yīng)度值進(jìn)行搜索。適應(yīng)度值越大的個(gè)體被選擇的概率越大。因此，適應(yīng)度函數(shù)的選取將直接影響個(gè)體的質(zhì)量，并同樣會(huì)影響收斂速度。

對(duì)于三參數(shù)威布爾分布，本文在d檢驗(yàn)(柯爾莫哥洛夫-斯米爾諾夫檢驗(yàn)，或k-s檢驗(yàn))[9-10]的基礎(chǔ)上經(jīng)過(guò)適當(dāng)?shù)奶幚碜鳛檫m應(yīng)度函數(shù)。

d檢驗(yàn)法是可靠性分析中常用的檢驗(yàn)方法，它比χ2檢驗(yàn)法精細(xì)，而且還適用于小樣本的情況。d檢驗(yàn)法是將n個(gè)試驗(yàn)數(shù)據(jù)按由小到大的次序排列，根據(jù)假設(shè)的分布，計(jì)算每個(gè)數(shù)據(jù)對(duì)應(yīng)的R0(ti)，將其與經(jīng)驗(yàn)分布函數(shù)Rn(ti)進(jìn)行比較，其中差值的最大絕對(duì)值即檢驗(yàn)統(tǒng)計(jì)量Dn的觀察值。將Dn與臨界值Dn,α進(jìn)行比較。若滿足式(5)，則接受原假設(shè)；否則，拒絕原假設(shè)。

(5)

式中：R0(t)為原假設(shè)分布函數(shù)，即三參數(shù)威布爾分布，函數(shù)中的參數(shù)對(duì)應(yīng)于種群中個(gè)體的染色體段。

由此可得三參數(shù)威布爾分布的適應(yīng)度函數(shù)的表達(dá)式如下

Fitness(t)=Dn,α-Dn

(6)

適應(yīng)度函數(shù)的取值越大，表明原假設(shè)分布函數(shù)越貼近經(jīng)驗(yàn)分布函數(shù)，即獲得的分布函數(shù)越精確，在遺傳迭代中被選擇的幾率越大。而且這個(gè)適應(yīng)度函數(shù)還有一個(gè)重要的優(yōu)勢(shì)：如果對(duì)于某個(gè)體使Fitness(t)>0，則可以判定該個(gè)體中的參數(shù)估計(jì)值通過(guò)d檢驗(yàn)，可以在遺傳迭代中直接獲得符合d檢驗(yàn)的參數(shù)估計(jì)值。

(4)選擇。在對(duì)個(gè)體適應(yīng)度進(jìn)行評(píng)價(jià)的基礎(chǔ)上把優(yōu)化的個(gè)體直接遺傳到下一代，或通過(guò)配對(duì)交叉產(chǎn)生新的個(gè)體，再遺傳到下一代。

(5)交叉。交叉操作是遺傳算法中最主要的操作。交叉的步驟為：首先對(duì)種群中個(gè)體進(jìn)行隨機(jī)配對(duì)；然后，在配對(duì)個(gè)體中隨機(jī)設(shè)定交叉處(一點(diǎn)或多點(diǎn))，配對(duì)個(gè)體彼此交換部分信息。

(6)變異。變異操作是按位(bit)進(jìn)行的，即把染色體某一位的內(nèi)容進(jìn)行編譯。編譯操作同樣也是隨機(jī)進(jìn)行的。一般而言，變異概率都取得較小。變異操作是十分微妙的遺傳操作，它需要和交叉操作配合使用，目的是挖掘群體中個(gè)體的多樣化，克服有可能限于局部解的弊病。

選擇、交叉、變異均屬遺傳算法的基本操作，方法繁多，但僅與收斂速度有關(guān)，在此不作贅述，請(qǐng)參看相關(guān)文獻(xiàn)。

2在可靠性分析中的應(yīng)用

計(jì)算機(jī)數(shù)控系統(tǒng)(Computer Numerical Control)是數(shù)控機(jī)床的核心部分，主要有計(jì)算機(jī)數(shù)控裝置(CNC)、伺服驅(qū)動(dòng)裝置、位置檢測(cè)裝置、可編程控制器PLC和接口電路組成。計(jì)算機(jī)數(shù)控系統(tǒng)的故障往往會(huì)帶來(lái)較大的危害度，其可靠性對(duì)于數(shù)控機(jī)床至關(guān)重要。表1給出了某系列數(shù)控機(jī)床數(shù)控系統(tǒng)出廠使用半年的故障數(shù)據(jù)。

圖2 給出了數(shù)控機(jī)床計(jì)算機(jī)數(shù)控系統(tǒng)故障間隔時(shí)間的散點(diǎn)圖，從圖中可以看出散點(diǎn)并不呈一條直線，具有明顯的拐點(diǎn)，符合三參數(shù)威布爾分布的特征。故而假設(shè)計(jì)算機(jī)數(shù)控系統(tǒng)的故障間隔時(shí)間服從三參數(shù)威布爾分布。按照?qǐng)D解法的步驟，在WPP圖上作兩條漸進(jìn)線L1和L2，見圖2。兩條漸近線表達(dá)式分別為：

(7)

圖2　計(jì)算機(jī)數(shù)控系統(tǒng)故障數(shù)據(jù)WPP及漸近線圖 Fig.2 WPP and asymptotic line graph of CNC system

序號(hào)tiRn(ti)xiyi11040.96734.6444-3.403521220.92064.8040-2.49173171———41710.82715.1417-1.661652570.78045.5491-1.394462730.73365.6095-1.17217282———82820.64025.6419-0.807493070.59355.7268-0.6505103740.54675.9243-0.5045113860.50005.9558-0.3665124940.45336.2025-0.2341135010.40656.2166-0.1053145150.35986.24420.0219155760.31316.35610.1495166740.26646.51320.2798176990.21966.54970.4160187010.17296.55250.5625198030.12626.68840.72762010360.07946.94310.92932111980.03277.08841.2297

由此可得三參數(shù)威布爾分布圖解法的初始解為：β=1.42，η=530.94，γ=99.48。

令寬松系數(shù)θ=0.5，可得三參數(shù)威布爾分布的搜索空間為

(8)

運(yùn)用遺傳算法進(jìn)行遺傳迭代可得三參數(shù)威布爾分布的參數(shù)估計(jì)值及適應(yīng)度見表2，其適應(yīng)度大于0且遠(yuǎn)大于圖解法得到的適應(yīng)度，說(shuō)明遺傳算法的參數(shù)估計(jì)值不但能通過(guò)d檢驗(yàn)，而且計(jì)算精度遠(yuǎn)大于圖解法。

表2　參數(shù)的估計(jì)值及其適應(yīng)度

因此，計(jì)算機(jī)數(shù)控系統(tǒng)的故障間隔時(shí)間服從三參數(shù)威布爾分布的假設(shè)成立，其可靠度函數(shù)為

(9)

圖3給出了計(jì)算機(jī)數(shù)控系統(tǒng)的可靠性函數(shù)曲線，圖中的散點(diǎn)即為故障時(shí)間點(diǎn)的經(jīng)驗(yàn)分布函數(shù)值。圖中的散點(diǎn)均勻的分布在曲線的兩側(cè)，說(shuō)明三參數(shù)威布爾分布及文中給出參數(shù)估計(jì)方法對(duì)于計(jì)算機(jī)數(shù)控系統(tǒng)的故障分析是適用的。

圖3　計(jì)算機(jī)數(shù)控系統(tǒng)的可靠度函數(shù)圖 Fig.3 Reliability graph of CNC system

3結(jié)論

首先揭示了三參數(shù)威布爾分布在可靠性分析中的特殊意義，然后針對(duì)其參數(shù)估計(jì)復(fù)雜的問(wèn)題，提出了一種綜合圖解法和遺傳算法的參數(shù)估計(jì)方法。該方法以圖解法的解為基礎(chǔ)形成搜索空間，用遺傳算法提高解的精度，從而得到更精確的參數(shù)估計(jì)值。隨后三參數(shù)威布爾分布應(yīng)用于計(jì)算機(jī)數(shù)控系統(tǒng)的可靠性分析中，驗(yàn)證了本文提出的參數(shù)估計(jì)方法的可行性。本文提出的參數(shù)估計(jì)方法對(duì)于其他復(fù)雜的威布爾分布形式(例如競(jìng)爭(zhēng)威布爾分布)的參數(shù)估計(jì)具有一定的啟示意義。

參考文獻(xiàn)

[1]張英芝，申桂香，吳甦，等. 隨機(jī)結(jié)尾數(shù)控機(jī)床三參數(shù)威布爾分布模型[J]. 吉林大學(xué)學(xué)報(bào)(工學(xué)版),2009,39(2):378-381.

ZHANG Ying-zhi, SHEN Gui-xiang, WU Su, et al. 3-parameter Weibull distribution for random truncated NC machine tool falt data[J]. Journal of Jilin University (Engineering and Technology Edition), 2009,39(2):378-381.

[2]Murthy D N P, Xie M, et al. Weibull models [M]. New York ,John Wiley,2004.

[3]楊謀存，聶宏. 三參數(shù)Weibull分布參數(shù)的極大似然估計(jì)數(shù)值解法[J]. 南京航空航天大學(xué)學(xué)報(bào), 2007,39(1):22-25.

YANG Mou-cun, NIE Hong. Advanced algorithm for maximum likelihood estimation of three parameter weibull distribution[J]. Journal of Nanjing University of Aeronautics & Astronautics, 2007,39(1):22-25.

[4]嚴(yán)曉東，馬翔，鄭榮躍. 三參數(shù)威布爾分布參數(shù)估計(jì)方法比較[J]. 寧波大學(xué)學(xué)報(bào)(理工版),2005,18(3):301-305.

YAN Xiao-dong, MA Xiang, ZHENG Rong-yue. Comparion of the parameters estimation methods for 3-parameter weibull distribution[J].Journal of Ningbo University (NSEE), 2005,18(3):301-305.

[5]Tsionas E G. Posterior analysis, prediction and reliability in three-parameter Weibull distributions [C]. Communications in Statistic-Theory and Methods, 2000,29(7):1435-1449.

[6]楊志忠，劉瑞元. 三參數(shù)Weibull分布參數(shù)估計(jì)求法改進(jìn)[J]. 工程數(shù)學(xué)學(xué)報(bào), 2004,21(2):281-284.

YANG Zhi-zhong, LIU Rui-yuan. Improvement of parameter estimation method for Three-parameter Weibull distribution[J]. Journal of Engineering Mathematics,2004,21(2):281-284.

[7]Teimouri M, Gupta A K. On the Three-parameter weibull distribution shape parameter estimaion[J]. Journal of Data Science,2013,11:403-414.

[8]盧安平,趙林,郭增偉,等.基于Monte Carlo法的極值分布類型及其參數(shù)估計(jì)方法比較[J]. 哈爾濱工業(yè)大學(xué)學(xué)報(bào), 2013,45(2):87-95.

LU An-ping, ZHAO Lin, GUO Zeng-wei, et al. A comparative study of extreme value distribution and parameterestimation based on the Monte Carlo method[J]. Journal of Harbin Institute of Technology, 2013,45(2):87-95.

[9]陸陽(yáng),王強(qiáng), 岳峰, 等.軟件可靠性關(guān)鍵參數(shù)判定的據(jù)估計(jì)方法[J]. 系統(tǒng)工程理論與實(shí)踐,2014,34(1): 248-254.

LU Yang,WANG Qiang, YUE Feng, et al. Key parameter judgment with moment estimation in software reliability prediction[J].Systems Engneering-Theory & Practice,2014,34(1): 248-254.

[10]Xing L, Boddu P, Sun Y, et al. Reliability analysis of static and dynamic fault-tolerant systems subject to probabilistic common-cause failures [J]. Proceedings of the Institution of Mechanical Engineers part 0-Journal of Risk and Reliability,2010,224(1):43-45.