1 Introduction

# The Weibull Birnbaum-Saunders Distribution: Properties and Applications

The Weibull Birnbaum-Saunders Distribution: Properties and Applications

Lazhar Benkhelifa

Laboratory of Applied Mathematics, Mohamed Khider University, Biskra, 07000, Algeria

Departement of Mathematics and Informatics, Larbi Ben M’Hidi University, Oum El Bouaghi, 04000, Algeria

l.benkhelifa@yahoo.fr

Abstract

This paper introduces a new four-parameter lifetime model called the Weibull Birnbaum-Saunders distribution. This new distribution represents a more flexible model for the lifetime data. Its failure rate function can be increasing, decreasing, upside-down bathtub shaped, bathtub-shaped or modified bathtub shaped depending on its parameters. Some structural properties of the proposed model are investigated including expansions for the cumulative and density functions, moments, generating function, mean deviations, order statistics and reliability. The maximum likelihood estimation method is used to estimate the model parameters and the observed information matrix is determined. The flexibility of the new model is shown by means of two real data sets.

Keywords: Birnbaum-Saunders distribution, Weibull-G class, moment, order statistic, maximum likelihood estimation, observed information matrix.

## 1 Introduction

The two-parameter Birnbaum–Saunders (BS) distribution, which was introduced by Birnbaum and Saunders [[3]], is a very important lifetime model and is widely used in reliability studies. This distribution, also known as the fatigue life distribution, was originally derived from a model that shows the total time that passes until that some type of cumulative damage, produced by the development and growth of a dominant crack, surpasses a threshold value and causes the material specimen to fail. Desmond [[10]] provided a more general derivation based on a biological model and strengthened the physical justification for the use of this distribution.

A random variable having the BS distribution with shape parameter and scale parameter , denoted by BS is defined by

 T=β⎡⎣αZ2+√(αZ2)2+1⎤⎦2,

where is a standard normal random variable. The cumulative distribution function (cdf) of is

 G(t)=Φ(v), t>0, (1)

where is the standard normal distribution function, and The probability density function (pdf) corresponding to (1) is given by

 g(t)=κ(α,β)t−3/2(t+β)exp{−τ(t/β)2α2}, t>0, (2)

where and The fractional moments of are given by (see [[20]])

 E(Tk)=βkI(α,β),

where

 I(α,β)=Kk+1/2(α−2)+Kk−1/2(α−2)2K1/2(α−2), (3)

and the function denotes the modified Bessel function of the third kind with representing its order and the argument. The parameter is the median of the BS distribution, because

Since the BS distribution was proposed, it has received much attention in the literature. This attention for the BS distribution is due to its many attractive properties and its relation to the normal distribution. For more details on the BS distribution, we refer to [[15]] and the references therein. The BS distribution has been used in several research areas such as engineering, environmental sciences, finance, and wind energy. However, it allows for upside-down hazard rates only (see [[14]]), hence cannot provide reasonable fits for modeling phenomenon with decreasing, increasing, modified bathtub shaped and bathtub-shaped failure rates which are common in reliability studies.

For this reason, several generalizations and extensions of the BS distribution have been proposed in the literature. For example, Cordeiro and Lemonte [[6]], using the beta-G class [[11]], proposed an extension of BS distribution named as the beta BS distribution. Saulo et al. [[22]], based on the work of Cordeiro and de Castro [[5]], defined the Kumaraswamy BS distribution. Lemonte [[16]], based on the scheme introduced by Marshall and Olkin [[17]], deifined the Marshall–Olkin extended BS distribution. Cordeiro et al. [[7]] adopted the McDonald-G class [[2]] to define the McDonald BS distribution. Cordeiro et al. [[8]] used the generator approach of Zografos and Balakrishnan [[28]] to introduce the gamma Birnbaum-Saunders distribution. In this paper, a new four-parameter extension for the BS distribution is proposed.

Recently, Bourguignon et al. [[4]] proposed an interesting method of adding a new parameter to an existing distribution. The resulting distribution, known as the Weibull-G distribution, gives more flexibility to model various types of data. Let be a continuous baseline distribution with density depends on a parameter vector and the Weibull cdf (for ) with positive parameters and . Bourguignon et al. [[4]] replaced the argument by where , and defined the cdf of their class by

 F(t;a,b,θ)=ab∫[G(t,θ)¯¯¯G(t,θ)]0wb−1e−awbdw=1−exp⎛⎝−a[G(t,θ)¯G(t,θ)]b⎞⎠. (4)

