1、http:/ - 1 - 威布尔分布参数估计在 EXCEL 中的实现方法研究 史景钊，花恒明，李祥付 河南农业大学机电工程学院，郑州（450002 ） E-mail: 摘 要： 三参数威布尔分布的参数估计比较复杂， 大多数估计方法都需要编程计算。 推导了 相关系数优化法求解三参数威布尔分布位置参数的公式，并把该公式利用 MS EXCEL 强大 的规划求解进行求解， 实现了位置参数估计的非编程求解， 同时利用图表功能实现了形状参 数和尺度参数的非编程求解。实例证明该法方便、快捷，可以满足估计精度。 关键 词： 可 靠性；威布尔分布；参数估什；EXCEL 中图 分类号：TB114.3 1. 引言
2、 威布尔分布是瑞典物理学家 Weibull W. 分析材料强度时在实际经验的基础上推导出来 的分布形式 1 ， 国内外大量研究表明， 用三参数威布尔分布比用对数正态分布往往能更准确 地描述结构疲劳寿命或腐蚀损伤的概率分布 2 ， 物理意义更加合理； 在以损耗为特征的机械 零件寿命评估中， 采用三参数威布尔分布比采用二参数威布尔分布拟合精度更高。 因此， 三 参数威布尔分布在强度与环境研究领域及机械零件磨损寿命评价中得到越来越广泛的应用。 在农业机械的强度设计中也经常要用到威布尔分布。 威布尔分布参数估计方法有很多， 国内外一直有人在进行相关研究 3-8 ， 现有几十种参 数估计方法， 但多数只
3、能用于形状参数和尺度参数的估计。 在众多的估计方法中， 能用于三 参数估计的并不多，见诸文献的有极大似然估计法、最大相关系数优化法、概率权重矩法、 灰色估计法、 图估计法等， 除图估计法外， 其他方法大都计算复杂， 应用不便， 即便是计 算 机水平发达的今天， 也只能通过 Matlab 或其他计 算机语言编程计算。EXCEL 提供了 超强的 数学运算 、 统计分析 等 实用程序 ， 利用它的 规 划求解功 能 可以快速 、 高效地求 解 三参数威 布尔分布的参数估计问题。 2. 三参数威布尔分布模型 威布尔分布的寿命分布函数由下式给出 (1) t exp 1 ) ( = m t t F 式中：
4、m 称 为形状参数，m0 ； 称 为尺度参数，0 ； 称为位置参数，也称最小寿 命， 表示产品在 以前不会失效， 对于产品寿命有 0 , =0 时退化 为二参数威布尔分布； t 是产品的工 作时间， t 。 当 m1 时，失效率是递增的， 适合于建模磨耗或老化失效。 设有 n 个产品进行寿命试验数据，按失效时间先后得到的寿命数据失效时间( 顺序统计 量) 为 n t t t 2 1 ，对应的累计失效概率( 经验分布函数) 为 ) ( ) ( ) ( 2 1 n t F t F t F 。 其中到第 i 个 产品失效时的累计失效概率 ) ( i t F 可用中位秩算法求得： http:/ - 2
5、 - 4 . 0 3 . 0 ) ( + n i t F i（2 ） 根据失效时间和累计失效概率即可用各种方法对其参数进行估计。 3. 最大相关系数优化法 对（1 ）式做 变形处理，并取两次自然对数得到： m t m t F ln ) ln( ) ( 1 1 ln ln = （3 ） 令： m B t x t F y ln ), ln( , ) ( 1 1 ln ln = = = 于是（3 ）式 化为线性方程 y=mx-B （4 ） 令： n i t F y t x i i i i , , 2 , 1 , ) ( 1 1 ln ln ); ln( = = = 则在 已知时可利用最小二乘法进行回
6、归求出 m 和 B ， 进而求出尺度参数 ) exp( m B = 。 x 与 y 间的相 关系数 R (x, y) 为： 1 2222 11 () (,) () () n ii i nn ii ii xy n x y Rxy x nx y ny = = = （5 ） 式中： 1 1 n i i x x n = = ， 1 1 n i i y y n = = 由(5) 式可知，相关系数 R (x, y) 是位 置参数 的函数，我们要求的 应是使 R (x, y) 最大 的那一个。 根据求极值的方法， 使 0 / ) , ( = d y x R d 的 即为所求， 0 / ) , ( = d y
7、 x R d 等 价于 0 / ) , ( 2 = d y x dR 。 若令： 2 1 () n ii i vx y n x y = = ， 2222 11 () () nn ii ii uxn xyn y = = 则有： 2 2 ) , ( u vu uv d y x dR = 即需满足： 0 = vu uv 因为： 111 2( ) ln( ) ln( ) nnn ii i i i iii dv d vx y n x yt y y t dd = = 11 2( ) nn i ii ii i yy xy n x y t = = http:/ - 3 - 2 22 2 111 22 11 1
8、 ( ) ln ( ) ln( ) 2( ) nnn iii iii nn i i ii i du d uy n ytt ddn xx yn y t = = = = 所以： 22 22 111111 2 ( ) ( ) ( ) ( ) nnnnnn ii ii ii i i iiiiii ii yy xx u vv u yn y x yn xy xn x x yn xy tt = = 因为 0|R| 1 ，所以 22 1 () 0 n i i yn y = ，否则 R= ； 1 () 0 n ii i xy n x y = ，否则 R=0 。 因此只有: 22 1111 ()()0 nnnn
9、ii ii i iiii ii yy xx xn x x yn xy tt = = （6 ） 满足（6 ）式 的 即为所求的位置参数。 4. 用 EXCEL 进行参数估计 （6 ）式所 表 示的方 程十 分复杂 ，解 该方程 一般 是通过 编程 ，用数 值解 法求出 ，然 后 求再用最小二乘法或其他方法求解形状参数和尺度参数。MS EXCEL 具有强大的统计和计 算功能，其“规划求解” 功 能更是求解 最优化问题 的强有力工 具，（6 ）式 所表示的方 程利用 EXCEL 的“ 规划求解” 功能可很容易解出，然后再利用其散点图的趋势线功能即可求出形状 参数和尺度参数。本文通过实例，就相关系数优
10、化法，用 EXCEL 进行求解。 例： 选取 5 台某产品进行可靠性试验，失效时间分别是 27,32,36,42,49 ， 已知产品寿命 服从威布尔分布，试估计分布参数。 1) 准备数据 表 按图 1 准备 数据表 在 A2 A6 单元格中输 入产品失效的顺序号 1 5 ； 在 B2 B6 单元格中输 入产品的失效时间 27 、32 、36 、42 、49 ； 图 1 位置参数估计数据准备表http:/ - 4 - 在 I8 单元 格中输入位置参数 的迭代初值， 初值可选择接近于第一个失效时间， 也 可用图估计法的估计值作为初值； 在 C2 单 元格中输入公式“=LN(B2-$I$8)” ，用
11、 填充柄填充 C3 C6 单元 格，C2 C6 单元格的值为 i x ，即 ) ln( i t ； 在 D2 单 元格中输入公式“=C2*C2” ，用填充柄填充 D3 D6 单元格，D2 D6 单元 格的值为为 2 i x ； 在 E2 单 元格中输入公式“=(A2-0.3)/5.4” ， 用填 充柄填充 E3 E6 单元格 ，E2 E6 单 元格的值为为 ) ( i t F ，这里 ) ( i t F 采用中位值算法，即 ) 4 . 0 ( ) 3 . 0 ( ) ( + = n i t F i ； 在 F2 单 元格中输入公式“=LN(LN(1/(1-E2)” ，用填充柄填充 F3 F6
12、单元格，F2 F6 单元格的值为为 i y ，即 ) ( 1 1 ln ln i t F ； 在 G2 单 元格中输入公式“=C2*F2” ， 用填充柄填充 G3 G6 单元格， G2 G6 单元格 的值为 i i y x ； 在 C7 单 元格中输入公式“=AVERAGE(C2:C6)” ，C7 单元格 的值为 x ； 在 C8 单 元格中输入公式“=SUM(D2:D6)” ，C8 单元格的值 为 = n i i x 2 ； 在 C9 单 元格中输入公式“=AVERAGE(F2:F6)” ，C9 单元格 的值为 y ； 在 C10 单 元格中输入公式“=SUM(G2:G6)” ，C10 单元
13、格的值 为 = n i i i y x ； 在 H2 单元格中输入公式“=($C$7-C2)/(B2-$I$8)” ，用填充柄填充 H3 H6 单元格， H2 H6 单元 格的值为 i i t x x ； 在 I2 单元 格中输入公式“=($C$9-F2)/(B2-$I$8)” ， 用填充柄填充 I3 I6 单元格，I2 I6 单元格的 值为 i i t y y ； 在 F8 单 元格中输入公式“=SUM(H2:H6)” ，F8 单元格的值为 = n i i i t x x 1 ； 在 F10 单 元格中输入公式“=SUM(I2:I6)” ，F10 单元格的值为 = n i i i t y y
14、 1 ； 在 I10 单 元格中输入公式“=(C8-5*C72)*F10-(C10-5*C7*C9)*F8” ，I10 单元格的值 就是（6 ）式 的左边。 其他文字仅用于说明，与求解关系不大，可以不填。 2 ）使用“ 规划求解” 功能估计位置参数 选择“ 工具 规划求解” 功能打开规划求解参数对话框，目标单元格设为$I$10 ，目 标值设为 0 ， 可变单元格设为$I$8:$B$2 。 对于产品 寿命有 1 0 t =0 ； http:/ - 5 - 单击“ 求解” 按钮，即可获得最大相关系数下的位置参数 =20.2395 ，如 图 2 所示， 此时可获得最大相关系数 R(x,y)=0.99
15、950878 ； 3 ）使用图表 功能求形状参数 m 和尺 度参数 插图散点 图，横坐标为 x i ，纵坐 标为 y i ； 在散点图 上添加趋势线，回归模型选择“ 线性” ，并选择“ 显示公式” 和“ 显示 R 2 值” ； EXCEL 自 动绘制回归直线， 并把结果显示在图上， 结果如图 3 所示。 其 中斜率 1.8486 即为形状参数 m ，而 5.5088 即为 m ln ，故 69 . 19 8486 . 1 5088 . 5 exp = = 5. 结语 （1 ）三 参数威布 尔 分布的参 数 估计较二 参 数威布尔 分 布的参数 估 计更为复 杂 ，用 MS EXCEL 解决三参
16、数威布尔分布的参数估计问题，实用方便。 （2 ）用“ 规划求解” 进行位置参数的估计时，要注意选择好迭代精度和初值，若精度和 初值选择不合适， 可能得不到最满意的解， 建议初值尽可能离第一个失效时间近一些， 精度 不低于 10 -7 。 （3 ）失效概 率有中位秩算法、平均秩算法，采用的算法不一样，估计结果也会稍有不 同。 图 3 回归结果 图 2 规划求解 结果 http:/ - 6 - 参考 文献 1 Weibull W. A statistical distribution function of wide applicabilityJ . Journal of Applied Micr
17、oelectron Reliability. 1951,28(4) ：613-617 2 Hallinan A J. A review of the Weibull DistributionJ. Journal of Quality Technology. 1993,25(2) ：85-93 3 Adatia A ，Chart L. K Robust. Estimators of the 3-parameter Weibull DistributionJ. IEEE Transactions on Reliability, 1985,34(3) 4 Garn G W. Moment Estim
18、ators for the 3-parameter Weibull DistributionJ. IEEE Transactions on Reliability ，1988, 37(3) 5 蒋卉，刘瑞 元. 用插值法求 威布尔分布位置参数估计J. 青海大学学报( 自然科学版) 2003,21(3) ：53-57 6 朱铭扬. 三参数威布尔分布的参数估计J. 江苏技术师范学院学报. 2006,12(6) ：31-34 7 胡恩平，罗兴柏，刘国庆. 三参数 Weibull 分布几种常用的参数估计方法J. 沈阳工业学院学报. 2000,19(3)：88-93 8 严晓东， 马翔， 郑荣跃等. 三
19、参数威布尔分布参数估计方法比较J. 宁 波大学学报( 理工版). 2005,18(3) ： 301-305 Study on Weibull Distribution Parameter Estimation in EXCEL SHI Jing-zhao, HUA Heng-Ming, LI Xiang-fu Mechanical and electrical Engineering college , Henan Agricultural University, Zhengzhou (450002) Abstract 3-parameter Weibull distribution of t
20、he parameter estimation more complicated, majority of the estimated methods need programming. Derivation of the correlation coefficient optimization method to solve the 3-parameter Weibull distribution parameters estimation of the formula and the formula to use MS EXCEL powerful Solver to solve a lo
21、cation parameter estimation of non-programming solution, while the use of the chart features a shape parameter and the scale parameter for solving the non-programming. Examples to prove that the method easily accessible, efficient and meet the estimated precision. Keywords: reliability; Weibull distribution; parameter estimation; Excel
