三参数威布尔分布参数估计方法比较_严晓东
第18卷第3期2005年9月
宁波大学学报(理工版)
JOURNAL OF N I N G BO UN I V ERSI TY (NSEE )
Vol . 18No . 3
Sep t . 2005
文章编号:100125132(2005) 0320301205
三参数威布尔分布参数估计方法比较
严晓东, 马 翔, 郑荣跃, 吴 亮
1
2
1
1
(1. 宁波大学建筑工程与环境学院, 浙江宁波315211; 2. 宁波房地产股份有限公司, 浙江宁波315000)
摘要:选择目前常用的5种三参数威布尔分布参数估计法(、双线性回
归法、相关系数优化法、灰色估计法) 进行比较研究, 7, 得出了这5种方法各自相应的适用范围, 指导意义.
关键词:; 双线性回归法; 概率权重矩法; 相关系数优化法;
中图分类号:O213. 2; T B114. 3 文献标识码:A
Comparis on of the Parameters Esti m ati on Methods
for 32Parameter W eibull D istributi on
Y AN Xiao 2dong , MA Xiang , Z HE NG Rong 2yue , WU L iang
1
2
1
1
(1. Faculty of A rchitectural, Civil Engieering and Envir on ment, N ingbo University, N ingbo 315211, China;
2. N ingbo Realestate Co . , L td . , N ingbo 315000, China )
Abstract:The W eibull distributi on is one of the most popular models in the field of reliability . Since 1939, many para meter esti m ati on methods have been p r oposed, including the p r obability weighted moment method , the maxi 2mum likelihood method, bilinear regressi on method, correlati on coefficient method and Grey esti m ati on method . The methods menti oned above are disscussed, p r ogra m s and exa mp les are given .
Key words:32para meterW eibull distributi on; the maxi m u m likelihood method; bilinear regressi on method;
p r obability weighted moments; correlati on coefficient method; grey esti m ati on method
CLC nu m ber:O213. 2; T B114. 3 D ocu m en t code:A
自W. W eibull 开发了三参数模型并应用它建模处理大量的失效数据以来, 威布尔分布模型成了可靠性领域里最广泛使用的模型
[2]
[1]
得它在曲线的拟合上富于弹性; 而它的可靠性基本函数都有封闭形式的解析表达式, 数学处理十分便利, 尤其是经过双对数变换后, 它能被线性化, 从而使得计算机图形处理和线性回归等技术能被方便而广泛地应用. 此外, 正态分布和指数分布等都可以看
. 威布尔分布可
包含递增和递减的失效率, 十分符合在实际应用中碰到的失效问题. 此外, 由于其形状参数的存在, 使
收稿日期:2004-10-15. 基金项目:浙江省教育厅资助项目(2002417) . 第一作者简介:严晓东(1967-) , 男, 汉族, 浙江宁波人, 讲师.
302宁波大学学报(理工版) 2005
成是威布尔分布的特例, 因此研究威布尔分布有广泛的实际意义
[3]
.
然而, 由于同时存在3个参数, 使得三参数威布尔分布参数估计有一定困难, 国内外一直有人在进行相关研究, 现有几十种参数估计方法中, 各自的适用范围和精度都不同, 给工程选用带来不便, 因此有必要对各种方法进行比较. 本文在综合分析基础上, 对目前比较常见的5种参数估计法(极大似然法、双线性回归法、概率权重矩法、相关系数优化法和灰色估计法) 进行了比较研究, 相关结论对在不同的场合用不同的威布尔分布数据类型进行拟合时, 正确地选用参数求解方法有一定的指导意义.
式中n 为试验时的样本容量, x i (i =1, 2, Λ, n ) 为样本的观测值, 如失效时间.
由此构成三参数W eibull 分布的似然方程为:
n
-(β-1) -i =1x -γγθ-γ9i
n
ββ
) +(β-1) ∑(x i -γi =1(θ-γ)
n
β(β-1)
) =0, β∑(x i -γ
(θ-γ) i =1
n
(θ-n ln -+∑ln (x i -γ) +
i =1β9
n β
) --(i =1
(6)
1 三参数威布尔分布的分布函数为
β
η]}, (x ≥γ) (1) F (x ) =1-exp {-[(x -γ) /
β和γ分别称为尺度参数、式中η、形状参数和位置参数, 且η, β>0, γ≥0.
威布尔分布的密度函数和失效率函数分别为
β
η) [(x -γ) /η]-1・f (x ) =(β/
exp {-[(x -γ) /η]},
β
η) [(x -γ) /η]-1, r (x ) =(β/
β
i =1
∑[(x i -γ) ln (x i -γ) ]
(θ-γ)
β
n
β
=0,
β=-+・
θθ-γ(θ-γ) (β+1) 9
∑(x i -γ) =0.
i =1n
β
(2) (3)
当γ=0时, 模型退化为二参数模型.
以上诸式中, 当β1时, 失效率是增的, 适合于建模磨耗或老化失效.
显然似然方程组为非线性方程组, 且结构十分复杂, 直接获得解析解得不可能的, 必须采用计算机
[5]
进行求解. 应用计算机实现极大似然估计实际上即对(6) 式给出的非线性方程组进行求解. 非线性方程组在计算机上求解大都采用Ne wt on 2Raphs on 迭代法.
基于Ne wt on 2Raphs on 迭代法的基本思想, 首先给出待估计参数γ, β, θ的初选值γ0, β0, θ0, 且令γΔγ=γΔβ=βΔθ=θ; β; θ, 然后在γ0+0+0+0, β0, θ0处将似然方程组的左端各项进行级数展开, 并作一阶近似, 则方程组(6) 转化为线性方程组.
γ, Δβ和Δθ的解后, 进求解线性方程组获得Δ
γ, Δβ和Δθ行判断, 如果Δ均小于给定的误差界限, γ, Δβ和Δθ即为估计值, 否则用γ则对应的Δ0+
Δγββ和θθ代替Δγ, Δβ和Δθ重新运算, 、0+Δ0+Δ
反复循环, 最终得γ, β和θ的估计值. 存在的问题是, 如果初选值选择不当将导致迭代过程发散, 而无法得到最终结果. 常采用图解法所获得的估计结果作为初选值.
[6]
2. 2双线性回归法
β1/
变换式(1) , 并令η=x 0, 得到如下2个方程:
(7) ln ln =β・ln (x -γ) -ln x 0,
1-F (x ) γ1/β
(ln ) (8) =-.
ηη1-F (x )
2 三参数威布尔分布的参数估计
2. 1极大似然法
极大似然法是一种十分有效和通用的参数估计方法. 其基本思想是选择待定参数使样本出现在观测值的领域内的概率最大, 并以这个值作为未知参数的点估计值.
先记
η=θ-γ,
改写式(2) , 再根据极大似然估计的基本原理成似然函数:
β-n βln (θ-γ) +(β-
1) ・ln L =n ln n n
β(γ) ∑ln x i --β∑(x i -r ) , i =1(θ-γ) i =1
[4]
(4) , 构
(5)
第3期 严晓东等:三参数威布尔分布参数估计方法比较
s
303
(8) 线性无关, 利用线性回归方易知上述两式(7) 、
程的最小二乘解法, 合并整理后得到
M 1, 0, 0M 1, 0, 1s
n i =1
n
∑x i , ∑x i (1-n
n
^β
^) ]ln ln [ln (x -γ
-
^2^2
[ln (x -γ) ]-[ln (x -γ) ]
n i =1n i =1
) , n n
(17)
M 1, 0, 3s
3
) . ∑x i (1-[8]
^) ]ln ln [ln (x -γ
,
22^^[ln (x -γ) ]-[ln (x -γ) ]
[x -x ]ln
x ln
2
2
2. 4相关系数优化法
(9)
1/β
改写式(1) 为1-F (x ) =-^
,
^1-βη
, (18)
^=x -γ
1-F (x 1-β
^-x ln
1F (x ^=β^・[ln (x -x 0
ln 1-F (x . (11)
(10) 两式可分别表示为一元函数关系, 式以上(9) 、
(11) 可表示为二元函数关系, 即
η, (19) (F x ) ]=(x -γ) -βln
:
Y ln[-ln (1-F (x ) ) ],X =ln (x -γ) ,
(20)
η, A =-βln
B =β.
(21) 可得Y =A +B X.
变量X 与Y 之间成线性关系, 根据观测数据(x i , F (x i ) ) , 通过式(20) 得到新的数据(X i , Y i ) , 于是由线性回归分析得到特定参数A 、B 和线性相关系数r 分别为
^=Y -B X, B ^=L /L, A X Y XX
(22)
r =L X Y /XX L YY ,
X =
^=f (γ^) , β1^=f (β^) , γ2^, γ^) . x 0=f 3(β
(12) (13) (14)
前两式表示两个曲线方程, 这两条曲线必定相交且仅交于一点, 若不相交则可断定数据不服从威布尔
(13) 式为超越方程组, 故先分布. 由于联立的(12) 、
按给定求解精度用迭代的方法求得β与γ的近似^. 估计算. 代入式(14) , 得到x 02. 3概率权重矩法(P WM )
[7]
n i =1
n i =1n
∑X i , Y 2
2
n
n i =1
∑Y i ,
n
L XX =∑X 1-n X , L YY =∑Y 1-n Y ,
i =1n
威布尔分布的概率权重矩式为
M 1, 0, k
22
ηΓ(1+)
γβ, 1+k (1+k ) 1+β
(15)
L X Y =∑X i Y i -n X Y .
i =1
由上式可见A 、B 和r 为γ的函数, 记为A (γ) 、
B (γ) 和r (γ) , 要求γ, 必须使线性相关系数r 的绝
式中Γ为Γ函数, k =0, 1, 3. 写出M 1, 0, 1、M 1, 0, 2、
M 1, 0, 3并联解得
γ=0, 通过数值计对值|r |取最大值, 由d |r (γ) |/d
2
4(M M -M )
γ=,
4M 1, 0, 3+M 1, 0, 0-4M 1, 0, 1
算方法可求得γ, 然后由式(22) 求得A 、B 和r , 再由式(20) 求得β和η的值. 2. 5灰色估计法
(16)
[9]
η,
M 1, 0, 0-2M 1, 0, 1
ln ) /M 1, 0, 1-2M 1, 0, 3β.
M -2M ln
2(M 1, 0, 1-2M 1, 0, 3)
M 1, 0, 0-γ
记R (x ) =1-F (x ) , 再对式(1) 两边取自然对数有
ln ln
R (x )
=βη
, (23)
观测样本的概率权重矩为
ln ln
变换得 x =ηe βR (x ) +γ,
304宁波大学学报(理工版) 2005
令 t i =-ln ln
R (x i )
,
48025, 48055, 48055, 48055, 48055, 48056, 51675, 52344, 52345, 52345, 52345, 52379, 55997, 56202, 57709, 57709, 57709, 57709, 63496
-t
则 x t i =ηe βi +γ.
-at
记η=c, 1/β=a, γ=b, 则公式为x =ce i +b . 这与
灰色模型G M (1, 1) 的结果完全一致. 因此可以用
G M (1, 1) 的直接建模法求得各参数.
算例4 子样数n =71, 某产品寿命试验数据3956. 42, 4004. 18, 4091. 61, 4355. 05, 4355. 40, 4376. 01, 4391. 79, 4487. 68, 4487. 68, 4736. 67, 4736. 67, 4939. 85, 4963. 62, 5220. 19, 5353. 41, 5372. 72, 5418. 04, 5444. 11, 5603. 17, 5698. 10, 5. 5843. 6, 6197. 41, 6249. 69, 279. 6740. 48, 6887. 65, , , 7209, 7209, 7209, 7366. 4, 581. 64, 7581. 64, 7581. 64, 7645. 59, 8246, 8599. 7, 8713. 97, 8936. 34, 9044. 22, 9197. 45, 9511. 73, 9754. 47, 9967. 45, 10136. 31, 10172. 88, 10172. 88, 10308. 04, 10395, 10609. 23, 10609. 23, 10788. 97, 10879. 97, 10971. 75, 11594. 41, 11990. 59, 12237. 31, 12400. 31, 12400. 31, 12550. 01, 13198. 73, 13947. 78, 15557. 12, 17646. 12, 19848. 23, 23199. 07.
3 算例与结果分析
为了验证方法对不同子样的适应程度, 进行了大量计算, 部分算例原始数据
[10]
如下:
算例1 子样数n =5, 试验数据(单位:MPa)
381, 395, 408, 算例2 子样数, 某一组铝合金试件的寿命试验数据(单位:千周)
350, 380, 400, 430, 450, 470, 480, 500, 520, 540, 550, 570, 600, 610, 630, 650, 670, 730, 770, 840
算例3 子样数n =31, 某电子产品的寿命试验数据
30926, 34554, 36381, 38423, 40103, 40501, 42200, 44392, 46092, 46125, 46175, 48025,
表1 参数估计结果及检验
Tab . 1 Results of param eter esti m a ti on and ver i f i ca ti on
算例1(n =5) 概率权重矩法极大似然法双线性回归法相关系数法灰色估计法算例2(n =20) 概率权重矩法极大似然法双线性回归法相关系数法灰色估计法算例3(n =31) 概率权重矩法极大似然法双线性回归法相关系数法灰色估计法
形状参数
-13. 79145. 22915. 234117. 600712. 6499
尺度参数
-642. 6586133. 1607133. 1500417. 899258. 4210
位置参数
1080. 7388284. 2496284. 21980. 2135158. 3513
K -S 检验值0. 2566000. 2105820. 2106000. 2094700. 198748K -S 检验值0. 0729850. 0669120. 0670120. 0873910. 059873K -S 检验值0. 1139070. 1137310. 1137400. 1146800. 102877
不等系数
0. 263360. 0426450. 0426450. 0449310. 005213
相关系数
0. 8470950. 9950690. 9950700. 9946830. 999984
形状参数
2. 12481. 99401. 99504. 52412. 0493
尺度参数
330. 2207316. 4807316. 5078608. 8870305. 3332
位置参数
264. 5440280. 8242280. 79840. 0001120288. 7339
不等系数
0. 0223310. 0207070. 0206990. 0571630. 009954
相关系数
0. 9991700. 9991080. 9991000. 9944050. 999950
形状参数
5. 95967. 00066. 99966. 79035. 9155
尺度参数
45489. 548553073. 795353070. 328951782. 535542713. 3705
位置参数
6271. 1115-1290. 9982-1291. 75620. 04108821. 1859
不等系数
0. 0789750. 0790190. 0790200. 0791800. 022323
相关系数
0. 9869060. 9869040. 9869070. 9868800. 999749
第3期续表
算例4(n =71) 概率权重矩法极大似然法双线性回归法相关系数法灰色估计法
严晓东等:三参数威布尔分布参数估计方法比较305
形状参数
1. 28701. 13861. 14162. 74961. 0372
尺度参数
5208. 145344946. 13584943. 26499426. 40584150. 3444
位置参数
3578. 95263821. 50413821. 23590. 59894368. 2277
K -S 检验值0. 0609950. 0582580. 0545910. 1254560. 056887
不等系数
0. 0358310. 0391930. 0391980. 0956010. 052517
相关系数
0. 9874010. 9969300. 9970000. 9875060. 991255
上述4个算例应用前述5种三参数威布尔分布
参数估计的计算结果见表1, 表中采用K -S 检验、不等系数及相关系数对计算结果拟合度进行检验. 表1和大量计算结果表明, 一般情况下, 5种方法都有相当的精度, 可在工程上应用. 如下结论:
(1) , 估计出威布尔分布3个参数, 同其他方法相比具有所需子样少, 拟合精度高的优点, 如最少3个数据即可进行参数估计, 并有一定精度; 其余3种方法离不开迭代求解, 尤其是最大似然法要用迭代法求解联立的3个超越方程, 相当复杂.
(2) 与通常的先估计位置参数再求得形状参数和尺度参数的做法不同, 灰色估计法先估计得到形状参数, 然后再进一步求得位置参数和尺度参数, 这为威布尔分布的参数估计方法提供了新的途径.
(3) 样本越大, 各种参数估计法所得拟合结果就越精确; 当样本量很小时, 除灰色估计法仍有较高精度外, 其他几种参数估计方法的拟合度都下降. 特别是样本数量小于10时, 用灰色估计法所得拟合结果在本文所讨论的5种参数估计方法中是精度最高的. 而概率权重矩法是这5种方法中拟合精度最差的.
(4) 样本数量属于大样本时, 用极大似然法、双线性法所得拟合结果在本文所讨论的5种参数估计法中是精度最高的.
(5) 当样本量属于中等时, 本文所讨论的5种方法所得的拟合精度相差不多, 可视情况选用.
参考文献:
[1]
W eibull W. A statistical distributi on functi on of wide app licability[J ].App lied M icr oelectr on, Re 2liab . , (617.
23]
A of W eibull distributi on [J ].Quality Technol ogy, 1993, 25(2) :85～93. J ing L ing, Jwo Pan . A ne w method f or selecti on of pop 2ulati on distributi on and para meter esti m ati on[J ].Relia 2bility Engineering and Syste m Safety, 1998, 60:247～255.
[4][5][6]
胡恩平. 三参数W eibull 分布几种常用的参数估计方法[J ].沈阳工业学院学报, 2000, 19(3) :88～94. 谷耀新. 三参数威布尔分布参数估计方法[J ].沈阳工业学院学报, 1997, 16(4) :53～57.
庄渭峰. 用微机实现威布尔参数的双线性回归最小二乘估计[J ].电子产品可靠性与环境试验, 1999, 5:
2～7.
[7]Green wood A J. Pr obability weighted moments:defini 2ti on and relati on and t o para meters of several distributi on exp ressable in inverse f or m [J ].W ater Res ources Re 2seach , 1979, 15(5) :1049～1054.
[8]R ichard A L, StephensM A. Esti m ati on and test -of -fit for the three para meter weibull distributi on [R ]. Technichal Report No . 1993.
476.
Stanf ord , California,
[9]郑荣跃, 严剑松, 威布尔分布参数估计新方法研究
[J ].机械强度, 2002, 24(4) :213～215.
[10] 郑荣跃. 灰色建模方法的改进及其在疲劳可靠性研
究中的应用[D ].长沙:国防科技大学硕士学位论文, 1989.
(责任编辑 洪明照)