Then, the Weibull-G density function is given by

 f(t;a,b,θ)=abg(t,θ)⎡⎣G(t,θ)b−1¯¯¯¯G(t,θ)b+1⎤⎦exp⎛⎜⎝−a⎡⎣G(t,θ)¯¯¯¯G(t,θ)⎤⎦b⎞⎟⎠. (5)

Some generalized distributions have been proposed under this methodology. Tahir et al. [[23], [24], [25]] defined the Weibull-Pareto, Weibull-Lomax and Weibull-Dagum distributions by taking to be the cdf of the Pareto, Lomax and Dagum distributions, respectively. More recently, Afify et al. [[1]] defined and studied the Weibull Fréchet distribution. In a similar way, we propose a new extension for the BS distribution called the Weibull BS (WBS) distribution, which has been applied to the modeling of fatigue failure times and reliability studies.

The rest of the paper is organized as follows. In Section 2, we introduce the WBS distribution and plot its density and failure rate functions. In Section 3, we provide a mixture representation for its density and cumulative distributions. Structural properties such as the ordinary moments, generating function, quantile function and simulation, mean deviations, the density of the order statistics and the reliability are derived in Section 4. In Section 5, we discuss maximum likelihood estimation of the WBS parameters and derive the observed information matrix. Two applications are presented in Section 6 to show the potentiality of the new distribution. Some concluding remarks are given in Section 7.

## 2 The WBS distribution

Substituting (1) in (4), the cdf of the WBS distribution can be written as

 F(t)=1−exp⎛⎝−a[Φ(v)1−Φ(v)]b⎞⎠. (6)

The pdf corresponding to (6) is given by

 f(t) =abκ(α,β)t−3/2(t+β)exp(−τ(t/β)2α2) ×[Φ(v)b−1{1−Φ(v)}b+1]exp⎛⎝−a[Φ(v)1−Φ(v)]b⎞⎠, (7)

where is a scale parameter and , and are positive shape parameters. It is clear that the BS distribution is not a special case of WBS distribution. If a random variable follows a WBS distribution with parameters and will be denoted by WBS The reliability and the failure rate function of are, respectively, given by

 R(t)=exp⎛⎝−a[Φ(v)1−Φ(v)]b⎞⎠,

and

 h(t)=abκ(α,β)t−3/2(t+β)[Φb−1(v){1−Φ(v)}b+1]exp{−τ(t/β)2α2}.

Plots of pdf and failure rate function of the WBS distribution for selected values of the parameters are given in Figure 1 and Figure 2, respectively. Figure 1 indicates that the WBS pdf can take various shapes such as symmetric, right-skewed and left-skewed depending on the parameter values. Figure 2 shows that the failure rate function of the WBS distribution can be increasing, decreasing, upside-down bathtub (unimodal) shaped, bathtub-shaped or modified bathtub shaped (unimodal shape followed by increasing) depending on the parameter values. So, the WBS distribution is quite flexible and can be used effectively in analyzing survival data.

Figure 1. Plots of the WBS pdf for some values of the parameters.
Figure 2. Plots of the WBS failure rate function for some values of the parameters.

## 3 Mixture representation

In this section, we derive expansions for the pdf and cdf of the WBS distribution. The pdf and cdf of the WBS distribution can be written as a linear combination of the pdf and cdf of exponentiated BS (EBS) distribution, respectively. A random variable having the EBS distribution with parameters , and , denoted by EBS, if its cdf and pdf are given by and , respectively, where is defined in (1) and is given in (2). The properties of exponentiated distributions have been studied by several authors. For example, Mudholkar and Srivastava [[19]] studied the exponentiated Weibull distribution, Gupta and Kundu [[13]] studied the exponentiated exponential distribution and Sarhan and Apaloo [[21]] proposed the exponentiated modified Weibull extension distribution.

The pdf of the WBS distribution (7) can be written as

 f(t)=abg(t)Φb−1(v)[1−Φ(v)]b+1exp⎛⎝−a[Φ(v)1−Φ(v)]b⎞⎠. (8)

Using the series expansion for the exponential function, we obtain

 exp⎛⎝−a[Φ(v)1−Φ(v)]b⎞⎠=∞∑k=0(−1)kakk!Φbk(v)[1−Φ(v)]bk. (9)

Substituting (9) in (8), we get

 f(t)=abg(t)∞∑k=0(−1)kakk!Φbk+b−1(v)[1−Φ(v)]−(bk+b+1).

Since , for and , then by using the binomial series expansion given by

 [1−Φ(v)]−(bk+b+1)=∞∑j=0Γ(bk+b+1+j)j!Γ(bk+b+1)Φj(v),

