# A New Finite Interval Lifetime Distribution Model for Fitting Bathtub-Shaped Failure Rate Curve

Authors: Xiaohong Wang; Chuang Yu; Yuxiang Li
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

This paper raised a new four-parameter fitting model to describe bathtub curve, which is widely used in research on components’ life analysis, then gave explanation of model parameters, and provided parameter estimation method as well as application examples utilizing some well-known lifetime data. By comparative analysis between the new model and some existing bathtub curve fitting model, we can find that the new fitting model is very convenient and its parameters are clear; moreover, this model is of universal applicability which is not only suitable for bathtub-shaped failure rate curves but also applicable for the constant, increasing, and decreasing failure rate curves.


---

## Body

## 1. Introduction

Typical failure rate curve is bathtub-shaped over product’s lifetime and has been widely accepted in the field of reliability and medicine. Bathtub curve can be divided into three phases: early failure, random failure, and wear-out failure phases. However, it is rather difficult to describe the model because the correlation between these three phases is not strong: both duration and intensity may be different; there may be even absence of one or two phases. In order to find a suitable bathtub curve fitting model, researchers have worked for decades. New bathtub curve fitting models were proposed one by one since the 1980s [1]. Nadarajah [2] and Jiang [3] had summarized some of these models, from which we can see that some models are on finite intervals and others are not. Models for finite interval have a good applicability which can fit the rapid increasing failure rate phase in wear-out phases well in short interval and can also fit the slowly increasing failure rate phase in long interval. On the other hand, models for infinite interval can be classified into addition type and extension type. The addition type [4, 5] is inconvenient in application because of its numerous parameters and sometimes the maximum likelihood estimation may not exist, while extension type [6–9] usually cannot fit the rapid increasing failure rate in wear-out phases well.Existing models for finite interval have two or three parameters; one parameter determines the interval, while the remaining one or two parameters describe the bathtub curve’s shape. Existing models could not match the bathtub curve very well because it has three phases which are generally not strongly related. Until now, it is difficult to find a bathtub-shaped model which can be widely accepted.Therefore, this paper put forward a four-parameter fitting model on finite interval for bathtub curve, gave explanation of model parameters in Section2, discussed parameter estimation method in Section 3, and then provided application examples utilizing this method in Section 4; in the end it summarized the model characteristic and application range in Section 5.


## 2. The Model

