- Email: [email protected]
A new method for selection of population distribution and parameter estimation
Reliability Engineeringand SystemSafety 60 (1998) 247-255 PIhS0951-8320(97)00171-3
ELSEVIER
© 1998 Elsevier Science Limited All rights reserved. Printed in Northern Ireland 0951-8320/98/$19.00
Technical note A new method for selection of population distribution and parameter estimation Jing Ling & Jwo Pan* Center for Automotive Structural Durability Simulation, Department of Mechanical Engineering and Applied Mechanics, The University of Michigan, Ann Arbor, MI 48109, USA
(Received 18 March 1997; accepted 3 December 1997)
An optimization method for parameter estimation is presented with the KolmogorovSmirnov distance used as the objective. A step-by-step implementation procedure is given. The method is demonstrated by estimating the parameters for three-parameter Weibull distributions from three different samples (with different sample sizes). A comparison of the proposed method and the usual methods such as the least-squares method, the matching moments method and the maximum likelihood method shows that more reasonable estimates of the parameters are given by the proposed optimization method. Then, the proposed method is successfully extended to estimate the parameters for the sum of two three-parameter Weibull distributions. Based on these findings, a new procedure for selection of population distribution and parameter estimation is presented. © 1998 Elsevier Science Limited.
1 INTRODUCTION
parameter estimation is proposed with the K o l m o g o r o v Smirnov distance used as the objective. This is in contrast to Ref. 9, in which minimizing the Kolmogorov-Smirnov distance is only used to estimate the location parameter of three-parameter Weibull distributions. Sometimes, two subpopulations may be mixed into one overall population, a sum model of the two distributions may be a suitable model to characterize the overall population. When three-parameter Weibull distributions are used, it is difficult to estimate the parameters for the model. Jiang and Kececioglu i i studied the behavior of a two-Weibull mixture on the Weibull probability papers and evaluated two graphical parameter estimation methods: the K a o Cran method 12,13 and the Jensen-Petersen method 14, for this kind of overall distribution. Jiang and Murthy 15 did similar research. However, all of these studies are limited to the mixed models of two-parameter Weibull distributions because of the limitation of the probability plotting method. More refined statistical methods are needed to estimate the parameters for the sum of three-parameter Weibull distributions ~5. Our proposed optimization method for parameter estimation is also applicable to the sum of Weibull distributions. In contrast to Refs l 1,15, the proposed optimization method can be used to estimate the parameters for the sum of three-parameter Weibull distributions. Based on
Parameter estimation plays a very important role in the evaluation of reliability. Much research has been done to find more accurate, practical methods for parameter estimation. The probability plotting method l is easiest and straightforward, but less precise than analytical methods. The probability plotting method is usually used to obtain the initial parameter estimation for analytical methods. The most commonly used analytical methods for parameter estimation are the least-squares method l, or the regression method 2 the matching moments method 3-5, and the maximum likelihood method 6,7. There are also other methods presented in literature, for example, the Bayes estimation method 8, the minimum-distance estimation method 9 and the robust estimation method 10 No matter what kind of parameter estimation method is used, the selected distribution with the estimated parameters must be tested to determine whether the distribution adequately and accurately represents the existing data of interest. The best distribution will be the one with the minimum difference test statistic value or the maximum likelihood test statistic value. In this paper, an optimization method for *To whom correspondence should be addressed. 247
Jing Ling, Jwo Pan
248
the proposed optimization method for parameter estimation, a new procedure for selection of population distribution is introduced. This procedure can easily be used to obtain the best population distribution with optimized parameters.
2 AN O P T I M I Z A T I O N M E T H O D OF P A R A M E T E R ESTIMATION No matter what kind of parameter estimation method is used, the selected distribution with the estimated parameters must be tested to determine whether it represents the existing data well. This analysis of the fitness of a selected distribution to the data may be accomplished by the Kolmogorov-Smirnov (K-S) test. The principle and the detailed implementation procedure of the K - S test are introduced in Ref. i. The K - S test is used to determine whether the maximum absolute difference between the observed probability of failure and the expected probability of failure is less than the allowable critical value corresponding to a specified significance level. If the maximum absolute difference is less than the allowable critical value, the selected distribution is accepted as a reasonable distribution. If the maximum absolute difference is equal to or greater than the allowable critical value, the selected distribution is rejected. Furthermore, the allowable critical value decreases with the increase of the significance level. Therefore, when a set of parameters of the selected distribution can be obtained by minimizing the maximum absolute difference between the observed probability of failure and the expected probability of failure, this set of parameters can be regarded as the optimized one for the selected distribution. The procedure to implement the above idea is as follows: 1. Arrange all data in ascending order, or rank them. 2. Calculate the observed probability of failure. There are three different rank definitions that can be used to calculate the observed probability of failure. The median rank as defined in eqn (1) is used here because of its wide acceptance.
Fo(Tg ) -
i-0.3 i = l, 2, n+0.4 . . . . .
n
(1)
where Fo(Ti) is the observed probability of failure as a function of Ti, which is the failure time of specimen i, and n is the total number of samples. 3. Select a suitable probability distribution from potential distributions such as exponential, normal or Weibuil distributions. 4. Determine the initial values of the parameters for the selected distribution. For example, for a three-parameter Weibull distribution, the initial values of the location, scale and shape parameters should be determined. This can be done by a graphical analysis. 5. Obtain a set of optimized parameters for the selected distribution by minimizing the maximum absolute difference between the observed probability of failure and the
expected probability of failure
Min
Max
i = 1,2 . . . . .
F~(Tg) - F o ( T i )
(2)
n
where Fe(Ti) is the expected probability of failure, which can be calculated from the selected distribution. Eqn (2) can be solved by the minimax algorithm introduced in Ref. 16. Function minimax in Matlab Toolbox 'Optimization' can also be used for this purpose. 6. Compare the maximum absolute difference between the observed probability of failure and the expected probability of failure with the allowable critical value, which is corresponding to a specified significance level. If the maximum absolute difference is less than the allowable critical value, this set of parameters is accepted as the optimized one for the selected distribution. Otherwise, the selected distribution may not be adequate to represent the data of interest. Other distributions may be selected to repeat the process from step 3 to step 6 until an acceptable distribution is found.
3 P A R A M E T E R E S T I M A T I O N FOR THREEP A R A M E T E R W E I B U L L DISTRIBUTIONS The above-proposed optimization method of parameter estimation can be used to estimate the population parameters for different probability distributions. Because of the difficulty of the parameter estimation for three-parameter Weibull distributions, the proposed optimization method is illustrated by the following three examples, in which three-parameter Weibull distributions are selected to represent the probability distributions. The first one is an example of small sample size. The second is an example of a larger sample size. The third is an example of grouped data.
3.1 Example 1 Six specimens were used to measure the ultimate strength of a material. The test results are listed in Table 1 in ascending order 17. If a three-parameter Weibull distribution is selected as a candidate distribution, the population parameters can be estimated by the proposed optimization method as follows. Since the data have been arranged in ascending order, the observed cumulative probability can be calculated according to the median rank in eqn (1). The results are listed in column 3 of Table 1. The cumulative distribution function (cdf) of the three-parameter Weibull distribution is
F(x)=I-exp
\
xa / J
(3)
where x0, xa and 13 are the location, scale and shape
249
P a r a m e t e r estimation
parameters of the three-parameter Weibull distribution, respectively. Substituting eqn (3) into eqn (2) yields
1
. . . .
i
,
'
T
'
'
'
0.8 Min
Max
1-expl-(Si-x°')~l
i=1,2 ......
L
"
n
FO~Si~
"x
~
xo k /j
0.6
(4) 0 . 4
Function m i n i m a x in Matlab is used to conduct the optimization process in the eqn (4). Eight different sets of values for x0, xa and /3 are used as the initial values to start the iteration processes for minimizing the maximum absolute difference between the observed cumulative probability and the expected cumulative probability. The iteration processes all converge to the point
j 0
/
,
300
,
I
,
400
,
,
• _Obseryed. --Optimization -- -- - Least Squares - - - - MatchingMoments ~
I
. . . .
500
I
600
i
i
J
i
~
i
L
700
i
i
i
i
800
i
i
900
i
i
b
1000
MaterialUltimateStrength(MPa)
Fig. 1. A comparison of the observed cumulative probability and the expected cumulative probability obtained by different methods.
x0 ----304.84 xa : 279.55 /3 = 2.1208. This shows the existence of an optimal solution. At the optimal point, the maximum absolute difference between the observed cumulative probability and the expected cumulative probability is 0.0043. From Ref. l, the allowable critical value corresponding to a significance level of 0.05 is 0.521. It is obvious that the three-parameter Weibull distribution is a good representation of the test results. The cdf of the ultimate strength of the material is F(S)= 1-
.,fl',~ .~¢"
0.2
F (
exp[-
\
279--~ - ]
(5)
J"
The parameters of the Weibull distribution can also be estimated by the usual methods such as the least-squares method, the matching moments method, and the maximum likelihood method. The estimated parameters by the leastsquares method and the matching moments method are listed in Table 2. When the maximum likelihood method is used to estimate the parameters of the Weibull distribution, there is no convergent solution. It can be seen from Table 2 that the proposed optimization method gives the minimum value for the maximum absolute difference
between the observed cumulative probability and the expected cumulative probability. Fig. 1 shows a comparison of the observed cumulative probability and the expected cumulative distribution curves obtained by the proposed optimization method, the leastsquares method, and the matching moments method. The expected cumulative distribution curves obtained by the proposed optimization method and the least-squares method are overlapping each other and they are better representations of the observed cumulative probability than that obtained by the matching moments method even though the estimated parameters by the matching moments method are acceptable. 3.2 E x a m p l e 2
Twenty specimens were used to measure the probability distribution of fatigue life of a material under a specified stress level. The obtained fatigue lives from fatigue tests are listed in Table 3 in ascending order 17. The observed cumulative probability calculated according to the median rank in eqn (1) are listed in column 3 of Table 3. If a three-parameter Weibull distribution is selected as a
Table 1. Test results of the material ultimate strength
Specimen No. i 1 2 3 4 5 6
Ultimate strength of specimen Si (MPa)
Observed cumulative probability Fo(S1)
Expected cumulative probability Fe(Si)
Absolute difference
408.0 464.0 514.0 566.0 626.0 708.0
0.1094 0.2656 0.4219 0.5781 0.7344 0.8906
0.1137 0.2613 0.4176 0.5792 0.7387 0.8863
0.0043 0.0043 0.0043 0.0011 0.0043 0.0043
IFe(Si) - Fo(Si)I
Table 2. Estimated parameters by different methods
Method
Xo
x~
~
MaxlF~ (Si) - Fo(Si)I
Optimization Least-squares Matching moments
304.84 309.91 299.92
279.55 274.04 279.44
2.1208 2.0681 2.4142
0.0043 0.0054 0.0314
Jing Ling, Jwo Pan
250
candidate distribution, the parameters of the three-parameter Weibull distribution can be estimated by the proposed optimization method similarly as in Example 1. The obtained results are listed in row 2 of Table 4. The parameters of the Weibull distribution can also be estimated by the least-squares method, the matching moments method and the maximum likelihood method. The estimated parameters are, respectively, listed in rows 3, 4 and 5 of Table 4. It can be seen that the proposed optimization method gives the minimum value for the maximum absolute difference between the observed probability of failure and the expected probability of failure. From Ref. l, the allowable critical value corresponding to significance level of 0.05 is 0.294. It shows that the three-parameter Weibull distribution can adequately represent the test results. The cdf obtained by the optimization method is T _ - 196.98" ]
F(T)=l-exp[-(
398.06
j
3.0413]. j
(6)
Fig. 2 shows a comparison of the observed probability of failure and the curves of the expected probability of failure obtained by the different methods. It can be seen that the cdf curve obtained by the optimization method is a very good representation of the observed failure probability.
1
,
,
,
,
i
. . . .
q
. . . .
1
. . . .
i
. . . .
r
.
,
~ ,
,
~-
,
0.8
,/
06 e~
_e o.4
E_ Sr 0.2
- -- --
jr' ~,~' r
0
- -
200
""-~-"~'"
300
''' 400
Least Squares
- Matching
Moments
----Maximum
~ a . ~ / * ' * '
....
Optimization -
J'
'''
. . . . .
Likelihood '
500 600 700 800 FatigueLife(1000 Cycles)
'''
~,,,
900
In order to determine the probability distribution of the failure life of an electronic component, 115 samples were tested to failure. The test results are listed in Table 5 as grouped data i. If a three-parameter Weibull distribution is used to represent the test data, its parameters can also be determined by the proposed optimization method.
Since the cumulative failures at the middle of each group can be approximated by i-I Hi U,= y , k + 7,
(7)
k = l
the observed failure probability at the middle of each group can be calculated by (8)
The reason for using T / - 250 as the time in eqn (8) is that the observed failure probability is taken at a group interval of 500 h. The calculated results are listed in column 5 of Table 5. Substituting the observed failure probability at the middle of each group and the cumulative distribution function of three-parameter Weibull distribution eqn (3) into
Table 3. Results of fatigue tests under a specified stress level
Specimen No. i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
1000
Fig. 2. A comparison of the observed probability of failure and the expected probability of failure obtained by different methods.
N i - 0.3 Fo(Ti - 250) -- 115 + 0
3.3 E x a m p l e 3
,
Fatigue life of specimen Ti( x 103 cycles)
Observed failure probability Fo(Ti)
Expected failure probability Fe(Ti)
Absolute difference
370 380 420 430 450 480 480 500 510 540 550 590 600 610 630 630 670 730 740 840
0.0343 0.0833 0.1324 0.1814 0.2304 0.2794 0.3284 0.3775 0.4265 0.4755 0.5245 0.5735 0.6225 0.6716 0.7206 0.7696 0.8186 0.8676 0.9167 0.9657
0.0383 0.0792 0.1146 0.1785 0.2273 0.2799 0.3074 0.3639 0.4216 0.4793 0.5079 0.5640 0.6435 0.6684 0.7156 0.7587 0.7976 0.8886 0.9295 0.9718
0.0420 0.0065 0.0254 0.0032 0.0076 0.0190 0.0300 0.0240 0.0444 0.0049 0.0240 0.0444 0.0235 0.0017 0.0046 0.0444 0.0031 0.0444 0.0069 0.0207
IFe(Ti)- Fo(Ti)L
Parameter
251
estimation
Table 4. Estimated parameters by different methods Method
Xo
x~
fl
MaxlFe(Ti) - Fo(Ti)I
Optimization Least Squares Matching Moments Maximum Likelihood
196.98 299.65 300.99 349.87
398.06 293.18 289.64 231.83
3.0413 2.0502 2.1380 1.6891
0.0444 0.0513 0.0569 0.0805
eqn (2), we have
'
i = 1,2,...,20
failure and the expected p d f obtained by the proposed optimization method. It can be seen that the component failure can be characterized by the three-parameter Weibull distribution with the pdf in eqn (10).
'
Xa
(9)
- Fo(Tg - 2 5 0 ) .
The parameters of the Weibull distribution can be determined by solving eqn (9). The estimated parameters are x0 = 2673.3 xa = 3856.9 /3 = 1.8469
4 P A R A M E T E R E S T I M A T I O N FOR THE SUM OF TWO WEIBULL DISTRIBUTIONS W h e n two subpopulations are mixed into one overall population, the sum model is referred to as the one where the cdf o f the overall population F ( x ) is characterized by the sum o f the cdfs of the two sub populations FI (x) and Fz(x). Now we consider three-parameter Weibull distributions as probability distributions, the cdfs of the two subpopulations are
exp[
Therefore, the probability density function (pdf) of the component life is 1.8469 ( T - - 2673.3x] °8469 f(T)--
3856.9\
3856.9
and
]
exp
× e x p E - ( T - 267.3)1.8469]
\ 3~K~ /
j.
(lO)
Fig. 3 shows a comparison o f the histogram of component
E_(x_x0qq \ ~ /
J
(12)
where x01, Xal and /31 are the location, scale and shape parameters of the first subpopulation, and x02, Xa2 and r2
Table 5. Grouped failure data and the associated median ranks Group number i
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
Time at group end point Ti (h) 3000 3500 4000 4500 5000 5500 6000 6500 7000 7500 8000 8500 9000 9500 10000 10500 11 000 11500 12 000 12 500
Failures in each group ni 2 5 8 10 12 13 12 11 10 8 7 5 4 3 2 1 0 1 0 1
Cumulative failures at the middle of each group Ni
Observed failure probability at the middle of each group Fo(Ti - 250)
1.0 4.5 11.0 20.0 31.0 43.5 56.0 67.5 78.0 87.0 94.5 100.5 105.0 108.5 111.0 112.5 113.0 113.5 114.0 114.5
0.0061 0.0364 0.0927 0.1707 0.2660 0.3744 0.4827 0.5823 0.6733 0.7513 0.8163 0.8683 0.9073 0.9376 0.9593 0.9723 0.9766 0.9809 0.9853 0.9896
252
Jing Ling, J w o P a n
are the location, scale and shape parameters of the second subpopulation. Then the cdf of the overall population is
The observed failure probability can be calculated according to the median rank in eqn (1) and the results are listed in column 3 of Table 6. Eqn (14) and the observed failure probability can be substituted into eqn (2) and the optimized parameters can be obtained by using function m i n i m a x in Matlab as follows
(13)
F(x) = rF l (x) + (1 -- r)F2(x)
where r and (1 - r) are the mixing weights for the two subpopulations, respectively. Eqn (13) can be rewritten as
0
X ~X01
(14)
X01 ~ X ~ X02
F(x) =
r{1-exp r{l-
I-
( ~ ) ~ ]
exp[-(~)~'l}+(1-r){1-
}
~:1
exp[-(~)
It can be seen from eqn (14) that there are seven parameters to be estimated to obtain the cdf of the overall population. They can be estimated by substituting eqn (14) into eqn (2) and solving eqn (2) using function m i n i m a x in Matlab. Here is a numerical example for demonstration of the proposed method.
}
X ~ X02
x01 = 3.4203 xal = 7.6764
/31 = 2.4547 r = 0.5767
4.1 E x a m p l e 4
x02 = 6.0717 The fatigue lives of 25 specimens of two different types are listed in Table 6 in ascending order. The distribution parameters can be estimated by the optimization method if the sum model is used to represent the data.
x,2 = 31.5442 /32 = 3.6680.
Table 6. Data for the sum model
Specimen number i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
Observed failure probability Fo(xi)
Expected failure probability Fe(xi)
Absolute difference
xi ( x 105 cycles)
Fatigue life of specimen 5.2426 6.3014 7.4485 8.1348 8.4369 9.0949 9.9299 10.258 10.997 11.367 11.695 12.631 13.714 15.504 19.891 23.785 26.467 28.345 31.098 34.850 35.689 37.225 41.059 44.892 49.839
0.0276 0.0669 0.1063 0.1457 0.1850 0.2244 0.2638 0.3031 0.3425 0.3819 0.4213 0.4606 0.5000 0.5394 0.5787 0.6181 0.6575 0.6969 0.7362 0.7756 0.8150 0.8543 0.8937 0.9331 0.9724
0.0167 0.0498 0.1071 0.1504 0.1711 0.2186 0.2809 0.3053 0.3582 0.3831 0.4042 0.4573 0.5051 0.5543 0.5958 0.6247 0.6541 0.6798 0.7240 0.7927 0.8086 0.8372 0.9019 0.9502 0.9848
0.0109 0.0171 0.0008 0.0047 0.0139 0.0058 0.0171 0.0022 0.0157 0.0012 0.0171 0.0033 0.0051 0.0149 0.0171 0.0066 0.0034 0.0171 0.0122 0.0171 0.0064 0.0171 0.0082 0.0171 0.0124
IFe(xi)- Fo(xi)l
Parameter estimation The corresponding maximum absolute difference between the observed probability of failure and the expected probability of failure is 0.0171, and is far less than the allowable critical value 0.270 corresponding to a significance level of 0.05. Substituting the values of the parameters into eqn (14), we obtain the cdf of specimen fatigue life as
253
Fig. 6 is a flowchart for selection of population distribution and the corresponding parameter estimation. At first, several potential distributions are selected based on a preliminary analysis of the sample data. According to the selected distributions, the population parameters are estimated from the sample data by the proposed optimization
x < 3.4203 {
F(x) =
1-
0.5767
0.5767
{
1 - exp
[ " ~X~ -- 3"4203"~ - - 2"4547exp
\
[ _ (x -- 3,4203"~2"4547]~ k. ~
J
7.6764 J
J J +0"4233
{
< x<
3.4203
. J
[ (x--6-0717"~3"66801}
1 - exp - \
~
6.0717
x__>6.0717
j
(15)
Fig. 4 shows a comparison of the theoretical cdf of failure and the observed failure probability. It can be seen from Fig. 4 that the sum model of the two three-parameter Weibull distributions is a good representation of the fatigue lives of specimens. The pdf of the fatigue life of the specimen is illustrated in Fig. 5. It is obvious that the specimens are from two different populations.
5 A N E W P R O C E D U R E FOR S E L E C T I O N OF POPULATION DISTRIBUTION AND PARAMETER ESTIMATION Traditionally, when we analyse the data, a potential distribution is selected, and its parameters are estimated based on the data. If the goodness-of-fit test shows that the selected population distribution with estimated parameters represents the data adequately, we stop here and accept this distribution. In this procedure, we are not sure if we get the best distribution or not. To get the best distribution with optimized parameters, a new procedure is proposed in the following.
method. For each distribution, there is a maximum absolute difference between the observed probability of failure and the expected probability of failure. The distribution with the minimum absolute difference between the observed probability of failure and the expected probability of failure is selected. Then, a Kolmogorov-Smirnov goodness-of-fit test is applied to this distribution to make sure that the distribution fits the sample data well. If the maximum absolute difference between the observed probability of failure and the expected probability of failure associated with this distribution is less than the allowable critical value corresponding to a specified significance level, this distribution is accepted as the best distribution for the sample data. The corresponding estimates of the population parameters are the optimized parameters. 5.1 Example 5 To demonstrate the procedure, we analyse the data in Table 6 according to the above procedure. First, a two-parameter exponential distribution, a normal distribution, a log normal
0.0003 . . . . . . . . . . . . . . . . . . . . . . .
I
0.00025 e.
"~
0.0002
0.00015
i~" 0.0001 51057 i
° 2
ComponentLife(1000Hours)
J
Fig. 3. A comparison of the histogram of component failure and the expected pdf obtained by the optimization method.
T
. . . .
i
. . . .
f
. . . .
1
.
.
.
.
.
.
.
.
0.8
.= i
0.6
•~
0.4
r.)
0.2 0
14
. . . .
ff¢¢ 0 10 ,
~.lr,
,
I
served Expected 20 30 40 50 60 SpecimenFatigueLife(100,000Cycles)
. . . .
I
,
i
h
I
b
i
i
~
i
b
~
,
I
,
,
,
,
Fig. 4. A comparison of the theoretical cdf of failure and the observed failure probability.
Jing Ling, Jwo Pan
254 0.08
....
i ....
~ ....
i ....
I ....
1
i ....
....
i ....
, ,.,
, .-I/'...." .'~""F ~:'k~''~'~''" ~..~. ..... ~-......
0.07 ~. ~
0.06 rr
0.05
~ ¢.".
0.6
"~
0.4
0.03
~'r/'
E
0.02
~ O.2
Lt
0 0
i
i ,
I0
20 30 40 50 Specimen Fatigue Life (100,000 Cycles)
(16)
where 3, and 3' are the failure rate and location parameter, and T is the failure time. The pdf of a normal distribution is (T2~o---2~)21
(17)
where # and a are the mean and standard deviation of
I
Sample Data
I
Potential Distributions Exponential Normal Lognormal Weibull Mixed-Weibull etc.
I °p'imize ]
@ Parameter Estimation
Best Distribution with Optimized Paramters
~
i , , 10
0
, k 20
.
.
-Sum of Two Weibulls . . . L , , , 30 40 50 60
Fig. 7. A comparison of the theoretical cdfs of failure and the observed failure probability for different distributions.
distribution, a three-parameter Weibull distribution and the sum model of two three-parameter Weibull distributions are selected as candidate distributions. The cdf of an exponential distribution is
f(T)-- v / ~ e x p [
'°~'
- - - W~bul]
Specimen Fatigue Life (100,000 Cycles)
Fig. 5. The PDF of specimen fatigue life.
F(T) = 1 - e - x(r- v)
• Observed ........ Exponential . . . . . Normal . . . . . Lognormal
~
/~t]
0 60
'
J .~'
0.01
I
~
e~
0.04
d:
0.8
No Suitable Distribution J
failure times. The pdf of a log normal distribution is 1
f ( T ) - T v / ~ o ~' exp
[
(logT~/~')21 -2~~ A
(18)
where/~' and o' are the mean and standard deviation of the logarithms of the failure times. The cdf of a three-parameter Weibull distribution is expressed in eqn (3). The sum model of two three-parameter Weibull distributions is expressed in eqn (14). The parameters in these distributions can be estimated by the proposed optimization method. The estimated parameters and the corresponding maximum absolute difference between the observed probability of failure and the expected probability of failure associated with each distribution are listed in Table 7. It can be seen from Table 7 that the maximum absolute difference between the observed probability of failure and the expected probability of failure associated with the sum model of two three-parameter Weibull distributions is the smallest among these different distributions. For the sample data in Table 6, the allowable critical value corresponding to significance level of 0.05 is 0.270, and far greater than 0.0171 which is the maximum absolute difference between the observed probability of failure and the expected probability of failure associated with the sum model of two three-parameter Weibull distributions. The cdf of the sum model of two three-parameter Weibull distributions is expressed in eqn (15). Fig. 7 shows a comparison of the theoretical cdfs of failure and the observed failure probability for different distributions. It is obvious that the sum model of two threeparameter Weibull distributions is the best representation of the sample data among these distributions.
6 CONCLUSIONS J
Fig. 6. A flowchart for selection of population distribution and parameter estimation.
An optimization method for parameter estimation is proposed. The method is illustrated by estimating the population parameters for three-parameter Weibull distributions. A comparison between the proposed method and the usual
Parameter estimation
255
Table 7. Results of the parameter estimation in example 5 Distribution
Parameter
MaxlFe(Ti) - Fo(Ti)l
Exponential Normal Log normal Weibull Sum model
X = 0.0583 3" = 4.1006 /x = 17.9416 tr = 12.7525 ~' = 6.2043 a' = 0.3594 x0 = 5.1552xa = 16.1626/~ = 0.8665 3'1 = 3.42037h = 7.6764/~1 = 2.4547 p = 0.5767 3'2 = 6.0717'/2 = 31.5442/32 = 3.6680
0.0710 0.1308 0.0724 0.0619 0.0171
methods such as the least-squares method, the matching moments method and the m a x i m u m likelihood method shows that more reasonable estimates of the parameters are given by the proposed optimization method. In addition, the proposed method can be used to estimate the parameters for the sum of two three-parameter Weibull distributions. A new procedure for selecting population distribution and parameter estimation is also presented.
8.
9.
10.
REFERENCES 1. Kececioglu, D. B., Reliability and Life Testing Handbook, Vol. 1. PTR Prentice Hall, 1993. 2. Li, Yong-Ming A general linear-regression analysis applied to the three-parameter Weibull distribution. IEEE Trans. Reliability, 1994, R-43, 255-263. 3. Weibull, W., Fatigue Testing and Analysis of Results. Pergamon Press, Oxford, 1961. 4. Newby, M. Properties of moment estimators for the threeparameter Weibull distribution. 1EEE Trans. Reliability, 1984, R-33, 192-195. 5. Cran, G. W. Moment estimators for the three-parameter Weibull distribution. IEEE Trans. Reliability, 1988, R-37, 360363. 6. Rockette, H., Antle, C. E. & Klimko, L. A. Maximum likelihood estimation with the Weibull Distribution. J. Am. Statist. Assoc., 1974, 69, 246-249. 7. Harter, H. L. & Moore, A. H. Maximum likelihood estimation
11.
12.
13. 14.
15.
16.
17.
of the parameters of gamma and Weibull populations from complete and from censored samples. Technometrics, 1965, 7, 639-643. Sinha, S. K. & Sloan, J. A. Bayes estimation of the parameters and reliability function of the three-parameter Weibull distribution. IEEE Trans. Reliability, 1988, R-37, 364-369. Hobbs, J. R., Moore, A. H. & Miller, R. M. Minimum-distance estimation of the parameters of the three-parameters Weibull distribution. IEEE Trans. Reliability, 1985, R-34, 495 -496. Adatia, A. & Chan, L. K. Robust estimators of the threeparameter Weibull distribution. 1EEE Trans. Reliability, 1985, R-34, 347-351. Jiang, S. & Kececioglu, D, B. Graphical representation of two mixed Weibull distributions. IEEE Trans. Reliability, 1992, 41, 241-247. Kao, J. H. K. A graphical estimation of mixed Weibull parameters in life testing of electron tubes. Technometrics, 1959, 1, 389-407. Cran, H. Graphical estimation methods for Weibull distributions. Microelectronics and Reliability, 1976, 15, 47-52. Jensen, F. and Petersen, An Engineering Approach to the Design and Analysis of Burn-In Procedures. John Wiley and Sons, Chichester, 1982. Jiang, R. & Murthy, D. N. P. Reliability modeling involving two Weibull distributions. Reliability Engineering and System Safety, 1995, 47, 187-198. Madsen, K. & Schjaer-Jacobsen, H. Algorithms for worstcase tolerance optimization. IEEE Trans. Circuits and Systems, 1979, CAS-26, 775-783. Gao, Z. T., Applied Statistics to Fatigue. Press of Industry of National Defense, Beijing, 1986.
ELSEVIER
© 1998 Elsevier Science Limited All rights reserved. Printed in Northern Ireland 0951-8320/98/$19.00
Technical note A new method for selection of population distribution and parameter estimation Jing Ling & Jwo Pan* Center for Automotive Structural Durability Simulation, Department of Mechanical Engineering and Applied Mechanics, The University of Michigan, Ann Arbor, MI 48109, USA
(Received 18 March 1997; accepted 3 December 1997)
An optimization method for parameter estimation is presented with the KolmogorovSmirnov distance used as the objective. A step-by-step implementation procedure is given. The method is demonstrated by estimating the parameters for three-parameter Weibull distributions from three different samples (with different sample sizes). A comparison of the proposed method and the usual methods such as the least-squares method, the matching moments method and the maximum likelihood method shows that more reasonable estimates of the parameters are given by the proposed optimization method. Then, the proposed method is successfully extended to estimate the parameters for the sum of two three-parameter Weibull distributions. Based on these findings, a new procedure for selection of population distribution and parameter estimation is presented. © 1998 Elsevier Science Limited.
1 INTRODUCTION
parameter estimation is proposed with the K o l m o g o r o v Smirnov distance used as the objective. This is in contrast to Ref. 9, in which minimizing the Kolmogorov-Smirnov distance is only used to estimate the location parameter of three-parameter Weibull distributions. Sometimes, two subpopulations may be mixed into one overall population, a sum model of the two distributions may be a suitable model to characterize the overall population. When three-parameter Weibull distributions are used, it is difficult to estimate the parameters for the model. Jiang and Kececioglu i i studied the behavior of a two-Weibull mixture on the Weibull probability papers and evaluated two graphical parameter estimation methods: the K a o Cran method 12,13 and the Jensen-Petersen method 14, for this kind of overall distribution. Jiang and Murthy 15 did similar research. However, all of these studies are limited to the mixed models of two-parameter Weibull distributions because of the limitation of the probability plotting method. More refined statistical methods are needed to estimate the parameters for the sum of three-parameter Weibull distributions ~5. Our proposed optimization method for parameter estimation is also applicable to the sum of Weibull distributions. In contrast to Refs l 1,15, the proposed optimization method can be used to estimate the parameters for the sum of three-parameter Weibull distributions. Based on
Parameter estimation plays a very important role in the evaluation of reliability. Much research has been done to find more accurate, practical methods for parameter estimation. The probability plotting method l is easiest and straightforward, but less precise than analytical methods. The probability plotting method is usually used to obtain the initial parameter estimation for analytical methods. The most commonly used analytical methods for parameter estimation are the least-squares method l, or the regression method 2 the matching moments method 3-5, and the maximum likelihood method 6,7. There are also other methods presented in literature, for example, the Bayes estimation method 8, the minimum-distance estimation method 9 and the robust estimation method 10 No matter what kind of parameter estimation method is used, the selected distribution with the estimated parameters must be tested to determine whether the distribution adequately and accurately represents the existing data of interest. The best distribution will be the one with the minimum difference test statistic value or the maximum likelihood test statistic value. In this paper, an optimization method for *To whom correspondence should be addressed. 247
Jing Ling, Jwo Pan
248
the proposed optimization method for parameter estimation, a new procedure for selection of population distribution is introduced. This procedure can easily be used to obtain the best population distribution with optimized parameters.
2 AN O P T I M I Z A T I O N M E T H O D OF P A R A M E T E R ESTIMATION No matter what kind of parameter estimation method is used, the selected distribution with the estimated parameters must be tested to determine whether it represents the existing data well. This analysis of the fitness of a selected distribution to the data may be accomplished by the Kolmogorov-Smirnov (K-S) test. The principle and the detailed implementation procedure of the K - S test are introduced in Ref. i. The K - S test is used to determine whether the maximum absolute difference between the observed probability of failure and the expected probability of failure is less than the allowable critical value corresponding to a specified significance level. If the maximum absolute difference is less than the allowable critical value, the selected distribution is accepted as a reasonable distribution. If the maximum absolute difference is equal to or greater than the allowable critical value, the selected distribution is rejected. Furthermore, the allowable critical value decreases with the increase of the significance level. Therefore, when a set of parameters of the selected distribution can be obtained by minimizing the maximum absolute difference between the observed probability of failure and the expected probability of failure, this set of parameters can be regarded as the optimized one for the selected distribution. The procedure to implement the above idea is as follows: 1. Arrange all data in ascending order, or rank them. 2. Calculate the observed probability of failure. There are three different rank definitions that can be used to calculate the observed probability of failure. The median rank as defined in eqn (1) is used here because of its wide acceptance.
Fo(Tg ) -
i-0.3 i = l, 2, n+0.4 . . . . .
n
(1)
where Fo(Ti) is the observed probability of failure as a function of Ti, which is the failure time of specimen i, and n is the total number of samples. 3. Select a suitable probability distribution from potential distributions such as exponential, normal or Weibuil distributions. 4. Determine the initial values of the parameters for the selected distribution. For example, for a three-parameter Weibull distribution, the initial values of the location, scale and shape parameters should be determined. This can be done by a graphical analysis. 5. Obtain a set of optimized parameters for the selected distribution by minimizing the maximum absolute difference between the observed probability of failure and the
expected probability of failure
Min
Max
i = 1,2 . . . . .
F~(Tg) - F o ( T i )
(2)
n
where Fe(Ti) is the expected probability of failure, which can be calculated from the selected distribution. Eqn (2) can be solved by the minimax algorithm introduced in Ref. 16. Function minimax in Matlab Toolbox 'Optimization' can also be used for this purpose. 6. Compare the maximum absolute difference between the observed probability of failure and the expected probability of failure with the allowable critical value, which is corresponding to a specified significance level. If the maximum absolute difference is less than the allowable critical value, this set of parameters is accepted as the optimized one for the selected distribution. Otherwise, the selected distribution may not be adequate to represent the data of interest. Other distributions may be selected to repeat the process from step 3 to step 6 until an acceptable distribution is found.
3 P A R A M E T E R E S T I M A T I O N FOR THREEP A R A M E T E R W E I B U L L DISTRIBUTIONS The above-proposed optimization method of parameter estimation can be used to estimate the population parameters for different probability distributions. Because of the difficulty of the parameter estimation for three-parameter Weibull distributions, the proposed optimization method is illustrated by the following three examples, in which three-parameter Weibull distributions are selected to represent the probability distributions. The first one is an example of small sample size. The second is an example of a larger sample size. The third is an example of grouped data.
3.1 Example 1 Six specimens were used to measure the ultimate strength of a material. The test results are listed in Table 1 in ascending order 17. If a three-parameter Weibull distribution is selected as a candidate distribution, the population parameters can be estimated by the proposed optimization method as follows. Since the data have been arranged in ascending order, the observed cumulative probability can be calculated according to the median rank in eqn (1). The results are listed in column 3 of Table 1. The cumulative distribution function (cdf) of the three-parameter Weibull distribution is
F(x)=I-exp
\
xa / J
(3)
where x0, xa and 13 are the location, scale and shape
249
P a r a m e t e r estimation
parameters of the three-parameter Weibull distribution, respectively. Substituting eqn (3) into eqn (2) yields
1
. . . .
i
,
'
T
'
'
'
0.8 Min
Max
1-expl-(Si-x°')~l
i=1,2 ......
L
"
n
FO~Si~
"x
~
xo k /j
0.6
(4) 0 . 4
Function m i n i m a x in Matlab is used to conduct the optimization process in the eqn (4). Eight different sets of values for x0, xa and /3 are used as the initial values to start the iteration processes for minimizing the maximum absolute difference between the observed cumulative probability and the expected cumulative probability. The iteration processes all converge to the point
j 0
/
,
300
,
I
,
400
,
,
• _Obseryed. --Optimization -- -- - Least Squares - - - - MatchingMoments ~
I
. . . .
500
I
600
i
i
J
i
~
i
L
700
i
i
i
i
800
i
i
900
i
i
b
1000
MaterialUltimateStrength(MPa)
Fig. 1. A comparison of the observed cumulative probability and the expected cumulative probability obtained by different methods.
x0 ----304.84 xa : 279.55 /3 = 2.1208. This shows the existence of an optimal solution. At the optimal point, the maximum absolute difference between the observed cumulative probability and the expected cumulative probability is 0.0043. From Ref. l, the allowable critical value corresponding to a significance level of 0.05 is 0.521. It is obvious that the three-parameter Weibull distribution is a good representation of the test results. The cdf of the ultimate strength of the material is F(S)= 1-
.,fl',~ .~¢"
0.2
F (
exp[-
\
279--~ - ]
(5)
J"
The parameters of the Weibull distribution can also be estimated by the usual methods such as the least-squares method, the matching moments method, and the maximum likelihood method. The estimated parameters by the leastsquares method and the matching moments method are listed in Table 2. When the maximum likelihood method is used to estimate the parameters of the Weibull distribution, there is no convergent solution. It can be seen from Table 2 that the proposed optimization method gives the minimum value for the maximum absolute difference
between the observed cumulative probability and the expected cumulative probability. Fig. 1 shows a comparison of the observed cumulative probability and the expected cumulative distribution curves obtained by the proposed optimization method, the leastsquares method, and the matching moments method. The expected cumulative distribution curves obtained by the proposed optimization method and the least-squares method are overlapping each other and they are better representations of the observed cumulative probability than that obtained by the matching moments method even though the estimated parameters by the matching moments method are acceptable. 3.2 E x a m p l e 2
Twenty specimens were used to measure the probability distribution of fatigue life of a material under a specified stress level. The obtained fatigue lives from fatigue tests are listed in Table 3 in ascending order 17. The observed cumulative probability calculated according to the median rank in eqn (1) are listed in column 3 of Table 3. If a three-parameter Weibull distribution is selected as a
Table 1. Test results of the material ultimate strength
Specimen No. i 1 2 3 4 5 6
Ultimate strength of specimen Si (MPa)
Observed cumulative probability Fo(S1)
Expected cumulative probability Fe(Si)
Absolute difference
408.0 464.0 514.0 566.0 626.0 708.0
0.1094 0.2656 0.4219 0.5781 0.7344 0.8906
0.1137 0.2613 0.4176 0.5792 0.7387 0.8863
0.0043 0.0043 0.0043 0.0011 0.0043 0.0043
IFe(Si) - Fo(Si)I
Table 2. Estimated parameters by different methods
Method
Xo
x~
~
MaxlF~ (Si) - Fo(Si)I
Optimization Least-squares Matching moments
304.84 309.91 299.92
279.55 274.04 279.44
2.1208 2.0681 2.4142
0.0043 0.0054 0.0314
Jing Ling, Jwo Pan
250
candidate distribution, the parameters of the three-parameter Weibull distribution can be estimated by the proposed optimization method similarly as in Example 1. The obtained results are listed in row 2 of Table 4. The parameters of the Weibull distribution can also be estimated by the least-squares method, the matching moments method and the maximum likelihood method. The estimated parameters are, respectively, listed in rows 3, 4 and 5 of Table 4. It can be seen that the proposed optimization method gives the minimum value for the maximum absolute difference between the observed probability of failure and the expected probability of failure. From Ref. l, the allowable critical value corresponding to significance level of 0.05 is 0.294. It shows that the three-parameter Weibull distribution can adequately represent the test results. The cdf obtained by the optimization method is T _ - 196.98" ]
F(T)=l-exp[-(
398.06
j
3.0413]. j
(6)
Fig. 2 shows a comparison of the observed probability of failure and the curves of the expected probability of failure obtained by the different methods. It can be seen that the cdf curve obtained by the optimization method is a very good representation of the observed failure probability.
1
,
,
,
,
i
. . . .
q
. . . .
1
. . . .
i
. . . .
r
.
,
~ ,
,
~-
,
0.8
,/
06 e~
_e o.4
E_ Sr 0.2
- -- --
jr' ~,~' r
0
- -
200
""-~-"~'"
300
''' 400
Least Squares
- Matching
Moments
----Maximum
~ a . ~ / * ' * '
....
Optimization -
J'
'''
. . . . .
Likelihood '
500 600 700 800 FatigueLife(1000 Cycles)
'''
~,,,
900
In order to determine the probability distribution of the failure life of an electronic component, 115 samples were tested to failure. The test results are listed in Table 5 as grouped data i. If a three-parameter Weibull distribution is used to represent the test data, its parameters can also be determined by the proposed optimization method.
Since the cumulative failures at the middle of each group can be approximated by i-I Hi U,= y , k + 7,
(7)
k = l
the observed failure probability at the middle of each group can be calculated by (8)
The reason for using T / - 250 as the time in eqn (8) is that the observed failure probability is taken at a group interval of 500 h. The calculated results are listed in column 5 of Table 5. Substituting the observed failure probability at the middle of each group and the cumulative distribution function of three-parameter Weibull distribution eqn (3) into
Table 3. Results of fatigue tests under a specified stress level
Specimen No. i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
1000
Fig. 2. A comparison of the observed probability of failure and the expected probability of failure obtained by different methods.
N i - 0.3 Fo(Ti - 250) -- 115 + 0
3.3 E x a m p l e 3
,
Fatigue life of specimen Ti( x 103 cycles)
Observed failure probability Fo(Ti)
Expected failure probability Fe(Ti)
Absolute difference
370 380 420 430 450 480 480 500 510 540 550 590 600 610 630 630 670 730 740 840
0.0343 0.0833 0.1324 0.1814 0.2304 0.2794 0.3284 0.3775 0.4265 0.4755 0.5245 0.5735 0.6225 0.6716 0.7206 0.7696 0.8186 0.8676 0.9167 0.9657
0.0383 0.0792 0.1146 0.1785 0.2273 0.2799 0.3074 0.3639 0.4216 0.4793 0.5079 0.5640 0.6435 0.6684 0.7156 0.7587 0.7976 0.8886 0.9295 0.9718
0.0420 0.0065 0.0254 0.0032 0.0076 0.0190 0.0300 0.0240 0.0444 0.0049 0.0240 0.0444 0.0235 0.0017 0.0046 0.0444 0.0031 0.0444 0.0069 0.0207
IFe(Ti)- Fo(Ti)L
Parameter
251
estimation
Table 4. Estimated parameters by different methods Method
Xo
x~
fl
MaxlFe(Ti) - Fo(Ti)I
Optimization Least Squares Matching Moments Maximum Likelihood
196.98 299.65 300.99 349.87
398.06 293.18 289.64 231.83
3.0413 2.0502 2.1380 1.6891
0.0444 0.0513 0.0569 0.0805
eqn (2), we have
'
i = 1,2,...,20
failure and the expected p d f obtained by the proposed optimization method. It can be seen that the component failure can be characterized by the three-parameter Weibull distribution with the pdf in eqn (10).
'
Xa
(9)
- Fo(Tg - 2 5 0 ) .
The parameters of the Weibull distribution can be determined by solving eqn (9). The estimated parameters are x0 = 2673.3 xa = 3856.9 /3 = 1.8469
4 P A R A M E T E R E S T I M A T I O N FOR THE SUM OF TWO WEIBULL DISTRIBUTIONS W h e n two subpopulations are mixed into one overall population, the sum model is referred to as the one where the cdf o f the overall population F ( x ) is characterized by the sum o f the cdfs of the two sub populations FI (x) and Fz(x). Now we consider three-parameter Weibull distributions as probability distributions, the cdfs of the two subpopulations are
exp[
Therefore, the probability density function (pdf) of the component life is 1.8469 ( T - - 2673.3x] °8469 f(T)--
3856.9\
3856.9
and
]
exp
× e x p E - ( T - 267.3)1.8469]
\ 3~K~ /
j.
(lO)
Fig. 3 shows a comparison o f the histogram of component
E_(x_x0qq \ ~ /
J
(12)
where x01, Xal and /31 are the location, scale and shape parameters of the first subpopulation, and x02, Xa2 and r2
Table 5. Grouped failure data and the associated median ranks Group number i
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
Time at group end point Ti (h) 3000 3500 4000 4500 5000 5500 6000 6500 7000 7500 8000 8500 9000 9500 10000 10500 11 000 11500 12 000 12 500
Failures in each group ni 2 5 8 10 12 13 12 11 10 8 7 5 4 3 2 1 0 1 0 1
Cumulative failures at the middle of each group Ni
Observed failure probability at the middle of each group Fo(Ti - 250)
1.0 4.5 11.0 20.0 31.0 43.5 56.0 67.5 78.0 87.0 94.5 100.5 105.0 108.5 111.0 112.5 113.0 113.5 114.0 114.5
0.0061 0.0364 0.0927 0.1707 0.2660 0.3744 0.4827 0.5823 0.6733 0.7513 0.8163 0.8683 0.9073 0.9376 0.9593 0.9723 0.9766 0.9809 0.9853 0.9896
252
Jing Ling, J w o P a n
are the location, scale and shape parameters of the second subpopulation. Then the cdf of the overall population is
The observed failure probability can be calculated according to the median rank in eqn (1) and the results are listed in column 3 of Table 6. Eqn (14) and the observed failure probability can be substituted into eqn (2) and the optimized parameters can be obtained by using function m i n i m a x in Matlab as follows
(13)
F(x) = rF l (x) + (1 -- r)F2(x)
where r and (1 - r) are the mixing weights for the two subpopulations, respectively. Eqn (13) can be rewritten as
0
X ~X01
(14)
X01 ~ X ~ X02
F(x) =
r{1-exp r{l-
I-
( ~ ) ~ ]
exp[-(~)~'l}+(1-r){1-
}
~:1
exp[-(~)
It can be seen from eqn (14) that there are seven parameters to be estimated to obtain the cdf of the overall population. They can be estimated by substituting eqn (14) into eqn (2) and solving eqn (2) using function m i n i m a x in Matlab. Here is a numerical example for demonstration of the proposed method.
}
X ~ X02
x01 = 3.4203 xal = 7.6764
/31 = 2.4547 r = 0.5767
4.1 E x a m p l e 4
x02 = 6.0717 The fatigue lives of 25 specimens of two different types are listed in Table 6 in ascending order. The distribution parameters can be estimated by the optimization method if the sum model is used to represent the data.
x,2 = 31.5442 /32 = 3.6680.
Table 6. Data for the sum model
Specimen number i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
Observed failure probability Fo(xi)
Expected failure probability Fe(xi)
Absolute difference
xi ( x 105 cycles)
Fatigue life of specimen 5.2426 6.3014 7.4485 8.1348 8.4369 9.0949 9.9299 10.258 10.997 11.367 11.695 12.631 13.714 15.504 19.891 23.785 26.467 28.345 31.098 34.850 35.689 37.225 41.059 44.892 49.839
0.0276 0.0669 0.1063 0.1457 0.1850 0.2244 0.2638 0.3031 0.3425 0.3819 0.4213 0.4606 0.5000 0.5394 0.5787 0.6181 0.6575 0.6969 0.7362 0.7756 0.8150 0.8543 0.8937 0.9331 0.9724
0.0167 0.0498 0.1071 0.1504 0.1711 0.2186 0.2809 0.3053 0.3582 0.3831 0.4042 0.4573 0.5051 0.5543 0.5958 0.6247 0.6541 0.6798 0.7240 0.7927 0.8086 0.8372 0.9019 0.9502 0.9848
0.0109 0.0171 0.0008 0.0047 0.0139 0.0058 0.0171 0.0022 0.0157 0.0012 0.0171 0.0033 0.0051 0.0149 0.0171 0.0066 0.0034 0.0171 0.0122 0.0171 0.0064 0.0171 0.0082 0.0171 0.0124
IFe(xi)- Fo(xi)l
Parameter estimation The corresponding maximum absolute difference between the observed probability of failure and the expected probability of failure is 0.0171, and is far less than the allowable critical value 0.270 corresponding to a significance level of 0.05. Substituting the values of the parameters into eqn (14), we obtain the cdf of specimen fatigue life as
253
Fig. 6 is a flowchart for selection of population distribution and the corresponding parameter estimation. At first, several potential distributions are selected based on a preliminary analysis of the sample data. According to the selected distributions, the population parameters are estimated from the sample data by the proposed optimization
x < 3.4203 {
F(x) =
1-
0.5767
0.5767
{
1 - exp
[ " ~X~ -- 3"4203"~ - - 2"4547exp
\
[ _ (x -- 3,4203"~2"4547]~ k. ~
J
7.6764 J
J J +0"4233
{
< x<
3.4203
. J
[ (x--6-0717"~3"66801}
1 - exp - \
~
6.0717
x__>6.0717
j
(15)
Fig. 4 shows a comparison of the theoretical cdf of failure and the observed failure probability. It can be seen from Fig. 4 that the sum model of the two three-parameter Weibull distributions is a good representation of the fatigue lives of specimens. The pdf of the fatigue life of the specimen is illustrated in Fig. 5. It is obvious that the specimens are from two different populations.
5 A N E W P R O C E D U R E FOR S E L E C T I O N OF POPULATION DISTRIBUTION AND PARAMETER ESTIMATION Traditionally, when we analyse the data, a potential distribution is selected, and its parameters are estimated based on the data. If the goodness-of-fit test shows that the selected population distribution with estimated parameters represents the data adequately, we stop here and accept this distribution. In this procedure, we are not sure if we get the best distribution or not. To get the best distribution with optimized parameters, a new procedure is proposed in the following.
method. For each distribution, there is a maximum absolute difference between the observed probability of failure and the expected probability of failure. The distribution with the minimum absolute difference between the observed probability of failure and the expected probability of failure is selected. Then, a Kolmogorov-Smirnov goodness-of-fit test is applied to this distribution to make sure that the distribution fits the sample data well. If the maximum absolute difference between the observed probability of failure and the expected probability of failure associated with this distribution is less than the allowable critical value corresponding to a specified significance level, this distribution is accepted as the best distribution for the sample data. The corresponding estimates of the population parameters are the optimized parameters. 5.1 Example 5 To demonstrate the procedure, we analyse the data in Table 6 according to the above procedure. First, a two-parameter exponential distribution, a normal distribution, a log normal
0.0003 . . . . . . . . . . . . . . . . . . . . . . .
I
0.00025 e.
"~
0.0002
0.00015
i~" 0.0001 51057 i
° 2
ComponentLife(1000Hours)
J
Fig. 3. A comparison of the histogram of component failure and the expected pdf obtained by the optimization method.
T
. . . .
i
. . . .
f
. . . .
1
.
.
.
.
.
.
.
.
0.8
.= i
0.6
•~
0.4
r.)
0.2 0
14
. . . .
ff¢¢ 0 10 ,
~.lr,
,
I
served Expected 20 30 40 50 60 SpecimenFatigueLife(100,000Cycles)
. . . .
I
,
i
h
I
b
i
i
~
i
b
~
,
I
,
,
,
,
Fig. 4. A comparison of the theoretical cdf of failure and the observed failure probability.
Jing Ling, Jwo Pan
254 0.08
....
i ....
~ ....
i ....
I ....
1
i ....
....
i ....
, ,.,
, .-I/'...." .'~""F ~:'k~''~'~''" ~..~. ..... ~-......
0.07 ~. ~
0.06 rr
0.05
~ ¢.".
0.6
"~
0.4
0.03
~'r/'
E
0.02
~ O.2
Lt
0 0
i
i ,
I0
20 30 40 50 Specimen Fatigue Life (100,000 Cycles)
(16)
where 3, and 3' are the failure rate and location parameter, and T is the failure time. The pdf of a normal distribution is (T2~o---2~)21
(17)
where # and a are the mean and standard deviation of
I
Sample Data
I
Potential Distributions Exponential Normal Lognormal Weibull Mixed-Weibull etc.
I °p'imize ]
@ Parameter Estimation
Best Distribution with Optimized Paramters
~
i , , 10
0
, k 20
.
.
-Sum of Two Weibulls . . . L , , , 30 40 50 60
Fig. 7. A comparison of the theoretical cdfs of failure and the observed failure probability for different distributions.
distribution, a three-parameter Weibull distribution and the sum model of two three-parameter Weibull distributions are selected as candidate distributions. The cdf of an exponential distribution is
f(T)-- v / ~ e x p [
'°~'
- - - W~bul]
Specimen Fatigue Life (100,000 Cycles)
Fig. 5. The PDF of specimen fatigue life.
F(T) = 1 - e - x(r- v)
• Observed ........ Exponential . . . . . Normal . . . . . Lognormal
~
/~t]
0 60
'
J .~'
0.01
I
~
e~
0.04
d:
0.8
No Suitable Distribution J
failure times. The pdf of a log normal distribution is 1
f ( T ) - T v / ~ o ~' exp
[
(logT~/~')21 -2~~ A
(18)
where/~' and o' are the mean and standard deviation of the logarithms of the failure times. The cdf of a three-parameter Weibull distribution is expressed in eqn (3). The sum model of two three-parameter Weibull distributions is expressed in eqn (14). The parameters in these distributions can be estimated by the proposed optimization method. The estimated parameters and the corresponding maximum absolute difference between the observed probability of failure and the expected probability of failure associated with each distribution are listed in Table 7. It can be seen from Table 7 that the maximum absolute difference between the observed probability of failure and the expected probability of failure associated with the sum model of two three-parameter Weibull distributions is the smallest among these different distributions. For the sample data in Table 6, the allowable critical value corresponding to significance level of 0.05 is 0.270, and far greater than 0.0171 which is the maximum absolute difference between the observed probability of failure and the expected probability of failure associated with the sum model of two three-parameter Weibull distributions. The cdf of the sum model of two three-parameter Weibull distributions is expressed in eqn (15). Fig. 7 shows a comparison of the theoretical cdfs of failure and the observed failure probability for different distributions. It is obvious that the sum model of two threeparameter Weibull distributions is the best representation of the sample data among these distributions.
6 CONCLUSIONS J
Fig. 6. A flowchart for selection of population distribution and parameter estimation.
An optimization method for parameter estimation is proposed. The method is illustrated by estimating the population parameters for three-parameter Weibull distributions. A comparison between the proposed method and the usual
Parameter estimation
255
Table 7. Results of the parameter estimation in example 5 Distribution
Parameter
MaxlFe(Ti) - Fo(Ti)l
Exponential Normal Log normal Weibull Sum model
X = 0.0583 3" = 4.1006 /x = 17.9416 tr = 12.7525 ~' = 6.2043 a' = 0.3594 x0 = 5.1552xa = 16.1626/~ = 0.8665 3'1 = 3.42037h = 7.6764/~1 = 2.4547 p = 0.5767 3'2 = 6.0717'/2 = 31.5442/32 = 3.6680
0.0710 0.1308 0.0724 0.0619 0.0171
methods such as the least-squares method, the matching moments method and the m a x i m u m likelihood method shows that more reasonable estimates of the parameters are given by the proposed optimization method. In addition, the proposed method can be used to estimate the parameters for the sum of two three-parameter Weibull distributions. A new procedure for selecting population distribution and parameter estimation is also presented.
8.
9.
10.
REFERENCES 1. Kececioglu, D. B., Reliability and Life Testing Handbook, Vol. 1. PTR Prentice Hall, 1993. 2. Li, Yong-Ming A general linear-regression analysis applied to the three-parameter Weibull distribution. IEEE Trans. Reliability, 1994, R-43, 255-263. 3. Weibull, W., Fatigue Testing and Analysis of Results. Pergamon Press, Oxford, 1961. 4. Newby, M. Properties of moment estimators for the threeparameter Weibull distribution. 1EEE Trans. Reliability, 1984, R-33, 192-195. 5. Cran, G. W. Moment estimators for the three-parameter Weibull distribution. IEEE Trans. Reliability, 1988, R-37, 360363. 6. Rockette, H., Antle, C. E. & Klimko, L. A. Maximum likelihood estimation with the Weibull Distribution. J. Am. Statist. Assoc., 1974, 69, 246-249. 7. Harter, H. L. & Moore, A. H. Maximum likelihood estimation
11.
12.
13. 14.
15.
16.
17.
of the parameters of gamma and Weibull populations from complete and from censored samples. Technometrics, 1965, 7, 639-643. Sinha, S. K. & Sloan, J. A. Bayes estimation of the parameters and reliability function of the three-parameter Weibull distribution. IEEE Trans. Reliability, 1988, R-37, 364-369. Hobbs, J. R., Moore, A. H. & Miller, R. M. Minimum-distance estimation of the parameters of the three-parameters Weibull distribution. IEEE Trans. Reliability, 1985, R-34, 495 -496. Adatia, A. & Chan, L. K. Robust estimators of the threeparameter Weibull distribution. 1EEE Trans. Reliability, 1985, R-34, 347-351. Jiang, S. & Kececioglu, D, B. Graphical representation of two mixed Weibull distributions. IEEE Trans. Reliability, 1992, 41, 241-247. Kao, J. H. K. A graphical estimation of mixed Weibull parameters in life testing of electron tubes. Technometrics, 1959, 1, 389-407. Cran, H. Graphical estimation methods for Weibull distributions. Microelectronics and Reliability, 1976, 15, 47-52. Jensen, F. and Petersen, An Engineering Approach to the Design and Analysis of Burn-In Procedures. John Wiley and Sons, Chichester, 1982. Jiang, R. & Murthy, D. N. P. Reliability modeling involving two Weibull distributions. Reliability Engineering and System Safety, 1995, 47, 187-198. Madsen, K. & Schjaer-Jacobsen, H. Algorithms for worstcase tolerance optimization. IEEE Trans. Circuits and Systems, 1979, CAS-26, 775-783. Gao, Z. T., Applied Statistics to Fatigue. Press of Industry of National Defense, Beijing, 1986.