where is the complete gamma function, we obtain

 f(t)=∞∑k,j=0wk,jhbk+b+j(t), (10)

where

 wk,j=(−1)kbak+1Γ(bk+b+1+j)k!j!(bk+b+j)Γ(bk+b+1),

and denotes the EBS density function. By integrating (10), we get

 F(t)=∞∑k,j=0wk,jΦbk+b+j(v). (11)

It easy clear that Equation (10) means that the pdf of the WBS distribution is a double linear mixture of the pdf of EBS distribution. Based on this equation, several structural properties of the WBS distribution can be obtained by knowing those of the EBS distribution. For example, the ordinary, inverse and factorial moments, generating function and characteristic function of the WBS distribution can be obtained directly from the EBS distribution.

If is a positive real non-integer, we can expand as

 Φbk+b+j(v)=∞∑r=0sr(bk+b+j)Φr(v), (12)

where

 Unknown environment '%

Thus, from equations (2), (10) and (12), we get

 f(t)=g(t)∞∑r=0drΦr(v), (13)

where

 dr=∞∑k,j=0wk,jsr(bk+b+j).

## 4 Some structural properties

In this section, we give some mathematical properties of the WBS distribution.

### 4.1 Moments

In this subsection, we derive the expression for th order moment of WBS distribution. The moments of some orders will help in determining the expected life time of adevice and also the dispersion, skewness and kurtosis in a given set of observations arising in reliability applications. The th moment of the WBS random variable can be derived from the probability weighted moments of the BS distribution. The probability weighted moments of the BS distribution are formally defined, for and non-negative integers, by

 τp,r=∫∞0tpg(t)Φr(v)dt. (14)

There are many softwares such as MAPLE, MATLAB and R that can be used to compute the integral (14) numerically. From [[6]], we have an alternative representation to compute that is

 τp,r =βp2rr∑j=1(rj)∞∑k1,…,kj=0A(k1,…,kj) ×2sj+j∑m=0(−1)m(2sj+jm)I(p+2sj+j−2m2,α), (15)

where ,

and is calculated from (3) in terms of the modified Bessel function of the third kind.

Therefore, the th moment of can be written from (13) as

 μ′s=∞∑r=0drτs,r, (16)

where is obtained from (15) and is given by (13). We can compute numerically the th moment in any symbolic software by taking in the sum a large number of summands instead of infinity.

### 4.2 Moment generating function

Let WBS. The moment generating function of , say , is an alternative specification of its probability distribution. Here, we provide two representations for . From (13), we obtain a first representation

 M(s)=∞∑r=0dr∫∞0estg(t)Φr(v)dt.

Expanding the exponential function, we can rewritte as

 M(s)=∞∑r=0∞∑p=0drτp,rp!sp.

The second representation for is based on the quantile expansion of the BS distribution. From (10), we have

 M(s)=∞∑k,j=0wk,j∫∞0estg(t)Φbk+b+j−1(v)dt,

where is the BS(,) pdf. By setting in above integral, we get

 M(s)=∞∑k,j=0wk,j∫10ubk+b+j−1exp(sQ(u))du,

where is the quantile function of the BS distribution and is given by (1). Using the exponential expansion, we get

 M(s)=∞∑k,j=0wk,j∞∑p=0spp!∫10ubk+b+j−1Qp(u)du. (17)

From [[6]], if the condition holds, we have the expansion for the quantile function of the BS distribution

 t=Q(u)=∞∑i=0ηi(u−1/2)i, (18)

where

 ηi=(2π)i/2∞∑j=0djej,i,

for for and the quantities can be determined from the recurrence equation

 ej,i=(ia0)−1i∑m=1(mj+m−i)amcj,i−m,

and . Here, the quantities are defined by (for ) and (for ), where the ’s are computed recursively from

 Missing or unrecognized delimiter for \left

The first constants are We have .

Inserting (18) in (17), we get

 M(s)=∞∑k,j=0wk,j∞∑p=0spp!∫10ubk+b+j−1(∞∑i=0ηi(u−1/2)i)pdu. (19)

From [[12], Sec. 0.314] for a power series raised to a positive integer we have

 (∞∑i=0ηi(u−1/2)i)p=∞∑i=0δp,i(u−1/2)i,

where the coefficients (for ) can be determined from the recurrence equation

 δp,i=(ia0)−1i∑m=1(mp+m−i)amδj,i−m,

and . Hence, comes directly from and, therefore, from . Then

 M(s)=∞∑k,j,p,i=0wk,jspp!δp,i∫10ub(k+1)+j−1(u−1/2)idu. (20)