Both Beta distribution and Mudholkar’s DTM model [10] can be classified into bathtub model on finite interval, and failure rate function of Beta distribution is given by(1)
h
t
=
t
p
-
1
·
1
-
t
q
-
1
B
p
,
q
-
B
t
p
,
q
.Ghitany [11] showed that h
(
t
) is bathtub-shaped if p
<
1, and failure rate function of DTM through data transformation is(2)
h
t
=
α
σ
α
·
t
α
-
1
1
-
θ
t
α
+
1
.Whenα
<
1, the failure rate curve is bathtub-shaped. It is obvious that the core of these two models to obtain bathtub failure rate curve is t
m
·
(
1
-
t
)
n. For Beta distribution, the shape of bathtub curve can be regulated flexibly by choosing suitable parameters, while the DTM model owns interval parameters and does not contain any integral. Combining the advantages of the two methods, we proposed a new four-parameter fitting model whose reliability function and failure rate function are(3)
R
t
=
exp
⁡
-
λ
·
t
/
η
α
1
-
t
/
η
β
,
(4)
h
t
=
λ
t
/
η
α
-
1
η
1
-
t
/
η
β
+
1
·
β
-
α
·
t
η
+
α
,where α
>
0, β
≥
0, η
>
0, and λ
>
0; we call this “finite interval distribution model used in reliability engineering,” FIRE for short. FIRE is equivalent to DTM when α
=
β, and the three-parameter submodel when λ
=
1 also has good properties, which is called 3-FIRE.For FIRE,λ is the multiplication factor for failure rate function, which only influences the failure rate axis scale; η is the interval and scale parameter; the model interval is (
0
,
η
) when β
≠
0 and it influences the axis scale; α and β are both shape parameters; specifically, FIRE is equivalent to Weibull distribution when β
=
0. Figure 1 shows ten failure rate curves when λ
=
η
=
1, as η
=
1, so the range of t is from 0 to 1, and the shape parameters of each curve are listed in Table 1.Table 1
Number
1
2
3
4
5
6
7
8
9
10
α
4
2
1.1
1
0.95
0.7
0.7
0.3
0.7
0.7
β
0.1
0.1
0.1
0.1
0.1
0.4
0.1
0.1
0.002
0Figure 1Failure rate curves under different shape parameters are summarized into 6 types as in Table2; constant C is the minimum failure rate in curve type 1. By making the derivative of the failure rate function equal to zero we can get the time when failure rate is the smallest:(5)
t
=
η
-
η
α
2
α
=
β
,
t
=
η
β
-
α
·
α
β
β
-
α
+
1
-
α
α
≠
β
.Table 2
Type
α
β
Interval
Failure rate
Significance
1
(
0,1
)
(
0
,
+
∞
)
(
0
,
η
)
+
∞
→
C
→
+
∞
Whole bathtub
2
1
(
0
,
+
∞
)
(
0
,
η
)
Increase fromλ
/
η to +
∞
No early failures
3
(
1
,
+
∞
)
(
0
,
+
∞
)
(
0
,
η
)
Increase from 0 to+
∞
Wear-out failures only
4
(
0,1
)
0
(
0
,
+
∞
)
Decrease from+
∞ to 0
Early failures only
5
1
0
(
0
,
+
∞
)
Identicallyλ
/
η
Random failures only
6
(
1
,
+
∞
)
0
(
0
,
+
∞
)
Increase from 0 to+
∞
Wear-out failures onlySubstituting formula (5) into (4) we can get the constant C.Both type 3 and type 6 curves can only describe wear-out failure, the failure rate increases to infinite on a finite interval for type 3 while type 6 does not. For some kind of special bathtub curve, FIRE could not fit the situation without wear-out failure, but we can use the curve whenα
<
1, β
→
0 (such as curve 9 in Figure 1) to fit this situation approximately. There is no clear demarcation point between random failure phase from the other two phases, the time of duration for random failure is mainly influenced by parameter β, the random failure phase gets shorter as β get bigger, and it can be considered that there is no random failure phase when parameter β is very big.Parameterλ is a multiplication factor for failure rate function, but for reliability function, it is also a shape parameter just as α or β. With regard to model with two shape parameters, the reliability curves under different parameters should have at most one intersection excluding the two endpoints. With one more shape parameter, FIRE could form double intersection curve families as is shown in Figure 2. Such curve families have similar failure rate in both ends but are quite different in the middle part, which reflects FIRE’s good applicability.Figure 2When only small data is provided, there is no need to choose such complex model; people can use 3-FIRE model.


## 3. Parameter Estimation

### 3.1. Determination of Interval Parameter

Interval parameterη means the product maximum life. Being affected by working stress or environmental stress, the structure or composition of products under specified conditions will be changed inevitably, which means there is a life limit for any product, and lifetime limit comes when failure rate increases to +
∞.The value ofη should be estimated in accordance with components’ malfunction states considering that the maximum life is determined by η:(6)
η
=
t
n
+
t
n
-
t
n
-
k
n
k
.In formula (6), n is the sample size, t
n is the failure time of the nth sample, and k is the failure sample number at t
n. This formula has considered the following factors:(a)
ensuring that the maximum life is larger than the existing maximum lifetime value by adding a positive value tot
n;
(b)
the larger failure gap between the last several samples resulting in larger maximum life;
(c)
the larger sample size resulting in larger maximum life.Because maximum life is an inherent nature of the product, the interval parameter should only relate to failure data but not to the model we select, so we choose the same interval parameter for different finite interval models later in this paper. It should be noted that formula (6) is just an empirical formula from data analysis; how to better estimate maximum life still needs to be further studied.


### 3.2. Maximum Likelihood Estimate

