- Email: [email protected]
A new bathtub curve model with a finite support
Reliability Engineering and System Safety 119 (2013) 44–51
Contents lists available at SciVerse ScienceDirect
Reliability Engineering and System Safety journal homepage: www.elsevier.com/locate/ress
A new bathtub curve model with a finite support R. Jiang n Faculty of Automotive and Mechanical Engineering, Changsha University of Science and Technology, Changsha, Hunan 410114, China
art ic l e i nf o
a b s t r a c t
Article history: Received 4 February 2013 Received in revised form 25 April 2013 Accepted 19 May 2013 Available online 25 May 2013
The failure rate with a bathtub shape usually increases very fast in the wear-out phase. In this case, the bathtub curve model with a finite support can better adapt the sharp change in failure rate. There are few models with the finite support. This paper presents such a model. However, the maximum likelihood estimator of the location parameter of such models sometimes converges to the largest observation of a dataset. An extended maximum spacing method is developed to estimate the parameters for the case where the maximum likelihood method fails. Three examples are included to illustrate the appropriateness of the proposed model and estimation method. & 2013 Elsevier Ltd. All rights reserved.
Keywords: Failure rate Bathtub curve Finite support Extended maximum spacing method
1. Introduction In many applications in reliability, the failure rate function can be of bathtub shape. The plot of failure rate with a bathtub shape is called the bathtub curve, and the model with bathtub failure rate is called the bathtub curve model. The bathtub curve is characterized by three phases (i.e., early use phase, normal use phase and wear-out phase), which correspond to infant mortality, random failure and wear-out failure, respectively. As such, the bathtub curve can explain the failure behavior of the population of nonrepairable components. A large number of models have been developed (e.g., see [1–4], and the literature cited therein). Typical methods to establish a bathtub curve model include [5]:
Seeking an appropriate model, Modifying the customary model by introducing additional parameters, and
Transforming the data to achieve compatibility with a well understood and convenient customary model.
Most of the bathtub curve models are derived from the second method by exponentiating, generalizing or extending some wellknown models such as the Weibull and gamma distributions, and the resulting models usually have an infinite support. There are few models with a finite support.
n
Tel.: +86 731 82309122; fax: +86 731 85258646. E-mail address: [email protected]
0951-8320/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ress.2013.05.019
Crevecoeur [6] describes the bathtub curve as “decreases at first, then is rather flat, and finally increases (usually without limit)”. Zhang et al. [7] mention that a good bathtub curve model should have a relatively wide flat portion. These imply that the bathtub failure rate curve can increase very fast during the wearout phase so that it has two relatively apparent change points (for the discussions on the change points, see [8–13]). For the case where the failure rate increases quickly in the wear-out phase, the bathtub curve model with a finite support can better adapt the sharp change in failure rate. However, a major problem of the models with the finite support is that the location parameter sometimes converges to the largest observation of a given dataset when using the maximum likelihood method (MLM) to estimate the parameters. When this problem occurs, we say that the MLM fails or the maximum likelihood estimators (MLE) of the parameters do not exist. The maximum spacing method (MSM, or maximum product of spacings method, see [14]) and its variants (e.g., see [15] and the literature cited therein) have been developed to solve this problem. However, the MSM is not applicable for tied data. In this paper, a new three-parameter bathtub curve model with a finite support is proposed; and an extended MSM (EMSM) is developed to estimate the parameters for the case where the MLM fails. Three examples are included to illustrate the appropriateness of the proposed model and parameter estimation method. The paper is organized as follows. Typical bathtub curve models are presented in Section 2. The proposed model is presented in Section 3, and the EMSM is presented in Section 4. The appropriateness of the proposed model is illustrated in Section 5. A discussion on the model support is presented in Section 6. The paper is concluded in Section 7.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
2. Some typical bathtub curve models It is not the focus of this paper to examine all the existing bathtub curve models. For the purpose of comparison and illustration, we give some typical bathtub curve models with different numbers of parameters and supports in this section. To make the notation simple, we use β and δ to denote the shape parameters, η and α to denote the scale parameters, γ to denote the location parameter, and λ to denote the dimensionless coefficient. As such, the models presented here can be slightly different from the original in the form. 2.1. Bathtub curve models with infinite support 2.1.1. Two-parameter model In terms of cumulative hazard function, Chen [16] proposes the following model: β
HðtÞ ¼ λðet −1Þ:
ð1Þ
The failure rate (rðtÞ ¼ H′ðtÞ) is bathtub-shaped when β∈ð0; 1Þ.
45
2.2.1. Two-parameter model The inverse Pareto distribution is given by β t FðtÞ ¼ ; t∈ð0; γÞ: γ The failure rate is bathtub-shaped when β∈ð0; 1Þ.
2.2.2. Three-parameter models The beta distribution defined in the interval [0, γ] is given by f ðtÞ ¼
Γðβ þ δÞ t β−1 ðγ−tÞδ−1 ; γ βþδ−1 ΓðβÞΓðδÞ
β; δ; γ 4 0; t∈ð0; γÞ:
The failure rate is bathtub-shaped when β∈ð0; 1Þ. The Weibull distribution truncated at γ (TWD) is given by FðtÞ ¼
HðtÞ ¼ ðt=ηÞβ et=α :
The failure rate is bathtub-shaped when β o1 (see [26]).
ð2Þ
The failure rate is bathtub-shaped when β∈ð0; 1Þ. Mudholkar and Srivastava [18] propose the exponentiated Weibull distribution (EWD) with the cdf given by ( " #)δ t β ; δ; η; β 40: ð3Þ FðtÞ ¼ 1−exp − η The model has bathtub failure rate when β 41 and βδ o 1 (see [19]). Xie et al. [20] introduce a scale parameter to the model (1) and call it the modified Weibull extension distribution (MWED) β
λ; β; η 4 0
δ
1−exp½−ðt=ηÞβ ; 1−exp½−ðγ=ηÞβ
t∈ð0; γÞ:
ð11Þ
3. Proposed bathtub curve model 3.1. Model In terms of the failure rate, the proposed model (FSM for short) is given by rðtÞ ¼
β 1 þ ; t þ η γ−t
β; η; γ 4 0; t o γ:;
ð12Þ
The derivation process of the model is shown in Appendix A.
ð4Þ
El-Gohary et al. [21] present the generalized Gompertz distribution (GGD) given by FðtÞ ¼ f1−exp½−λðet=η −1Þg ;
ð9Þ
where Γð:Þ is the gamma function. The failure rate is bathtubshaped when β∈ð0; 1Þ. Mudholkar et al. [5] propose a data-transformation model (DTM) given by t=η β HðtÞ ¼ ; t∈ð0; γÞ: ð10Þ 1−t=γ
2.1.2. Three-parameter models Lai et al. [17] propose the modified Weibull distribution (MWD) given by
HðtÞ ¼ λðeðt=ηÞ −1Þ;
ð8Þ
λ; β; η 40:
ð5Þ
3.2. Characterization of failure rate function The failure rate at t ¼ 0 and γ are given respectively by rð0Þ ¼
When δ∈ð0; 1Þ, the model has the bathtub failure rate.
β 1 þ ; η γ
rðγÞ ¼ ∞:
ð13Þ
From (12), we have 2.1.3. Four-parameter models The two-fold Weibull competing risk model (e.g., see [22–24]) is given by β1 −1 β2 −1 β t β t rðtÞ ¼ 1 þ 2 : ð6Þ η1 η1 η2 η2 The failure rate is bathtub-shaped when β1 o 1 oβ2 . Sarhan and Apaloo [25] present the exponentiated modified Weibull extension given by FðtÞ ¼ f1−exp½−HðtÞgδ
ð7Þ
where HðtÞ is given by (4). The model includes (4) (when δ ¼ 1) and (5) (when β ¼ 1) as its special cases. 2.2. Bathtub curve models with finite support This category of models has a location parameter as the life upper limit.
r′ðtÞ ¼ −
β ðt þ ηÞ2
þ
1 ðγ−tÞ2
:
ð14Þ
From (14) we have r 0 ð0Þ ¼ −
β 1 þ ; η2 γ 2
r 0 ðγÞ ¼ ∞:
ð15Þ
From (15), we have
When γ ≤η=pffiffiffiβ, r′ð0Þ≥0 so that the failure rate is increasing, and When γ 4η=pffiffiffiβ, r′ð0Þ o 0 so that the failure rate is bathtubshaped. Let r m denote the minimum of rðtÞ, which is achieved at t ¼ t m . Letting r 0 ðtÞ ¼ 0, we have pffiffiffi pffiffiffi ð1 þ βÞ2 βγ−η pffiffiffi ; r m ¼ rðt m Þ ¼ : ð16Þ tm ¼ ηþγ 1þ β
46
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
3.3. Characterization of density function
maximizing the following log spacing function:
From (12), we have ð1 þ t=ηÞβ HðtÞ ¼ ln : 1−t=γ
lnðLS Þ ¼ ∑ ln½Fðt i Þ−Fðt i−1 Þ
nþ1
ð26Þ
i¼1
ð17Þ
The reliability and density functions are given respectively by 1−t=γ γ−t 1 RðtÞ ¼ ; f ðtÞ ¼ 1 þ β : ð18Þ β η þ t γð1 þ t=ηÞβ ð1 þ t=ηÞ
where Fðt 0 Þ ¼ 0 and Fðt nþ1 Þ ¼ 1. In practice, censored and tied observations occur and the spacing function given by (26) is not applicable for these cases. Cheng and Traylor [27] deal with the handling of censored sample and tied data. We develop a new approach to handle the tied data as follows.
From the second relation of (18), we have 4.2. Extended MSM
0
f ðtÞ ηþγ β ¼− o 0: − f ðtÞ ðη þ tÞ2 þ βðη þ tÞðγ−tÞ η þ t
ð19Þ
This implies that the pdf is decreasing with f ð0Þ ¼
1 β þ ; γ η
f ðγÞ ¼
1 1 o : γ γð1 þ γ=ηÞβ
ð20Þ
Titterington [28] interprets the MSM as a maximum grouped likelihood estimator with one observation in each group. The EMSM is based on this idea. Let di (≥1) denote the number of t i 's. Consider the tied data given by n
t 1 ðd1 Þ ot 2 ðd2 Þ o::: ot n ðdn Þ; N ¼ ∑ di :
ð27Þ
i¼1
3.4. Moments The k-th origin moment is given by Z γ Z γ mk ¼ t k f ðtÞdt ¼ k t k−1 RðtÞdt: 0
ð21Þ
0
Let ρ ¼ η=γ and z ¼ 1 þ t=η. Using (18) into (21) and after some simplification, we have Z 1þ1=ρ mk ¼ kηk ½ðz−1Þk−1 −ρðz−1Þk z−β dz: ð22Þ 1
ð23Þ
Letting k ¼ 2 yields the second-order origin moment given by ( m2 ¼ s2 þ μ2 ¼ 2η2
1þρ 1−β
lnðLi Þ ¼
di fln½Fðt i Þ−Fðt i−1 Þ þ ln Fðt iþ1 Þ−Fðt i Þ g: 2
ð28Þ
The overall grouped log-likelihood function is given by n
lnðLE Þ ¼ ∑ lnðLi Þ:
ð29Þ
i¼1
Letting k ¼ 1 yields the mean given by " # ρη 1 2−β 2−β 1þ −1 : − μ ¼ m1 ¼ ð1−βÞð2−βÞ ρ ρ
−
When d1 ¼ d2 ¼ ::: ¼ dn ¼ 1, (27) reduces into (25). Assume that t i 's are rough observations so that di =2 of them can fall in the interval (t i−1 ; t i ) and the other di =2 can fall in the interval (t i ; t iþ1 ). As such, the grouped log-likelihood function of the t i 's is given by
"
−ρ 3−β
1þ
1 ρ
" 1þ
# " # 1 3−β 1 þ 2ρ 1 2−β −1 þ −1 1þ ρ 2−β ρ
4.3. Example 1
#)
1−β −1
ð24Þ
where s2 is the life variance. In a similar way, the other moments can be obtained.
4. Extended maximum spacing method For a given dataset, if the MLEs of the parameters exist (this usually occurs when the largest observation is right-censored as illustrated later), we use the MLM to estimate the parameters; otherwise, we use the EMSM (to be presented in Section 4.2) to estimate the parameters. 4.1. Maximum spacing method The basis of the MSM is as follows. The spacings in the data are the differences between the values of the cdf at consecutive observations. The spacings follow the uniform distribution on [0, 1], and hence their geometric mean can achieve the maximum when the spacings are as equal as possible. Consider an ordered and complete dataset t 1 o t 2 o::: ot n :
The parameters are estimated by maximizing lnðLE Þ. It is noted that (29) is slightly different from (26) for the special case d1 ¼ d2 ¼ ⋯ ¼ dn ¼ 1, i.e., the coefficient of both ln½Fðt 1 Þ and ln½1−Fðt n Þ is 0.5 in (29) and is 1 in (26). The appropriateness of the EMSM is illustrated as follows:
ð25Þ
Assume that the data can be fitted by the distribution Fðt; θÞ, where θ is the parameter set. The parameters are estimated by
The data shown in Table 1 come from Jiang and Murthy [29], which are randomly generated from the Weibull distribution with the shape, scale and location parameters (β, η, γ)¼(2.5, 30, 10). For this example, the MLEs exist and are shown in the 3rd row of Table 2; the parameters estimated from the MSM and EMSM are shown in the 4th and 5th rows of Table 2, respectively. Since the true values of the parameters are known, the performance of the estimation methods can be evaluated by the relative errors in the shape parameter (denoted as εβ ) and mean life (denoted as εμ ), which are shown in the last two columns of Table 2. As seen, the EMSM is an improvement to the MSM. This will be further confirmed in Example 3.
5. Illustrations In this section, we present two examples (with small and moderate sample sizes) to illustrate the appropriateness of proposed model and estimation method. We also compare the performance of two categories of bathtub curve models. Each of the models for comparison has three parameters. The models with Table 1 Data for Example 1. 17.8 32.3 43.4
21.3 33.5 44.5
23.8 34.9 47.0
25.9 36.6 48.8
27.4 38.5 52.5
29.4 39.7 61.4
30.6 41.2
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
Table 2 Estimated parameters for Example 1.
True MLM MSM EMSM
β
η
γ
μ
2.5
30 25.08
10 14.30
30.29 27.84
9.96 11.88
36.62 36.51 36.79 36.55
2.1652 2.2923 2.3766
εβ (%)
lnðLÞ
−75.367 −75.919 −75.446
εμ (%)
0.29 0.47 0.17
Table 3 Data for Example 2. 2 106 275
10 143 293
13 147 300+
23 173 300+
23 181 300+
28 212 300+
30 245 300+
65 247 300+
80 261 300+
88 266 300+
Table 4 Estimated parameters and goodness-of-fit measures for Example 2. Model
β
η
δ,α
MWD (2) EWD (3) MWED (4) GGD (5) BETA (9) DTM (10) TWD (11)
0.6537
633.56
408.28
5.4197 0.6888
447.96 186.14 186.15
0.1370
FSM (12)
γ,λ
0.6165 0.6190
653.25
0.4219 0.2223 354.30 423.90
0.7334 0.06674
111678 9.5118
459.39 452.35
0.6086 0.5433
SSE
p
−141.869
0.0289
−141.736 −142.033 −141.774 −141.408 −141.422 −141.680 −141.356
0.0386 0.0356 0.0262 0.0195 0.0188 0.0345 0.0288
0.2061 0.2116 0.1647 0.2597 0.4081 0.3787 0.2226 0.2211
lnðLÞ
the infinite support are Models (2)–(5); and the models with the finite support are Models (9)–(11). 5.1. Goodness-of-fit measures For a given dataset, we evaluate the goodness-of-fit of a fitted model using three measures. The first measure is the loglikelihood value lnðLÞ (not lnðLE Þ for the EMSM). The second measure is the sum of squared errors given by n
SSE ¼ ∑ di ½Fðt i Þ−F i 2
It is noted that the scale parameter of the truncated Weibull model given by (11) is very large. In this case, the model can be approximated by FðtÞ≈
13.39 8.31 4.94
ð31Þ
5.2.3. Comparison of models The models can be ranked based on each of lnðLÞ, SSE and p, and the results are shown in the 2–4th columns of Table 5. As seen, the proposed model is the best in terms of lnðLÞ; the model given by (10) is the best in terms of SSE; and the beta model given by (9) is the best in terms of p. This confirms the appropriateness of the proposed model. The median empirical CDF and the CDFs of the fitted models with the infinite and finite supports are shown in Figs. 1 and 2, respectively. As seen, the median empirical CDF is inverse S-shaped; and the CDFs of the models with the finite support are more similar to the median empirical CDF than those of the models with the infinite support in shape. Table 5 Rank numbers of the fitted models. Model
Rank no. for Example 2 Rank no. for Example 3 Average of ranks
MWD (2) EWD (3) MWED (4) GGD (5) BETA (9) DTM (10) TWD (11) FSM (12)
lnðLÞ
SSE
p
lnðLÞ
SSE
p
7 5 8 6 2 3 4 1
5 8 7 3 2 1 6 4
7 6 8 3 1 2 4 5
7 5 8 6 4 3 2 1
7 6 8 5 4 2 3 1
7 6 8 5 4 1 3 2
150
200
6.7 6.0 7.8 4.7 2.8 2.0 3.7 2.3
1
ð30Þ
F(t)
0.8 0.6 0.4 0.2 0
0
50
100
250
300
350
400
450
t Fig. 1. Median empirical CDF and the CDFs of the fitted models with the infinite support for Example 2.
5.2. Example 2
1
(9)
(11)
0.8
F(t)
5.2.1. Data The failure and running times of 30 devices shown in Table 3 come from Meeker and Escobar [33]. The sign “þ” indicates the running times. It is noted that the largest observation is rightcensored and there are two failures at t ¼ 23.
β ðt=ηÞβ t ¼ : γ ðγ=ηÞβ
This is the model given by (8), implying that the truncated Weibull model reduces into the inverse Pareto distribution in this case.
i¼1
where Fðt i Þ is the cdf of the fitted model and F i is the median empirical cdf at t i . Microsoft Excel provides a standard function to evaluate the median of the beta distribution (i.e., Betainv ð0:5; m; N−m þ 1Þ, where m ¼ Σ ij ¼ 1 dj ), and hence we do not need to use its approximation (e.g., see [30,31]). The third measure is the χ 2 statistic and corresponding p value. The critical value of p is usually in the range 0.1–0.3 (see, for example, [32]). Since the degree of freedom for computing the p value must at least be 1, the number of data intervals has to be larger than or equal to 5 for the three-parameter models under consideration.
47
0.6
(12)
0.4
(10)
0.2
5.2.2. Parameter estimation The MLEs exist for this example. To compute the p value, we divide the dataset into five intervals with an interval length of 75, and the last interval is half-open. The estimated parameter and the values of lnðLÞ, SSE and p are shown in Table 4.
0
0
50
100
150
200
250
300
350
400
450
t Fig. 2. Median empirical CDF and the CDFs of the fitted models with the finite support for Example 2.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
The bathtub curves of the models with the infinite and finite supports are shown in Figs. 3 and 4, respectively. As seen, the failure rates of the models with the finite support have faster increase rate in the wear-out phase than those of the models with the infinite support. 5.2.4. Some results derived from model (12) The proposed model is mathematically tractable, and hence it is easy to derive some important characteristics (e.g., moments, initial and minimal failure rates) of the devices. The results are shown in the second row of Table 6. The density function obtained from the fitted model is displayed in Fig. 5. As seen, the pdf decreases and the rate of decrease is also decreasing with time.
Table 6 Characteristics derived from Model (12). Example
μ
s
s=μ
rð0Þ,10−3
tm
rðt m Þ,10−3
2 3
191.40 173.41
138.74 135.87
0.7248 0.7836
9.23 9.04
85.31 99.65
3.43 3.98
0.01 0.008 0.006 f(t)
48
Example 3
0.004 0.002
5.3. Example 3
0
5.3.1. Data The data shown in Table 7 deal with failure times of electronic devices reported in Wang [34]. It is noted that the largest observation value is a failure observation. 5.3.2. Parameter estimation For the four models with the finite support, the MLEs exist for the model (10), and the MLM fails for the other three models. We use the EMSM to estimate their parameters. To further confirm the appropriateness of the EMSM, we fix the value of γ obtained from the EMSM and re-estimate the other two parameters using the MLM. To compute the p value, we divide the data into five intervals so that the degree of freedom equals 1; the interval length is taken as 81 so that the number of data in each interval is not smaller than 3. The estimated parameters and corresponding values of goodness-of-fit measures are shown in Table 8. To further confirm the appropriateness of the EMSM, we look at the goodness-of-fit measures of the fitted models (9), (11) and (12). Though the values of lnðLÞ obtained from the EMSM are slightly smaller than those obtained from the MLM (with γ being the same), the values of SSE are smaller and the values of p are larger than the corresponding values obtained from the MLM. This further validates the EMSM. 0.015 (3) (5) 0.01
0
100
200
300
400
500
t Fig. 5. Density functions derived from Model (12) for Examples 2 and 3.
Table 7 Data for Example 3. 5 98 245
11 122 293
21 145 321
31 165 330
46 195 350
75 224 420
5.3.3. Comparison of models The rank numbers of the fitted models are shown in the 5th– 7th columns of Table 5. As seen, the model proposed in this paper is the best in terms of both lnðLÞ and SSE; and the model given by (10) is the best in terms of p. This further confirms the appropriateness of the proposed model. The median empirical CDF and the CDFs of the fitted models with the infinite and finite supports are shown in Figs. 6 and 7, respectively. As seen, the CDFs of the models with the finite support are closer to the median empirical CDF than those of the models with the infinite support. The bathtub curves of the models with the infinite and finite supports are shown in Figs. 8 and 9, respectively. As seen, the failure rates of the models with the finite support have longer flat portion and faster increase rate in the wear-out phase than those of the models with the infinite support.
r(t)
(2)
5.3.4. Some results derived from model (12) Some important characteristics of the electronic devices are shown in the last row of Table 6. The density function is also displayed in Fig. 5.
(4)
0.005
0
Example 2
0
50
100
150
200
250
300
350
400
450
t
5.4. Summary
Fig. 3. Bathtub curves of the models with the infinite support for Example 2.
To integrate the results from the two examples, we calculate the average of the six rank numbers for each model and the results are shown in the last column of Table 5. In terms of the average rank number, the proposed model is the second best and the models with the finite support have smaller averages than the models with the infinite support.
0.015 (12)
(11)
(10)
0.01
r(t)
(9)
0.005
6. Discussion on the model support 0
0
50
100
150
200
250
300
350
400
450
t Fig. 4. Bathtub curves of the models with the finite support for Example 2.
The possible concerns with the model with the finite support are physical reasons and possible consequences of using such a model. We discuss these issues as follows.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
49
Table 8 Estimated parameters and goodness-of-fit measures for Example 3. β
MWD (2) EWD (3) MWED (4) GGD (5) BETA (9), EMSM BETA (9), MLM DTM (10) TWD (11), EMSM TWD (11), MLM FSM (12), EMSM FSM (12), MLM
0.6468 6.2271 0.7522
η
δ,α 665.71 373.00 134.04 144.95
0.7383 0.7793 0.6871 0.7549 0.8049 0.12073 0.1304
0.6012 1.2081 1.2781
410.34 1202.20 729.34 17.88 20.97
SSE
p
0.3413 0.1831 449.54 449.54 480.55 440.38 440.38 437.89 437.89
−108.933 −108.397 −109.117 -108.695 −108.162 −108.143 −108.145 -108.069 −108.055 −107.833 −107.830
0.0235 0.0231 0.0341 0.0184 0.0100 0.0135 0.0080 0.0084 0.0108 0.0062 0.0066
0.4175 0.4542 0.3471 0.4930 0.6375 0.5950 0.7446 0.6435 0.6160 0.6801 0.6750
5
0.8
4
0.6
3
0.4
0
0
50
100
150
200
250
300
350
400
450
Fig. 6. Median empirical CDF and the CDFs of the fitted models with the infinite support for Example 3.
1 0.8 0.6 0.4 0.2 0
Truncated Weibull
Weibull
2 1
t
F(t)
lnðLÞ
1
0.2
0
50
100
150
200
250
300
350
400
450
t Fig. 7. Median empirical CDF and the CDFs of the fitted models with the finite support for Example 3.
0.05 0.04
(3)
0.03
r(t)
γ,λ
276.88 0.1234
r(t)
F(t)
Model
(5) (2)
0.02
(4)
0.01 0
0
50
100
150
200
250
300
350
400
450
t Fig. 8. Bathtub curves of the models with the infinite support for Example 3. 0.05
(9)
0.04 (10)
r(t)
0.03 0.02 0.01 (12) 0
(11) 0
50
100
150
200
250
300
350
400
450
t
Fig. 9. Bathtub curves of the models with the finite support for Example 3.
0
t2 0
5
10
15 t
20
25
30
Fig. 10. Failure rate functions of the Weibull and truncated-Weibull distributions.
From a perspective of engineering, the model with the finite support is more appropriate than the model with the infinite support since the lifetime of a component is always finite; from a mathematical perspective, the model with the infinite support is more attractive since it avoids the problem to specify the upper bound of the lifetime. As such, the lifetime distribution with the infinite support can be viewed as an approximate representation of the lifetime. The selection between the model with the finite support and the model with the infinite support is problem-specific. In the context of bathtub failure rate modeling, we have shown above that the model with the finite support is preferred due to the favorite property that the failure rate can increase more quickly in the wear-out phase. When using the finite support model, the location parameter γ should be viewed as a point estimate of the upper bound of lifetime. When γ is appropriately estimated, the probability of event T≥γ should be small, but this does not imply that T≥γ is impossible. The value of γ provides useful information about the upper bound of lifetime, and can be applied in reliability and maintenance analysis. In the context of preventive maintenance, the component is at best preventively replaced at or before the second change point of the bathtub curve. The method to specify this change point is out of the scope of this paper. There is not a big gap between the model with the finite support and the model with infinite support. To illustrate, we consider the Weibull distribution Fðt; β; ηÞ with (β, η)¼ (3.5, 10.0). It can be well approximated by the truncated Weibull distribution given by (11) with γ ¼ 21:17 and the same values of β and η. This is because FðγÞ ¼ 1−10−6 so that the plots of their pdfs or cdfs are almost overlapped in (0, γ). However, the failure rate of the truncated-Weibull distribution quickly increases after some time point (e.g., t 2 ¼ 19:24, see Fig. 10). Generally, a distribution with the infinite support can be well approximated by its truncated distribution with the finite support as long as γ is appropriately specified; and obtaining an estimate of
50
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
the lifetime upper bound is more practical than assuming that the lifetime upper bound is infinite.
ðε1 þ ε2 þ ε1 ε2 Þr 0 ðtÞ≈−r 2 ðtÞ:
ðA:4Þ
Let β ¼ ε1 þ ε2 þ ε1 ε2 , which is positive and small. From (A.4) and after some simplifications, we have
7. Conclusions In this paper we developed a new three-parameter bathtub curve model with a finite support. Two examples have been used to illustrate its appropriateness. The new model serves as a good alternative for modeling bathtub-shaped failure rate data. An extended maximum spacing method has been developed to estimate the parameters of the models with a location parameter. The appropriateness of the method has been illustrated by two examples. We have compared the performances of the bathtub curve models with infinite and finite supports, and the results support the use of the models with the finite support. The other findings have been 1. The model given by (10) and the beta distribution are good alternatives for fitting the data with a bathtub failure rate. 2. The upper truncated Weibull model reduces into the inverse Pareto distribution when the scale parameter is very large. The results obtained in this paper are useful for the researchers who are attempting to develop new bathtub curve models; and for the practitioners who need to choose an appropriate bathtub curve model for fitting a given dataset that is considered to have a bathtub failure rate. A topic for the future study is to develop an appropriate method to specify the change points of the bathtub curve.
Acknowledgment The author would like to thank the two referees for their helpful comments and suggestions which have greatly enhanced the clarity of the paper. The research was supported by the National Natural Science Foundation (No. 71071026).
Appendix A. Derivation of the model (12) We use an approximate approach to derive the model (12) as follows. The starting point is the failure rate function given by rðtÞ ¼ f ðtÞ=RðtÞ:
ðA:1Þ
From (A.1), we have r′ðtÞ ¼ r 2 ðtÞ þ rðtÞ
Using (A.3) to (A.2) we have
f ′ðtÞ : f ðtÞ
ðA:2Þ
For a bathtub failure rate, rðtÞ has a minimum value r m at t ¼ t m . Using t m as a partition point, the bathtub curve can be divided into a decreasing part and an increasing part. We first consider the decreasing part. When t is small, we have:
rðtÞ 4 f ðtÞ but rðtÞ is close to f ðtÞ, and jr′ðtÞj o jf ′ðtÞj but r′ðtÞ is close to f ′ðtÞ. Let ε1 and ε2 be small positive numbers. The above two groups of relations can be written as follows: rðtÞ ≈1 þ ε1 ; f ′ðtÞ≈ð1 þ ε2 Þr′ðtÞ: f ðtÞ
ðA:3Þ
rðtÞ≈
β ; η 4 0: tþη
ðA:5Þ
We now consider the increasing part. For large t (t m ≪t o γ), we 0 assume that f ðtÞ decreases slowly so that f ðtÞ=f ðtÞ o0 and ≈0. In this case, (A.2) can be written as below r′ðtÞ≈r 2 ðtÞ:
ðA:6Þ
From (A.6) and after some simplifications, we have rðtÞ≈
1 γ−t
ðA:7Þ
where γ is a positive constant and can be viewed as a life upper limit. It is noted that (A.7) is actually the failure rate of the uniform distribution defined over (0; γ), which is an increasing function of t. Eq. (A.5) represents the early failure mode and Eq. (A.7) represents the wear-out failure mode. As such, the total failure rate is given by rðtÞ ¼
β 1 þ ; β; η; γ 40; t o γ: t þ η γ−t
ðA:8Þ
This is the model (12). References [1] Rajarshi S, Rajarshi MB. Bathtub distributions—a review. Communications in Statistics—Theory and Methods 1988;17(8):2597–621. [2] Lai CD, Xie M. Stochastic ageing and dependence for reliability. New York: Springer; 2006 [Chapter 3]. [3] Pham H, Lai CD. On recent generalizations of the Weibull distribution. IEEE Transactions on Reliability 2007;56:454–8. [4] Bebbington M, Lai CD, Zitikis R. Bathtub curves in reliability and beyond. Australian and New Zealand Journal of Statistics 2007;49:251–65. [5] Mudholkar GS, Asubonteng KO, Hutson AD. Transformation of the bathtub failure rate data in reliability for using Weibull-model analysis. Statistical Methodology 2009;6:622–33. [6] Crevecoeur GU. A model for the integrity assessment of ageing repairable systems. IEEE Transactions on Reliability 1993;42(1):148–55. [7] Zhang T, Xie M, Tang LC, Ng SH. Reliability and modeling of systems integrated with firmware and hardware. International Journal of Reliability, Quality and Safety Engineering 2005;12(3):227–39. [8] Clarotti C, Lannoy A, Odin S, Procaccia H. Detection of equipment aging and determination of the efficiency of a corrective measure. Reliability Engineering and System Safety 2004;84(1):57–64. [9] Bebbington M, Lai CD, Zitikis R. Useful periods for lifetime distributions with bathtub shaped failure rate functions. IEEE Transactions on Reliability 2006;55 (2):245–51. [10] Bebbington M, Lai CD, Zitikis R. Estimating the turning point of a bathtub shaped failure distribution. Journal of Statistical Planning and Inference 2008;138(4):1157–66. [11] Kim J. A change-point estimator with the hazard ratio. Journal of the Korean Statistical Society 2009;38(4):377–82. [12] Jiang R. Determination of two change points of a bathtub failure rate curve. In: Proceedings of the third international Asia conference on industrial engineering and management innovation. Beijing; 2012, p. 1116–20. [13] Jiang R. Specification of change points of failure rate or intensity function: a non-parametric approach. In: Proceedings of the IEEE international conference on industrial engineering and engineering management. Hong Kong; 2012, p. 747–51. [14] Cheng RC, Amin NAK. Estimating parameters in continuous univariate distributions with a shifted origin. Journal of Royal Statistical Society B 1983;45:394–403. [15] Ekström M. Alternatives to maximum likelihood estimation based on spacings and the Kullback–Leibler divergence. Journal of Statistical Planning and Inference 2008;138(6):1778–91. [16] Chen Z. A new two-parameter lifetime distribution with bathtub shape or increasing failure rate function. Statistics and Probability Letters 2000;49: 155–61. [17] Lai CD, Xie M. Murthy DNP. A modified Weibull distribution. IEEE Transactions on Reliability 2003;52(1):33–7. [18] Mudholkar GS, Srivastava DK. Exponentiated Weibull family for analyzing bathtub failure-rate data. IEEE Transactions on Reliability 1993;42(2):299–302.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
[19] Jiang R, Murthy DNP. The exponentiated Weibull family: a graphical approach. IEEE Transactions on Reliability 1999;48(1):68–72. [20] Xie M, Tang Y, Goh TNA. Modified Weibull extension with bathtub-shaped failure rate function. Reliability Engineering and System Safety 2002;76: 279–85. [21] El-Gohary A, Alshamrani A, Al-Otaibi A. The generalized gompertz distribution. Applied Mathematical Modeling 2013;37(1-2):13–24. [22] Jiang R, Murthy DNP. Reliability modeling involving two Weibull distributions. Reliability Engineering and System Safety 1995;47:187–98. [23] Jiang R, Murthy DNP. Parametric study of competing risk model involving two Weibull distributions. International Journal of Reliability, Quality and Safety Engineering 1997;4(1):17–34. [24] Xie M, Lai CD. Reliability analysis using an additive Weibull model with bathtub-shaped failure rate function. Reliability Engineering and System Safety 1996;52(1):87–93. [25] Sarhan AM, Apaloo J. Exponentiated modified Weibull extension distribution. Reliability Engineering and System Safety 2013;112(4):137–44. [26] Zhang T, Xie M. On the upper truncated Weibull distribution and its reliability implications. Reliability Engineering and System Safety 2011;96:194–200.
51
[27] Cheng RCH, Traylor L. Non-regular maximum likelihood problems. Journal of the Royal Statistical Society: Series B 1995;57(1):3–44. [28] Titterington DM. Comment on estimating parameters in continuous univariate distribution. Journal of the Royal Statistical Society: Series B 1985;47:115–6. [29] Jiang R, Murthy DNP. Comments on: a general linear-regression analysis applied to the 3 parameter Weibull model. IEEE Transactions on Reliability 1997;46(3):389–93. [30] Hyndman RJ, Fan Y. Sample quantiles in statistical packages. The Amerian Statistician 1996;50:361–5. [31] Folland C, Anderson C. Estimating changing extremes using empirical ranking methods. Journal of Climate 2002;15:2954–60. [32] Blischke WR, Murthy DNP. Reliability: modeling, prediction, and optimization. New York: John Wiley & Sons; 2000. [33] Meeker WQ, Escobar LA. Statistical methods for reliability data. New York: John Wiley; 1998. [34] Wang FK. A new model with bathtub-shaped failure rate using an additive Burr XII distribution. Reliability Engineering and System Safety 2000;70(3):305–12.
Contents lists available at SciVerse ScienceDirect
Reliability Engineering and System Safety journal homepage: www.elsevier.com/locate/ress
A new bathtub curve model with a finite support R. Jiang n Faculty of Automotive and Mechanical Engineering, Changsha University of Science and Technology, Changsha, Hunan 410114, China
art ic l e i nf o
a b s t r a c t
Article history: Received 4 February 2013 Received in revised form 25 April 2013 Accepted 19 May 2013 Available online 25 May 2013
The failure rate with a bathtub shape usually increases very fast in the wear-out phase. In this case, the bathtub curve model with a finite support can better adapt the sharp change in failure rate. There are few models with the finite support. This paper presents such a model. However, the maximum likelihood estimator of the location parameter of such models sometimes converges to the largest observation of a dataset. An extended maximum spacing method is developed to estimate the parameters for the case where the maximum likelihood method fails. Three examples are included to illustrate the appropriateness of the proposed model and estimation method. & 2013 Elsevier Ltd. All rights reserved.
Keywords: Failure rate Bathtub curve Finite support Extended maximum spacing method
1. Introduction In many applications in reliability, the failure rate function can be of bathtub shape. The plot of failure rate with a bathtub shape is called the bathtub curve, and the model with bathtub failure rate is called the bathtub curve model. The bathtub curve is characterized by three phases (i.e., early use phase, normal use phase and wear-out phase), which correspond to infant mortality, random failure and wear-out failure, respectively. As such, the bathtub curve can explain the failure behavior of the population of nonrepairable components. A large number of models have been developed (e.g., see [1–4], and the literature cited therein). Typical methods to establish a bathtub curve model include [5]:
Seeking an appropriate model, Modifying the customary model by introducing additional parameters, and
Transforming the data to achieve compatibility with a well understood and convenient customary model.
Most of the bathtub curve models are derived from the second method by exponentiating, generalizing or extending some wellknown models such as the Weibull and gamma distributions, and the resulting models usually have an infinite support. There are few models with a finite support.
n
Tel.: +86 731 82309122; fax: +86 731 85258646. E-mail address: [email protected]
0951-8320/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ress.2013.05.019
Crevecoeur [6] describes the bathtub curve as “decreases at first, then is rather flat, and finally increases (usually without limit)”. Zhang et al. [7] mention that a good bathtub curve model should have a relatively wide flat portion. These imply that the bathtub failure rate curve can increase very fast during the wearout phase so that it has two relatively apparent change points (for the discussions on the change points, see [8–13]). For the case where the failure rate increases quickly in the wear-out phase, the bathtub curve model with a finite support can better adapt the sharp change in failure rate. However, a major problem of the models with the finite support is that the location parameter sometimes converges to the largest observation of a given dataset when using the maximum likelihood method (MLM) to estimate the parameters. When this problem occurs, we say that the MLM fails or the maximum likelihood estimators (MLE) of the parameters do not exist. The maximum spacing method (MSM, or maximum product of spacings method, see [14]) and its variants (e.g., see [15] and the literature cited therein) have been developed to solve this problem. However, the MSM is not applicable for tied data. In this paper, a new three-parameter bathtub curve model with a finite support is proposed; and an extended MSM (EMSM) is developed to estimate the parameters for the case where the MLM fails. Three examples are included to illustrate the appropriateness of the proposed model and parameter estimation method. The paper is organized as follows. Typical bathtub curve models are presented in Section 2. The proposed model is presented in Section 3, and the EMSM is presented in Section 4. The appropriateness of the proposed model is illustrated in Section 5. A discussion on the model support is presented in Section 6. The paper is concluded in Section 7.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
2. Some typical bathtub curve models It is not the focus of this paper to examine all the existing bathtub curve models. For the purpose of comparison and illustration, we give some typical bathtub curve models with different numbers of parameters and supports in this section. To make the notation simple, we use β and δ to denote the shape parameters, η and α to denote the scale parameters, γ to denote the location parameter, and λ to denote the dimensionless coefficient. As such, the models presented here can be slightly different from the original in the form. 2.1. Bathtub curve models with infinite support 2.1.1. Two-parameter model In terms of cumulative hazard function, Chen [16] proposes the following model: β
HðtÞ ¼ λðet −1Þ:
ð1Þ
The failure rate (rðtÞ ¼ H′ðtÞ) is bathtub-shaped when β∈ð0; 1Þ.
45
2.2.1. Two-parameter model The inverse Pareto distribution is given by β t FðtÞ ¼ ; t∈ð0; γÞ: γ The failure rate is bathtub-shaped when β∈ð0; 1Þ.
2.2.2. Three-parameter models The beta distribution defined in the interval [0, γ] is given by f ðtÞ ¼
Γðβ þ δÞ t β−1 ðγ−tÞδ−1 ; γ βþδ−1 ΓðβÞΓðδÞ
β; δ; γ 4 0; t∈ð0; γÞ:
The failure rate is bathtub-shaped when β∈ð0; 1Þ. The Weibull distribution truncated at γ (TWD) is given by FðtÞ ¼
HðtÞ ¼ ðt=ηÞβ et=α :
The failure rate is bathtub-shaped when β o1 (see [26]).
ð2Þ
The failure rate is bathtub-shaped when β∈ð0; 1Þ. Mudholkar and Srivastava [18] propose the exponentiated Weibull distribution (EWD) with the cdf given by ( " #)δ t β ; δ; η; β 40: ð3Þ FðtÞ ¼ 1−exp − η The model has bathtub failure rate when β 41 and βδ o 1 (see [19]). Xie et al. [20] introduce a scale parameter to the model (1) and call it the modified Weibull extension distribution (MWED) β
λ; β; η 4 0
δ
1−exp½−ðt=ηÞβ ; 1−exp½−ðγ=ηÞβ
t∈ð0; γÞ:
ð11Þ
3. Proposed bathtub curve model 3.1. Model In terms of the failure rate, the proposed model (FSM for short) is given by rðtÞ ¼
β 1 þ ; t þ η γ−t
β; η; γ 4 0; t o γ:;
ð12Þ
The derivation process of the model is shown in Appendix A.
ð4Þ
El-Gohary et al. [21] present the generalized Gompertz distribution (GGD) given by FðtÞ ¼ f1−exp½−λðet=η −1Þg ;
ð9Þ
where Γð:Þ is the gamma function. The failure rate is bathtubshaped when β∈ð0; 1Þ. Mudholkar et al. [5] propose a data-transformation model (DTM) given by t=η β HðtÞ ¼ ; t∈ð0; γÞ: ð10Þ 1−t=γ
2.1.2. Three-parameter models Lai et al. [17] propose the modified Weibull distribution (MWD) given by
HðtÞ ¼ λðeðt=ηÞ −1Þ;
ð8Þ
λ; β; η 40:
ð5Þ
3.2. Characterization of failure rate function The failure rate at t ¼ 0 and γ are given respectively by rð0Þ ¼
When δ∈ð0; 1Þ, the model has the bathtub failure rate.
β 1 þ ; η γ
rðγÞ ¼ ∞:
ð13Þ
From (12), we have 2.1.3. Four-parameter models The two-fold Weibull competing risk model (e.g., see [22–24]) is given by β1 −1 β2 −1 β t β t rðtÞ ¼ 1 þ 2 : ð6Þ η1 η1 η2 η2 The failure rate is bathtub-shaped when β1 o 1 oβ2 . Sarhan and Apaloo [25] present the exponentiated modified Weibull extension given by FðtÞ ¼ f1−exp½−HðtÞgδ
ð7Þ
where HðtÞ is given by (4). The model includes (4) (when δ ¼ 1) and (5) (when β ¼ 1) as its special cases. 2.2. Bathtub curve models with finite support This category of models has a location parameter as the life upper limit.
r′ðtÞ ¼ −
β ðt þ ηÞ2
þ
1 ðγ−tÞ2
:
ð14Þ
From (14) we have r 0 ð0Þ ¼ −
β 1 þ ; η2 γ 2
r 0 ðγÞ ¼ ∞:
ð15Þ
From (15), we have
When γ ≤η=pffiffiffiβ, r′ð0Þ≥0 so that the failure rate is increasing, and When γ 4η=pffiffiffiβ, r′ð0Þ o 0 so that the failure rate is bathtubshaped. Let r m denote the minimum of rðtÞ, which is achieved at t ¼ t m . Letting r 0 ðtÞ ¼ 0, we have pffiffiffi pffiffiffi ð1 þ βÞ2 βγ−η pffiffiffi ; r m ¼ rðt m Þ ¼ : ð16Þ tm ¼ ηþγ 1þ β
46
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
3.3. Characterization of density function
maximizing the following log spacing function:
From (12), we have ð1 þ t=ηÞβ HðtÞ ¼ ln : 1−t=γ
lnðLS Þ ¼ ∑ ln½Fðt i Þ−Fðt i−1 Þ
nþ1
ð26Þ
i¼1
ð17Þ
The reliability and density functions are given respectively by 1−t=γ γ−t 1 RðtÞ ¼ ; f ðtÞ ¼ 1 þ β : ð18Þ β η þ t γð1 þ t=ηÞβ ð1 þ t=ηÞ
where Fðt 0 Þ ¼ 0 and Fðt nþ1 Þ ¼ 1. In practice, censored and tied observations occur and the spacing function given by (26) is not applicable for these cases. Cheng and Traylor [27] deal with the handling of censored sample and tied data. We develop a new approach to handle the tied data as follows.
From the second relation of (18), we have 4.2. Extended MSM
0
f ðtÞ ηþγ β ¼− o 0: − f ðtÞ ðη þ tÞ2 þ βðη þ tÞðγ−tÞ η þ t
ð19Þ
This implies that the pdf is decreasing with f ð0Þ ¼
1 β þ ; γ η
f ðγÞ ¼
1 1 o : γ γð1 þ γ=ηÞβ
ð20Þ
Titterington [28] interprets the MSM as a maximum grouped likelihood estimator with one observation in each group. The EMSM is based on this idea. Let di (≥1) denote the number of t i 's. Consider the tied data given by n
t 1 ðd1 Þ ot 2 ðd2 Þ o::: ot n ðdn Þ; N ¼ ∑ di :
ð27Þ
i¼1
3.4. Moments The k-th origin moment is given by Z γ Z γ mk ¼ t k f ðtÞdt ¼ k t k−1 RðtÞdt: 0
ð21Þ
0
Let ρ ¼ η=γ and z ¼ 1 þ t=η. Using (18) into (21) and after some simplification, we have Z 1þ1=ρ mk ¼ kηk ½ðz−1Þk−1 −ρðz−1Þk z−β dz: ð22Þ 1
ð23Þ
Letting k ¼ 2 yields the second-order origin moment given by ( m2 ¼ s2 þ μ2 ¼ 2η2
1þρ 1−β
lnðLi Þ ¼
di fln½Fðt i Þ−Fðt i−1 Þ þ ln Fðt iþ1 Þ−Fðt i Þ g: 2
ð28Þ
The overall grouped log-likelihood function is given by n
lnðLE Þ ¼ ∑ lnðLi Þ:
ð29Þ
i¼1
Letting k ¼ 1 yields the mean given by " # ρη 1 2−β 2−β 1þ −1 : − μ ¼ m1 ¼ ð1−βÞð2−βÞ ρ ρ
−
When d1 ¼ d2 ¼ ::: ¼ dn ¼ 1, (27) reduces into (25). Assume that t i 's are rough observations so that di =2 of them can fall in the interval (t i−1 ; t i ) and the other di =2 can fall in the interval (t i ; t iþ1 ). As such, the grouped log-likelihood function of the t i 's is given by
"
−ρ 3−β
1þ
1 ρ
" 1þ
# " # 1 3−β 1 þ 2ρ 1 2−β −1 þ −1 1þ ρ 2−β ρ
4.3. Example 1
#)
1−β −1
ð24Þ
where s2 is the life variance. In a similar way, the other moments can be obtained.
4. Extended maximum spacing method For a given dataset, if the MLEs of the parameters exist (this usually occurs when the largest observation is right-censored as illustrated later), we use the MLM to estimate the parameters; otherwise, we use the EMSM (to be presented in Section 4.2) to estimate the parameters. 4.1. Maximum spacing method The basis of the MSM is as follows. The spacings in the data are the differences between the values of the cdf at consecutive observations. The spacings follow the uniform distribution on [0, 1], and hence their geometric mean can achieve the maximum when the spacings are as equal as possible. Consider an ordered and complete dataset t 1 o t 2 o::: ot n :
The parameters are estimated by maximizing lnðLE Þ. It is noted that (29) is slightly different from (26) for the special case d1 ¼ d2 ¼ ⋯ ¼ dn ¼ 1, i.e., the coefficient of both ln½Fðt 1 Þ and ln½1−Fðt n Þ is 0.5 in (29) and is 1 in (26). The appropriateness of the EMSM is illustrated as follows:
ð25Þ
Assume that the data can be fitted by the distribution Fðt; θÞ, where θ is the parameter set. The parameters are estimated by
The data shown in Table 1 come from Jiang and Murthy [29], which are randomly generated from the Weibull distribution with the shape, scale and location parameters (β, η, γ)¼(2.5, 30, 10). For this example, the MLEs exist and are shown in the 3rd row of Table 2; the parameters estimated from the MSM and EMSM are shown in the 4th and 5th rows of Table 2, respectively. Since the true values of the parameters are known, the performance of the estimation methods can be evaluated by the relative errors in the shape parameter (denoted as εβ ) and mean life (denoted as εμ ), which are shown in the last two columns of Table 2. As seen, the EMSM is an improvement to the MSM. This will be further confirmed in Example 3.
5. Illustrations In this section, we present two examples (with small and moderate sample sizes) to illustrate the appropriateness of proposed model and estimation method. We also compare the performance of two categories of bathtub curve models. Each of the models for comparison has three parameters. The models with Table 1 Data for Example 1. 17.8 32.3 43.4
21.3 33.5 44.5
23.8 34.9 47.0
25.9 36.6 48.8
27.4 38.5 52.5
29.4 39.7 61.4
30.6 41.2
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
Table 2 Estimated parameters for Example 1.
True MLM MSM EMSM
β
η
γ
μ
2.5
30 25.08
10 14.30
30.29 27.84
9.96 11.88
36.62 36.51 36.79 36.55
2.1652 2.2923 2.3766
εβ (%)
lnðLÞ
−75.367 −75.919 −75.446
εμ (%)
0.29 0.47 0.17
Table 3 Data for Example 2. 2 106 275
10 143 293
13 147 300+
23 173 300+
23 181 300+
28 212 300+
30 245 300+
65 247 300+
80 261 300+
88 266 300+
Table 4 Estimated parameters and goodness-of-fit measures for Example 2. Model
β
η
δ,α
MWD (2) EWD (3) MWED (4) GGD (5) BETA (9) DTM (10) TWD (11)
0.6537
633.56
408.28
5.4197 0.6888
447.96 186.14 186.15
0.1370
FSM (12)
γ,λ
0.6165 0.6190
653.25
0.4219 0.2223 354.30 423.90
0.7334 0.06674
111678 9.5118
459.39 452.35
0.6086 0.5433
SSE
p
−141.869
0.0289
−141.736 −142.033 −141.774 −141.408 −141.422 −141.680 −141.356
0.0386 0.0356 0.0262 0.0195 0.0188 0.0345 0.0288
0.2061 0.2116 0.1647 0.2597 0.4081 0.3787 0.2226 0.2211
lnðLÞ
the infinite support are Models (2)–(5); and the models with the finite support are Models (9)–(11). 5.1. Goodness-of-fit measures For a given dataset, we evaluate the goodness-of-fit of a fitted model using three measures. The first measure is the loglikelihood value lnðLÞ (not lnðLE Þ for the EMSM). The second measure is the sum of squared errors given by n
SSE ¼ ∑ di ½Fðt i Þ−F i 2
It is noted that the scale parameter of the truncated Weibull model given by (11) is very large. In this case, the model can be approximated by FðtÞ≈
13.39 8.31 4.94
ð31Þ
5.2.3. Comparison of models The models can be ranked based on each of lnðLÞ, SSE and p, and the results are shown in the 2–4th columns of Table 5. As seen, the proposed model is the best in terms of lnðLÞ; the model given by (10) is the best in terms of SSE; and the beta model given by (9) is the best in terms of p. This confirms the appropriateness of the proposed model. The median empirical CDF and the CDFs of the fitted models with the infinite and finite supports are shown in Figs. 1 and 2, respectively. As seen, the median empirical CDF is inverse S-shaped; and the CDFs of the models with the finite support are more similar to the median empirical CDF than those of the models with the infinite support in shape. Table 5 Rank numbers of the fitted models. Model
Rank no. for Example 2 Rank no. for Example 3 Average of ranks
MWD (2) EWD (3) MWED (4) GGD (5) BETA (9) DTM (10) TWD (11) FSM (12)
lnðLÞ
SSE
p
lnðLÞ
SSE
p
7 5 8 6 2 3 4 1
5 8 7 3 2 1 6 4
7 6 8 3 1 2 4 5
7 5 8 6 4 3 2 1
7 6 8 5 4 2 3 1
7 6 8 5 4 1 3 2
150
200
6.7 6.0 7.8 4.7 2.8 2.0 3.7 2.3
1
ð30Þ
F(t)
0.8 0.6 0.4 0.2 0
0
50
100
250
300
350
400
450
t Fig. 1. Median empirical CDF and the CDFs of the fitted models with the infinite support for Example 2.
5.2. Example 2
1
(9)
(11)
0.8
F(t)
5.2.1. Data The failure and running times of 30 devices shown in Table 3 come from Meeker and Escobar [33]. The sign “þ” indicates the running times. It is noted that the largest observation is rightcensored and there are two failures at t ¼ 23.
β ðt=ηÞβ t ¼ : γ ðγ=ηÞβ
This is the model given by (8), implying that the truncated Weibull model reduces into the inverse Pareto distribution in this case.
i¼1
where Fðt i Þ is the cdf of the fitted model and F i is the median empirical cdf at t i . Microsoft Excel provides a standard function to evaluate the median of the beta distribution (i.e., Betainv ð0:5; m; N−m þ 1Þ, where m ¼ Σ ij ¼ 1 dj ), and hence we do not need to use its approximation (e.g., see [30,31]). The third measure is the χ 2 statistic and corresponding p value. The critical value of p is usually in the range 0.1–0.3 (see, for example, [32]). Since the degree of freedom for computing the p value must at least be 1, the number of data intervals has to be larger than or equal to 5 for the three-parameter models under consideration.
47
0.6
(12)
0.4
(10)
0.2
5.2.2. Parameter estimation The MLEs exist for this example. To compute the p value, we divide the dataset into five intervals with an interval length of 75, and the last interval is half-open. The estimated parameter and the values of lnðLÞ, SSE and p are shown in Table 4.
0
0
50
100
150
200
250
300
350
400
450
t Fig. 2. Median empirical CDF and the CDFs of the fitted models with the finite support for Example 2.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
The bathtub curves of the models with the infinite and finite supports are shown in Figs. 3 and 4, respectively. As seen, the failure rates of the models with the finite support have faster increase rate in the wear-out phase than those of the models with the infinite support. 5.2.4. Some results derived from model (12) The proposed model is mathematically tractable, and hence it is easy to derive some important characteristics (e.g., moments, initial and minimal failure rates) of the devices. The results are shown in the second row of Table 6. The density function obtained from the fitted model is displayed in Fig. 5. As seen, the pdf decreases and the rate of decrease is also decreasing with time.
Table 6 Characteristics derived from Model (12). Example
μ
s
s=μ
rð0Þ,10−3
tm
rðt m Þ,10−3
2 3
191.40 173.41
138.74 135.87
0.7248 0.7836
9.23 9.04
85.31 99.65
3.43 3.98
0.01 0.008 0.006 f(t)
48
Example 3
0.004 0.002
5.3. Example 3
0
5.3.1. Data The data shown in Table 7 deal with failure times of electronic devices reported in Wang [34]. It is noted that the largest observation value is a failure observation. 5.3.2. Parameter estimation For the four models with the finite support, the MLEs exist for the model (10), and the MLM fails for the other three models. We use the EMSM to estimate their parameters. To further confirm the appropriateness of the EMSM, we fix the value of γ obtained from the EMSM and re-estimate the other two parameters using the MLM. To compute the p value, we divide the data into five intervals so that the degree of freedom equals 1; the interval length is taken as 81 so that the number of data in each interval is not smaller than 3. The estimated parameters and corresponding values of goodness-of-fit measures are shown in Table 8. To further confirm the appropriateness of the EMSM, we look at the goodness-of-fit measures of the fitted models (9), (11) and (12). Though the values of lnðLÞ obtained from the EMSM are slightly smaller than those obtained from the MLM (with γ being the same), the values of SSE are smaller and the values of p are larger than the corresponding values obtained from the MLM. This further validates the EMSM. 0.015 (3) (5) 0.01
0
100
200
300
400
500
t Fig. 5. Density functions derived from Model (12) for Examples 2 and 3.
Table 7 Data for Example 3. 5 98 245
11 122 293
21 145 321
31 165 330
46 195 350
75 224 420
5.3.3. Comparison of models The rank numbers of the fitted models are shown in the 5th– 7th columns of Table 5. As seen, the model proposed in this paper is the best in terms of both lnðLÞ and SSE; and the model given by (10) is the best in terms of p. This further confirms the appropriateness of the proposed model. The median empirical CDF and the CDFs of the fitted models with the infinite and finite supports are shown in Figs. 6 and 7, respectively. As seen, the CDFs of the models with the finite support are closer to the median empirical CDF than those of the models with the infinite support. The bathtub curves of the models with the infinite and finite supports are shown in Figs. 8 and 9, respectively. As seen, the failure rates of the models with the finite support have longer flat portion and faster increase rate in the wear-out phase than those of the models with the infinite support.
r(t)
(2)
5.3.4. Some results derived from model (12) Some important characteristics of the electronic devices are shown in the last row of Table 6. The density function is also displayed in Fig. 5.
(4)
0.005
0
Example 2
0
50
100
150
200
250
300
350
400
450
t
5.4. Summary
Fig. 3. Bathtub curves of the models with the infinite support for Example 2.
To integrate the results from the two examples, we calculate the average of the six rank numbers for each model and the results are shown in the last column of Table 5. In terms of the average rank number, the proposed model is the second best and the models with the finite support have smaller averages than the models with the infinite support.
0.015 (12)
(11)
(10)
0.01
r(t)
(9)
0.005
6. Discussion on the model support 0
0
50
100
150
200
250
300
350
400
450
t Fig. 4. Bathtub curves of the models with the finite support for Example 2.
The possible concerns with the model with the finite support are physical reasons and possible consequences of using such a model. We discuss these issues as follows.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
49
Table 8 Estimated parameters and goodness-of-fit measures for Example 3. β
MWD (2) EWD (3) MWED (4) GGD (5) BETA (9), EMSM BETA (9), MLM DTM (10) TWD (11), EMSM TWD (11), MLM FSM (12), EMSM FSM (12), MLM
0.6468 6.2271 0.7522
η
δ,α 665.71 373.00 134.04 144.95
0.7383 0.7793 0.6871 0.7549 0.8049 0.12073 0.1304
0.6012 1.2081 1.2781
410.34 1202.20 729.34 17.88 20.97
SSE
p
0.3413 0.1831 449.54 449.54 480.55 440.38 440.38 437.89 437.89
−108.933 −108.397 −109.117 -108.695 −108.162 −108.143 −108.145 -108.069 −108.055 −107.833 −107.830
0.0235 0.0231 0.0341 0.0184 0.0100 0.0135 0.0080 0.0084 0.0108 0.0062 0.0066
0.4175 0.4542 0.3471 0.4930 0.6375 0.5950 0.7446 0.6435 0.6160 0.6801 0.6750
5
0.8
4
0.6
3
0.4
0
0
50
100
150
200
250
300
350
400
450
Fig. 6. Median empirical CDF and the CDFs of the fitted models with the infinite support for Example 3.
1 0.8 0.6 0.4 0.2 0
Truncated Weibull
Weibull
2 1
t
F(t)
lnðLÞ
1
0.2
0
50
100
150
200
250
300
350
400
450
t Fig. 7. Median empirical CDF and the CDFs of the fitted models with the finite support for Example 3.
0.05 0.04
(3)
0.03
r(t)
γ,λ
276.88 0.1234
r(t)
F(t)
Model
(5) (2)
0.02
(4)
0.01 0
0
50
100
150
200
250
300
350
400
450
t Fig. 8. Bathtub curves of the models with the infinite support for Example 3. 0.05
(9)
0.04 (10)
r(t)
0.03 0.02 0.01 (12) 0
(11) 0
50
100
150
200
250
300
350
400
450
t
Fig. 9. Bathtub curves of the models with the finite support for Example 3.
0
t2 0
5
10
15 t
20
25
30
Fig. 10. Failure rate functions of the Weibull and truncated-Weibull distributions.
From a perspective of engineering, the model with the finite support is more appropriate than the model with the infinite support since the lifetime of a component is always finite; from a mathematical perspective, the model with the infinite support is more attractive since it avoids the problem to specify the upper bound of the lifetime. As such, the lifetime distribution with the infinite support can be viewed as an approximate representation of the lifetime. The selection between the model with the finite support and the model with the infinite support is problem-specific. In the context of bathtub failure rate modeling, we have shown above that the model with the finite support is preferred due to the favorite property that the failure rate can increase more quickly in the wear-out phase. When using the finite support model, the location parameter γ should be viewed as a point estimate of the upper bound of lifetime. When γ is appropriately estimated, the probability of event T≥γ should be small, but this does not imply that T≥γ is impossible. The value of γ provides useful information about the upper bound of lifetime, and can be applied in reliability and maintenance analysis. In the context of preventive maintenance, the component is at best preventively replaced at or before the second change point of the bathtub curve. The method to specify this change point is out of the scope of this paper. There is not a big gap between the model with the finite support and the model with infinite support. To illustrate, we consider the Weibull distribution Fðt; β; ηÞ with (β, η)¼ (3.5, 10.0). It can be well approximated by the truncated Weibull distribution given by (11) with γ ¼ 21:17 and the same values of β and η. This is because FðγÞ ¼ 1−10−6 so that the plots of their pdfs or cdfs are almost overlapped in (0, γ). However, the failure rate of the truncated-Weibull distribution quickly increases after some time point (e.g., t 2 ¼ 19:24, see Fig. 10). Generally, a distribution with the infinite support can be well approximated by its truncated distribution with the finite support as long as γ is appropriately specified; and obtaining an estimate of
50
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
the lifetime upper bound is more practical than assuming that the lifetime upper bound is infinite.
ðε1 þ ε2 þ ε1 ε2 Þr 0 ðtÞ≈−r 2 ðtÞ:
ðA:4Þ
Let β ¼ ε1 þ ε2 þ ε1 ε2 , which is positive and small. From (A.4) and after some simplifications, we have
7. Conclusions In this paper we developed a new three-parameter bathtub curve model with a finite support. Two examples have been used to illustrate its appropriateness. The new model serves as a good alternative for modeling bathtub-shaped failure rate data. An extended maximum spacing method has been developed to estimate the parameters of the models with a location parameter. The appropriateness of the method has been illustrated by two examples. We have compared the performances of the bathtub curve models with infinite and finite supports, and the results support the use of the models with the finite support. The other findings have been 1. The model given by (10) and the beta distribution are good alternatives for fitting the data with a bathtub failure rate. 2. The upper truncated Weibull model reduces into the inverse Pareto distribution when the scale parameter is very large. The results obtained in this paper are useful for the researchers who are attempting to develop new bathtub curve models; and for the practitioners who need to choose an appropriate bathtub curve model for fitting a given dataset that is considered to have a bathtub failure rate. A topic for the future study is to develop an appropriate method to specify the change points of the bathtub curve.
Acknowledgment The author would like to thank the two referees for their helpful comments and suggestions which have greatly enhanced the clarity of the paper. The research was supported by the National Natural Science Foundation (No. 71071026).
Appendix A. Derivation of the model (12) We use an approximate approach to derive the model (12) as follows. The starting point is the failure rate function given by rðtÞ ¼ f ðtÞ=RðtÞ:
ðA:1Þ
From (A.1), we have r′ðtÞ ¼ r 2 ðtÞ þ rðtÞ
Using (A.3) to (A.2) we have
f ′ðtÞ : f ðtÞ
ðA:2Þ
For a bathtub failure rate, rðtÞ has a minimum value r m at t ¼ t m . Using t m as a partition point, the bathtub curve can be divided into a decreasing part and an increasing part. We first consider the decreasing part. When t is small, we have:
rðtÞ 4 f ðtÞ but rðtÞ is close to f ðtÞ, and jr′ðtÞj o jf ′ðtÞj but r′ðtÞ is close to f ′ðtÞ. Let ε1 and ε2 be small positive numbers. The above two groups of relations can be written as follows: rðtÞ ≈1 þ ε1 ; f ′ðtÞ≈ð1 þ ε2 Þr′ðtÞ: f ðtÞ
ðA:3Þ
rðtÞ≈
β ; η 4 0: tþη
ðA:5Þ
We now consider the increasing part. For large t (t m ≪t o γ), we 0 assume that f ðtÞ decreases slowly so that f ðtÞ=f ðtÞ o0 and ≈0. In this case, (A.2) can be written as below r′ðtÞ≈r 2 ðtÞ:
ðA:6Þ
From (A.6) and after some simplifications, we have rðtÞ≈
1 γ−t
ðA:7Þ
where γ is a positive constant and can be viewed as a life upper limit. It is noted that (A.7) is actually the failure rate of the uniform distribution defined over (0; γ), which is an increasing function of t. Eq. (A.5) represents the early failure mode and Eq. (A.7) represents the wear-out failure mode. As such, the total failure rate is given by rðtÞ ¼
β 1 þ ; β; η; γ 40; t o γ: t þ η γ−t
ðA:8Þ
This is the model (12). References [1] Rajarshi S, Rajarshi MB. Bathtub distributions—a review. Communications in Statistics—Theory and Methods 1988;17(8):2597–621. [2] Lai CD, Xie M. Stochastic ageing and dependence for reliability. New York: Springer; 2006 [Chapter 3]. [3] Pham H, Lai CD. On recent generalizations of the Weibull distribution. IEEE Transactions on Reliability 2007;56:454–8. [4] Bebbington M, Lai CD, Zitikis R. Bathtub curves in reliability and beyond. Australian and New Zealand Journal of Statistics 2007;49:251–65. [5] Mudholkar GS, Asubonteng KO, Hutson AD. Transformation of the bathtub failure rate data in reliability for using Weibull-model analysis. Statistical Methodology 2009;6:622–33. [6] Crevecoeur GU. A model for the integrity assessment of ageing repairable systems. IEEE Transactions on Reliability 1993;42(1):148–55. [7] Zhang T, Xie M, Tang LC, Ng SH. Reliability and modeling of systems integrated with firmware and hardware. International Journal of Reliability, Quality and Safety Engineering 2005;12(3):227–39. [8] Clarotti C, Lannoy A, Odin S, Procaccia H. Detection of equipment aging and determination of the efficiency of a corrective measure. Reliability Engineering and System Safety 2004;84(1):57–64. [9] Bebbington M, Lai CD, Zitikis R. Useful periods for lifetime distributions with bathtub shaped failure rate functions. IEEE Transactions on Reliability 2006;55 (2):245–51. [10] Bebbington M, Lai CD, Zitikis R. Estimating the turning point of a bathtub shaped failure distribution. Journal of Statistical Planning and Inference 2008;138(4):1157–66. [11] Kim J. A change-point estimator with the hazard ratio. Journal of the Korean Statistical Society 2009;38(4):377–82. [12] Jiang R. Determination of two change points of a bathtub failure rate curve. In: Proceedings of the third international Asia conference on industrial engineering and management innovation. Beijing; 2012, p. 1116–20. [13] Jiang R. Specification of change points of failure rate or intensity function: a non-parametric approach. In: Proceedings of the IEEE international conference on industrial engineering and engineering management. Hong Kong; 2012, p. 747–51. [14] Cheng RC, Amin NAK. Estimating parameters in continuous univariate distributions with a shifted origin. Journal of Royal Statistical Society B 1983;45:394–403. [15] Ekström M. Alternatives to maximum likelihood estimation based on spacings and the Kullback–Leibler divergence. Journal of Statistical Planning and Inference 2008;138(6):1778–91. [16] Chen Z. A new two-parameter lifetime distribution with bathtub shape or increasing failure rate function. Statistics and Probability Letters 2000;49: 155–61. [17] Lai CD, Xie M. Murthy DNP. A modified Weibull distribution. IEEE Transactions on Reliability 2003;52(1):33–7. [18] Mudholkar GS, Srivastava DK. Exponentiated Weibull family for analyzing bathtub failure-rate data. IEEE Transactions on Reliability 1993;42(2):299–302.
R. Jiang / Reliability Engineering and System Safety 119 (2013) 44–51
[19] Jiang R, Murthy DNP. The exponentiated Weibull family: a graphical approach. IEEE Transactions on Reliability 1999;48(1):68–72. [20] Xie M, Tang Y, Goh TNA. Modified Weibull extension with bathtub-shaped failure rate function. Reliability Engineering and System Safety 2002;76: 279–85. [21] El-Gohary A, Alshamrani A, Al-Otaibi A. The generalized gompertz distribution. Applied Mathematical Modeling 2013;37(1-2):13–24. [22] Jiang R, Murthy DNP. Reliability modeling involving two Weibull distributions. Reliability Engineering and System Safety 1995;47:187–98. [23] Jiang R, Murthy DNP. Parametric study of competing risk model involving two Weibull distributions. International Journal of Reliability, Quality and Safety Engineering 1997;4(1):17–34. [24] Xie M, Lai CD. Reliability analysis using an additive Weibull model with bathtub-shaped failure rate function. Reliability Engineering and System Safety 1996;52(1):87–93. [25] Sarhan AM, Apaloo J. Exponentiated modified Weibull extension distribution. Reliability Engineering and System Safety 2013;112(4):137–44. [26] Zhang T, Xie M. On the upper truncated Weibull distribution and its reliability implications. Reliability Engineering and System Safety 2011;96:194–200.
51
[27] Cheng RCH, Traylor L. Non-regular maximum likelihood problems. Journal of the Royal Statistical Society: Series B 1995;57(1):3–44. [28] Titterington DM. Comment on estimating parameters in continuous univariate distribution. Journal of the Royal Statistical Society: Series B 1985;47:115–6. [29] Jiang R, Murthy DNP. Comments on: a general linear-regression analysis applied to the 3 parameter Weibull model. IEEE Transactions on Reliability 1997;46(3):389–93. [30] Hyndman RJ, Fan Y. Sample quantiles in statistical packages. The Amerian Statistician 1996;50:361–5. [31] Folland C, Anderson C. Estimating changing extremes using empirical ranking methods. Journal of Climate 2002;15:2954–60. [32] Blischke WR, Murthy DNP. Reliability: modeling, prediction, and optimization. New York: John Wiley & Sons; 2000. [33] Meeker WQ, Escobar LA. Statistical methods for reliability data. New York: John Wiley; 1998. [34] Wang FK. A new model with bathtub-shaped failure rate using an additive Burr XII distribution. Reliability Engineering and System Safety 2000;70(3):305–12.