Therefore, using the binomial expansion in (20), we obtain

 M(s)=∞∑k,j,i,p=0i∑l=0(−1)i−l(il)wk,jδp,ip!(bk+b+j+l)2i−lsp.

### 4.3 Quantile function and simulation

In this subsection, we give an expression for WBS quantile function, in terms of the BS quantile function . The BS quantile function is straightforward computed from the standard normal quantile function . We have (see , [[6]])

 QBS(u)=β2(αΦ−1(u)+√4+[αΦ−1(u)]2)2.

Then, by inverting , we obtain

 QWBS(u)=β2(αΦ−1(p)+√4+[αΦ−1(p)]2)2,

where

 p=(−1aln(1−u))1b1+(−1aln(1−u))1b.

Therefore, it is easy to simulate the WBS distribution. Let be a continuous uniform variable on the unit interval . Thus, using the inverse transformation method, the random variable given by

 T=QWBS(U)=β2(αΦ−1(P)+√4+[αΦ−1(P)]2)2, (21)

where

 P=(−1aln(1−U))1b1+(−1aln(1−U))1b,

has the WBS distribution. Equation (21) may be used to generate random numbers from the WBS distribution when the parameters are known.

We plot the exact and the empirical cdf of WBS distribution in Figure 3 using a pseudo random sample of size1000 to check the correctness of the procedure for simulating a data set from WBS distribution. The histograms for two generated data sets and the exact WBS density plots from two simulated data sets for some parameter values are given in Figure 4. These plots indicate that the simulated values are consistent with the WBS distribution.

Figure 3. Comparison of exact and empirical cdf of the WBS distribution to simulate random numbers.

Figure 4. Plots of the WBS densities for simulated data sets: (a) , and ; (b) , and .

### 4.4 Mean deviations

Let be a random variable having the WBS() distribution. The mean deviations of about the mean and about the median can be used as measures of spread in a population. They are given by

 δ1=E(∣∣T−μ′1∣∣)=∫∞0∣∣T−μ′1∣∣f(t)dt

and

 δ2=E(|T−m|)=∫∞0|t−m|f(t)dt,

respectively, where the mean is calculated from (16) and the median is given by The measures and   can be expressed as

 δ1=2μ′1F(μ′1)−J(μ′1) \ \ \ and \ \ \ \ δ2=E(|T−m|)=μ′1−2J(m),

where . From (13), can be written as

 J(q)=∞∑r=0drφ(q,r), (22)

where

 φ(q,r)=∫q0tg(t)Φr(v)dt.

From [[6]], we have

 φ(q,r) =κ(α,β)2rr∑j=1(rj)∞∑k1,…,kj=0β−(2sj+j)/2A(k1,…,kj)2sj+j∑m=0(−β)m Missing or unrecognized delimiter for \right

where and are defined in (15). Consider

 D(p,q)=∫q0tpexp{−(t/β+β/t)2α2}dt=βp+1∫q/β0upexp{−u+u−12α2}du.

From [[26]], we can write

 D(p,q)=2βp+1Kp+1(α−2)−qp+1Kp+1(q2α2β,β2α2q),

where, is the incomplete Bessel function with order and arguments and . Then, we obtain

 φ(q,r) =κ(α,β)2rr∑j=1(rj)∞∑k1,…,kj=0β−(2sj+j)/2A(k1,…,kj)2sj+j∑m=0(−β)m ×(2sj+jm){D(2sj+j−2m+12,q)+βD(2sj+j−2m−12,q)},

which can be calculated from the function . Hence, we can use this expression for to compute . From (22), we obtain the Bonferroni and Lorenz curves defined by and , respectively. These curves have applications in reliability.

### 4.5 Reliability

In the stress-strength modelling, is a measure of component reliability when it is subjected to random stress and has strength . The component fails at the instant that the stress applied to it exceeds the strength and the component will function satisfactorily whenever . The parameter is referred to as the reliability parameter. This type of functional can be of practical importance in many applications. In this Section, we derive the reliability when and have independent WBS and WBS distributions. The pdf of and the cdf of can be obtained from (10) and (11) as

 f1(t)=g(t)∞∑k,j=0w∗1k,jΦb1k+b1+j−1,

and

 F2(t)=∞∑l,m=0w2l,mΦb2l+b2+m(v)

respectively, where

 w∗1k,j=(−1)kb1ak+11Γ(b1k+b1+1+j)k!j!Γ(b1k+b1+1).

and

 w2l,m=(−1)lb2al+12Γ(b2l+b2+1+m)l!m!(b2l+b2+m)Γ(b2l+b2+1).