Standard statistical techniques such as method of maximum likelihood can also be used in this case. The likelihood equations, given the complete or censored failure data set, can be written down and solved which besidesn is equal to one.Letx
=
t
/
η after determination of interval parameter; the failure rate function becomes(7)
f
t
=
R
t
·
h
t
=
λ
·
x
α
-
1
η
·
1
-
x
β
+
1
·
β
-
α
·
x
+
α
·
exp
⁡
-
λ
·
x
α
1
-
x
β
.Taking logarithm of the likelihood functionL
=
∏
i
=
1
n
f
(
t
i
),(8)
ln
⁡
L
α
,
β
,
λ
=
n
·
ln
⁡
λ
+
α
-
1
·
∑
i
=
1
n
ln
⁡
x
i
-
β
+
1
·
∑
i
=
1
n
ln
⁡
⁡
1
-
x
i
+
∑
i
=
1
n
ln
⁡
⁡
β
-
α
x
i
+
α
-
λ
·
∑
i
=
1
n
x
i
α
1
-
x
i
β
-
n
·
ln
⁡
η
.By taking the partial derivation of the above function we can obtain(9)
∂
ln
⁡
L
∂
λ
=
n
λ
-
∑
i
=
1
n
x
i
α
1
-
x
i
β
,
∂
ln
⁡
L
∂
α
=
∑
i
=
1
n
1
-
x
i
β
-
α
x
i
+
α
+
ln
⁡
x
i
-
λ
·
x
i
α
·
ln
⁡
x
i
1
-
x
i
β
,
∂
ln
⁡
L
∂
β
=
∑
i
=
1
n
λ
·
x
i
α
·
ln
⁡
1
-
x
i
1
-
x
i
β
x
i
β
-
α
x
i
+
α
a
a
a
a
a
a
-
ln
⁡
1
-
x
i
+
λ
·
x
i
α
·
ln
⁡
1
-
x
i
1
-
x
i
β
.Whenn
=
1, the equations have no solution.Whenn is not equal to one, it is difficult to solve the parameters through the method of making the partial derivative equal to zero directly; in this case numerical methods could be used, such as finding combinations of parameters which can get the maximum ln
⁡
L by searching algorithm and then substituting it into the partial derivation function to verify if the solution is correct.


## 3.1. Determination of Interval Parameter

Interval parameterη means the product maximum life. Being affected by working stress or environmental stress, the structure or composition of products under specified conditions will be changed inevitably, which means there is a life limit for any product, and lifetime limit comes when failure rate increases to +
∞.The value ofη should be estimated in accordance with components’ malfunction states considering that the maximum life is determined by η:(6)
η
=
t
n
+
t
n
-
t
n
-
k
n
k
.In formula (6), n is the sample size, t
n is the failure time of the nth sample, and k is the failure sample number at t
n. This formula has considered the following factors:(a)
ensuring that the maximum life is larger than the existing maximum lifetime value by adding a positive value tot
n;
(b)
the larger failure gap between the last several samples resulting in larger maximum life;
(c)
the larger sample size resulting in larger maximum life.Because maximum life is an inherent nature of the product, the interval parameter should only relate to failure data but not to the model we select, so we choose the same interval parameter for different finite interval models later in this paper. It should be noted that formula (6) is just an empirical formula from data analysis; how to better estimate maximum life still needs to be further studied.


## 3.2. Maximum Likelihood Estimate

Standard statistical techniques such as method of maximum likelihood can also be used in this case. The likelihood equations, given the complete or censored failure data set, can be written down and solved which besidesn is equal to one.Letx
=
t
/
η after determination of interval parameter; the failure rate function becomes(7)
f
t
=
R
t
·
h
t
=
λ
·
x
α
-
1
η
·
1
-
x
β
+
1
·
β
-
α
·
x
+
α
·
exp
⁡
-
λ
·
x
α
1
-
x
β
.Taking logarithm of the likelihood functionL
=
∏
i
=
1
n
f
(
t
i
),(8)
ln
⁡
L
α
,
β
,
λ
=
n
·
ln
⁡
λ
+
α
-
1
·
∑
i
=
1
n
ln
⁡
x
i
-
β
+
1
·
∑
i
=
1
n
ln
⁡
⁡
1
-
x
i
+
∑
i
=
1
n
ln
⁡
⁡
β
-
α
x
i
+
α
-
λ
·
∑
i
=
1
n
x
i
α
1
-
x
i
β
-
n
·
ln
⁡
η
.By taking the partial derivation of the above function we can obtain(9)
∂
ln
⁡
L
∂
λ
=
n
λ
-
∑
i
=
1
n
x
i
α
1
-
x
i
β
,
∂
ln
⁡
L
∂
α
=
∑
i
=
1
n
1
-
x
i
β
-
α
x
i
+
α
+
ln
⁡
x
i
-
λ
·
x
i
α
·
ln
⁡
x
i
1
-
x
i
β
,
∂
ln
⁡
L
∂
β
=
∑
i
=
1
n
λ
·
x
i
α
·
ln
⁡
1
-
x
i
1
-
x
i
β
x
i
β
-
α
x
i
+
α
a
a
a
a
a
a
-
ln
⁡
1
-
x
i
+
λ
·
x
i
α
·
ln
⁡
1
-
x
i
1
-
x
i
β
.Whenn
=
1, the equations have no solution.Whenn is not equal to one, it is difficult to solve the parameters through the method of making the partial derivative equal to zero directly; in this case numerical methods could be used, such as finding combinations of parameters which can get the maximum ln
⁡
L by searching algorithm and then substituting it into the partial derivation function to verify if the solution is correct.