We have

 R=∫∞0f1(t)F2(t)dt.

Then

 R=∞∑k,j,l,m=0w∗1k,jw2l,m∫∞0g(t)Φd∗(v)dt,

where

 d∗=b1(k+1)+b2(l+1)+j+m−1.

From (12), we can write

 Φd∗(v)=∞∑r=0sr(δ)Φr(v),

and then, we get

 R=∞∑k,j,l,m=0w∗1k,jw2l,m∞∑r=0sr(δ)τ0,r,

where is obtained from (15).

### 4.6 Order statistics

In this section, the distribution of the th order statistic for the WBS distribution are presented. The order statistics play an important role in reliability and life testing. Let be a simple random sample from WBS distribution with cdf and pdf as in (7) and (8), respectively. Let denote the order statistics obtained from this sample. In reliability literature, the th order statistic, say , denotes the lifetime of an –out–of– system which consists of independent and identically components.  The pdf of is given by

 fi,n(t)=n!(n−i)!(n−1)!n−i∑j=0(−1)j(n−ij)f(t)F(t)i+j−1. (23)

From (7), we have

 Fi+j−1(t)=⎡⎣1−exp⎛⎝−a[Φ(v)1−Φ(v)]b⎞⎠⎤⎦i+j−1.

Using the binomial series expansion, we get

 Fi+j−1(t)=∞∑k=0(−1)k(i+j−1k)exp⎛⎝−ka[Φ(v)1−Φ(v)]b⎞⎠. (24)

Inserting (8) and (24) in (23), we obtain

 fi,n(t) Missing or unrecognized delimiter for \left ×n−k∑j=0∞∑k=0(−1)k+j(i+j−1k)(n−ij)exp⎛⎝−a(k+1)[Φ(v)1−Φ(v)]b⎞⎠.

Using the power series for the exponential function, we have

 exp⎛⎝−a(k+1)[Φ(v)1−Φ(v)]b⎞⎠=∞∑l=0(−1)l(ak+a)ll!Φ(v)bl[1−Φ(v)]bl.

Then

 fi,n(t) =g(t)n!ab(n−i)!(n−1)!n−k∑j=0∞∑k=0(−1)k+j(i+j−1k)(n−ij) ×∞∑l=0(−1)l(ak+a)ll!Φ(v)bl+b−1[1−Φ(v)]−(bl+b+1).

Since , for and , we have

 [1−Φ(v)]−(bl+b+1)=∞∑m=0Γ(bl+b+1+m)m!Γ(bl+b+1)Φm(v).

Therefore, the pdf of the th order statistic for WBS distribution is

 fi,n(t)=∞∑l,m=0ϑl,mhbl+b+m(t), (25)

where

 ϑl,m=∞∑k=0n−i∑j=0(−1)k+j+ln!bak+1(k+1)lΓ(bl+b+1+m)l!m!(n−i)!(n−1)!(bl+b+m)Γ(bl+b+1)(i+j−1k)(n−ij),

and is the EBS density function with power parameter . Equation (25) means that the density function of the WBS order statistics is a linear mixture of the EBS densities. Then, we can easily obtain the mathematical properties for . For example, the th moment of is

 E(Tpi,n)=∞∑l,m=0ϑl,m(bl+l+m)τp,(bl+l+m−1).

## 5 Estimation

In this section, we consider estimation of the unknown parameters of the WBS distribution by the method of maximum likelihood. Let be observed values of , independent random variables having the WBS distribution with unknown parameter vector . The total log-likelihood function for , is given by

 ℓ =ℓ(ξ)=nlog(a)+log(b)+log[κ(α,β)]−32n∑i=1log(ti)+n∑i=1log(ti+β) −12α2n∑i=1τ(tiβ)+(b−1)n∑i=1log[Φ(vi)]−(b+1)n∑i=1log[1−Φ(vi)] −an∑i=1[Φ(vi)1−Φ(vi)]b.

Then, the components of the unit score vector are given by

 ∂ℓ∂α =−nα(1+2α2)+1α3n∑i=1(tiβ+βti)−(b−1)αn∑i=1viϕ(vi)Φ(vi) −(b+1)αn∑i=1viϕ(vi)1−Φ(vi)+abαn∑i=1viϕ(vi)Φ(vi)b−1[1−Φ(vi)]b+1, ∂ℓ∂β =−n2β+n∑i=11ti+β+12βα2n∑i=1(tiβ−βti)−(b−1)2βαn∑i=1τ(√ti/β)ϕ(vi)Φ(vi) Missing or unrecognized delimiter for \right ∂ℓ∂a =n