## 4. Example and Analysis

For the sake of comparative analysis, some well-known failure data are fit both by FIRE and other models below:(a)
Mudholkar’s DTM model [10], which is three-parameter model on finite interval is the submodel of FIRE when α
=
β;
(b)
Jiang’s FSM model [3], which is a three-parameter model on finite interval;
(c)
Almalki’s NMW model [4], whose submodels include additive Weibull model, modified Weibull model, and S–Z modified Weibull model; it is a five-parameter model on infinite interval;
(d)
Xie et al.’s MWE model [12], whose submodels include Chen’s model; it is a three-parameter model on infinite interval;
(e)
Wang’s ABXII model [5], which is a six-parameter model on infinite interval;
(f)
Lemonte’s ENH model [7], which is a new three-parameter family of exponential-type distributions on infinite interval.Example 1.
Aarset’s 50 failure data [13]: 0.1, 0.2, 1, 1, 1, 1, 1, 2, 3, 6, 7, 11, 12, 18, 18, 18, 18, 18, 21, 32, 36, 40, 45, 46, 47, 50, 55, 60, 63, 63, 67, 67, 67, 67, 72, 75, 79, 82, 82, 83, 84, 84, 84, 85, 85, 85, 85, 85, 86, 86.Interval parameter is obtained by formula (6):(10)
η
=
86
+
86
-
85
50
×
2
=
86.01
.For DTM,θ
=
1
/
η
=
0.011627 and γ
=
η
=
86.01 for FSM.Other parameters are estimated by MLE, and ln(L), AIC value, and P value are utilized to evaluate the fitting results as in Table 3.Table 3
Example1.
Model
Parameter
ln
⁡
(
L
)
AIC
P value
FIRE
α
=
0.5384
,
β
=
0.1967
,
λ
=
0.8499
,
η
=
86.01
−202.18
412.36
0.79028
3-FIRE
α
=
0.5709
,
β
=
0.1763
,
η
=
86.01
−202.49
410.99
0.8639
DTM
α
=
0.2762
,
σ
=
688.414
,
θ
=
0.011627
−207.29
420.58
0.36173
FSM
β
=
0.093
,
γ
=
86.01
,
η
=
16.00
−221.93
449.85
0.03575
NMW
α
=
0.071
,
β
=
7.015
×
1
0
-
8
,
γ
=
0.016
,
θ
=
0.595
,
λ
=
0.197
−212.92
435.83
0.87463
MWE
α
=
13.8597
,
β
=
0.5890
,
λ
=
0.007557
−231.65
469.29
0.17501
ABXII
c
1
=
0.7016
,
s
1
=
424555
,
k
1
=
492.92
c
2
=
82.2343
,
s
2
=
92.43
,
k
2
=
1075.4
−206.10
424.20
0.79963
ENH
α
=
0.8200
,
β
=
0.7966
,
λ
=
0.0252
241.36
488.73
0.06385Example 2.
Wang’s 18 data [5]: 5, 11, 21, 31, 46, 75, 98, 122, 145, 165, 195, 224, 245, 293, 321, 330, 350, 420. Consider(11)
η
=
420
+
420
-
350
18
×
1
=
423.8889
.Example 3.
Coetzee’s 127 data [14]: 0.01, 0.01, 0.01, 0.01, 0.01, 0.02, 0.02, 0.02, 0.02, 0.03, 0.04, 0.06, 0.08, 0.1, 0.1, 0.12, 0.12, 0.12, 0.13, 0.14, 0.15, 0.15, 0.15, 0.16, 0.16, 0.17, 0.18, 0.18, 0.19, 0.2, 0.21, 0.22, 0.23, 0.25, 0.26, 0.28, 0.28, 0.3, 0.32, 0.34, 0.36, 0.38, 0.39, 0.41, 0.41, 0.42, 0.43, 0.44, 0.44, 0.45, 0.45, 0.5, 0.53, 0.56, 0.58, 0.58, 0.61, 0.62, 0.62, 0.62, 0.64, 0.66, 0.7, 0.7, 0.7, 0.72, 0.77, 0.78, 0.78, 0.8, 0.82, 0.83, 0.85, 0.86, 0.96, 0.97, 0.98, 0.99, 1.05, 1.06, 1.07, 1.18, 1.35, 1.36, 1.42, 1.55, 1.59, 1.65, 1.73, 1.77, 1.79, 1.8, 1.91, 2.09, 2.14, 2.15, 2.15, 2.31, 2.33, 2.36, 2.43, 2.45, 2.5, 2.51, 2.58, 2.64, 2.68, 3.08, 3.94, 4.12, 4.33, 4.42, 4.53, 4.88, 4.97, 5.11, 5.32, 5.55, 6.63, 6.89, 7.62, 11.41, 11.76, 11.85, 12.36, 13.22. Consider(12)
η
=
13.22
+
13.22
-
12.36
127
×
1
=
13.2268
.Examples1 and 2 are two sets of data which are most widely used in the study of bathtub curve model. Wang [5] has illustrated that the failure rate corresponding to these two datasets is bathtub-shaped through TTT transforming. The failure rate grows faster in wear-out phase in Example 1 while there is almost no wear failure in Example 3 (Table 5). According to the fitting results, FIRE has greater advantages for data in Example 1 and also works well for Examples 2 and 3.As mentioned earlier, extension-type model on infinite interval (e.g., MWE) is difficult to describe the rapid increase (Example1) failure rate, and the maximum likelihood estimation of addition type model on finite interval may not necessarily exist; NMW and ABXII are models acquired by adding the failure rate function using this method. In Example 2, the maximum likelihood estimates of both methods do not exist because when the NMW’s β
→
0
+ and when the ABXII’s s
2
→
42
0
+, the likelihood function L
→
+
∞.According to the parameters in Table4 we can get the density curve in Figure 3; from its partial enlarged figure (Figure 4) we can see that the density function of the fault value exceeded 6000 similar to the NMW; a partial enlarged figure (Figure 4) shows that the density function of the fault value exceeds 6000; this estimation result is obviously pointless. In contrast, by choosing the interval parameters, model on finite interval has avoided the situation that maximum likelihood estimation may not exist, as the FIRE density curve which is shown in Figure 3.Table 4
Example2.
Model
Parameter
ln
⁡
(
L
)
AIC
P value
FIRE
α
=
0.8810
,
β
=
0.1999
,
λ
=
1.6530
,
η
=
423.8889
−107.46
222.91
0.99970
3-FIRE
α
=
0.7113
,
β
=
0.3179
,
η
=
423.8889
−108.24
222.49
0.7251
DTM
σ
=
465.5967
,
θ
=
0.002359
,
α
=
0.4364
−110.15
226.29
0.99678
FSM
β
=
0.1207
,
γ
=
423.8889
,
η
=
20.4219
−107.43
220.87
0.999998
NMW
α
=
0.0178
,
β
=
1
×
1
0
-
200
,
γ
=
1.078
,
θ
=
0.782
,
λ
=
1.298
−107.38
224.76
—
MWE
α
=
0.2148
,
β
=
0.2465
,
λ
=
0.02180
−109.74
225.49
0.99345
ABXII
c
1
=
0.998
,
s
1
=
4.599
×
1
0
9
,
k
1
=
2.28
×
1
0
7
c
2
=
6.1
×
1
0
7
,
s
2
=
420
+
9.54
×
1
0
-
7
,
k
2
=
8.36
×
1
0
5
−94.66
201.33
—
ENH
α
=
0.5575
,
β
=
1.3087
,
λ
=
0.0191
111.8
229.60
0.96739Table 5
Example3.
Model
ln
⁡
(
L
)
AIC
FIRE
−173.62
355.23
3-FIRE
−231.75
469.50
DTM
−217.50
440.99
FSM
−173.63
353.26Figure 3Figure 4Comparison of the four models on finite interval showed that FIRE fit the 3 datasets the best judging by ln(L), and the AIC of FIRE is sometimes higher because of the one more parameter.The fitting results of the other models are not stable; therefore, FIRE with 4 parameters has an advantage over the other three methods considering model applicability.


## 5. Conclusions

In this paper, we proposed a new four-parameter fitting model on finite interval for bathtub curve, which can not only fit bathtub-shaped failure rate curve well, but is also applicable for the constant, increasing, and decreasing failure rate curves. This model has good applicability in life data fitting and has advantages of moderate parameters number, clear parameters significance, and convenience parameters estimation and is going to be widely used in reliability and life analysis.



---

*Source: 954327-2015-03-26.xml*