A New Life Time Distribution: Burr III Modified Weibull Distribution And Its Application in Burn in Process Текст научной статьи по специальности «Математика»

A New Life Time Distribution: Burr III Modified Weibull Distribution And Its Application in Burn in Process

Deepthy G S1, Nicy Sebastian1

department of Statistics, St. Thomas College (Autonomous) Thrissur, Kerala - 680 001 (India)

deepthygs@gmail.com,nicycms@gmail.com

Abstract

In burn-in analysis, models with a bathtub-shaped hazard rate and a bimodal density function are inevitable.This work focusses on a new five parameter distribution called Burr III Modified Weibull distribution which can be used to design burn-in procedures and preventative maintenance for incurable devices. The statistical properties such as quantile function, hazard rate function and order statistics have been discussed.The model parameters are estimated using the maximum likelihood estimation technique, and the performance of the proposed model is evaluated using the simulation technique. Finally, a real data set is presented to demonstrate the model's utility and its application in the burn-in process.

Keywords: Burr III distribution, Modified Weibull distribution, maximum likelihood estimation.

1. Introduction

In lifetime analysis, one is often interested to know about the reasons for the failure of a system or component. The different ways through which a system or product may fail is known as competing risk. Some normal reasons for the failures are material defects, imperfection of manufacturing process, wear out processes etc. Here comes the importance of burn-in process in reliability. Burn-in is an important technique for detecting and eliminating early fault, thereby increasing system reliability. It is used to find the defective units before they reach the customers. Items that survive burn-in period ensure the product quality. Burn-in must have a high failure rate in the early stages of life in order to be successful. This characteristic is mainly found in a class of life time distributions with bathtub-shaped failure rates. In the modelling of the lifetime of a system, distributions with bathtub-shaped failure rates are preferred. The bathtub curve has three phases: an infant mortality phase (decreasing failure rate), a normal life phase (constant failure rate), and a wear-out phase (increasing failure rate). The necessary condition in distribution of lifetime of system is bathtub hazard rate and sufficient condition is bimodal density function. That is, a model with bathtub hazard rate function and a bimodal density function allows designing burn-in procedures and preventive maintenance of malfunctioning devices. The goal of this study is to create a new reliability distribution that meets both the necessary and sufficient conditions of the lifetime of system. The important works related to this are reliability enhancement through optimal burn-in by Way K [22], a new model in relation to burn-in of components and its consequences are discussed by M^ltoft [8] , the burn-in problems by Lawrence [12] and Chandrasekaran [5]. Park [17] learned about the impact of burn-in on the average residual life.

In reliability analysis, lifetime distriburions has an importance because of its ability to explain the nature of data and its properties. For analysing lifetime data, one often consider Weibull distribution. But one of the disadvantage of Weibull distribution is that it possess only monotone hazard rates not non-monotone hazard rates. Hence it cannot be used in the case of lifetime data possessing bathtub hazard rates. To overcome this problem, several new models were

proposed such as three parameter modified Weibull distribution introduced by Lai et al. [11], a four-parameter generalisation of the Weibull distribution that can simulate a bathtubshaped hazard rate function by Carrasco et al. [4], a modification is given to weibull distribution which exhibit non monotone behavior by Sarhan et al. [19], beta modified Weibull distribution introduced by Silva et al. [20], extended flexible Weibull distribution by Bebbington et.al. [1]. Khan [15] introduced another flexible distribution called modified beta Weibull distribution etc. Among which modified Weibull proposed by Lai et al. [11] seems to be more flexible than other proposed ones. Another important life time distribution that we consider in our study is Burr III distribution. In 1942, Burr introduced twelve types of cumulative distribution function based on the Pearson system of distributions among which Burr XII and Burr III are commonly used. The Burr type III distribution is used in a variety of domains, including survival and reliability research, forestry, and environmental studies etc. BurrIII is the inverse of BurrXII distribution. Burr III distribution is also known as dagum distribution studies on the wealth of distribition [6] and as kappa distribution in the meteorological literature [14]. Many modification of Burr III have been already intoduced such as Burr and Cislak [3], Johnson et al. [9] etc. This distribution has a wide range of applications in statistical modelling, including forestry by Gove et al. [7], meteorological field Mielke [14] and actuarial literature Kleiber [10] etc.

Using the concept of competing risk model, a new five parameter distribution called Burr III Modified Weibull (BIIIMW) is proposed in this paper. If we model lifetime of units subjecting to two risk as series system, then the life time of the observed unit is the minimum of the individual potential lifetimes associated with each risk (see [23]). That is the realibility function of the BIIIMW model is the product of the reliability function of Burr III and Modified Weibull distribution. This model exhibits both bathtub hazard function and bimodal density function, which are commonly present when dealing with survival and lifetime data, and it can also be utilised in burn-in procedures.

The cumulative distribution function (cdf) of the modified Weibull (MW) distribution proposed by Lai is given by,

FMW(x; a, ft, A) = 1 - e-axftA; x > 0, a > 0, ft > 0, A > 0. (1)

where a is the scale parameter, ft is the shape parameter and A is the accelerating factor in the time of imperfection and a factor of fragility in the individual's survival as time increases. The cumulative distribution (cdf) of the Burr III distribution is given by,

GB(x; c,k) = (1 + x-c)-k; x > 0, k > 0, c > 0. (2)

where c and k are shape parameters.

2. Proposed Distribution

The concept of competing risk model is used to create the new Burr III Modified Weibull distribution (BIIIMW). Realibility function of the model is the product of reliability functions of Burr III and Modified Weibull distributions. The realibility function, cumulative distribution

function and probability density function of the Burr III modified Weibull distribution (BIIIMW) are respectively given by,

r(x, c, k, a, ft, A) = e-axfteAx (1 - (1 + x-c)"k) x > 0, c > 0, k > 0, a > 0, ft > 0, A > 0,

F(x, c, k, a, ft, A) = 1 - e-axfteAx (1 - (1 + x-c)-k) , x > 0, c > 0, k > 0, a > 0, ft > 0, A > 0

Deepthy G S, Nicy Sebastian RT&A, No 1 (67)

On Burr III Modified Weibull Distribution Volume 17, March 2022

and

e—axP eAx

ckx

-c-1 (1 + x-c)-k-1 + (1 - (1 + x-c)-k) axp-1 eAx(p + Ax)

f (x c,K a p A) = ^ for x > 0, c > 0, k > 0, a > 0, p > 0, A > 0.

0, otherwise.

(4)

Where a is the scale parameter. p, c, k are the shape parameters and A is the accelerating factor in the imperfection time and a factor of fragility increases.

Plot of the density function of BIIIMW distribution is given in Figure 1. The figures shows that

Figure 1: Probability density functions of the BIIIMW distribution for different parameteric values.

the the BIIIMW distribution can be decresing, positively skewed, negatively skewed, unimodal and bimodal for selected values of parameters.

The particular case of the BIIIMW distribution is included in Table 1.

Table 1: Special cases of BIIIMW distribution

Parameters Distribution

a = 0 BurrIII

c=1, a = 0 Inverse Lomax

A = 0 BurrIII-Weibull

A = 0, p = 1 BurrIII-Exponential

A = 0, p = 2 BurrIII-Rayleigh

c=1 Inverse Lomax-Modified Weibull

c=1, A=0, p=1 Inverse Lomax-Exponential

3. Reliability Analysis

The hazard and reverse hazard function of the BIIIMW distribution is obtained respectively, as follows

-c-In I >v-e\-k-l , (t it I ,v-c\~k\ \nkxl

h(x)

ckx

(1 + x-c)-k-1 + (l - (1 + x-c)-k^j ax?-1 eAx(? + Ax) (l - (1 + x-c)-k)

and

t(x)

e-ax? eAx [ckx-

c1

(1 + x-c)-k-1 + (1 - (1 + x-c)-k^j ax?-1 eAx(ft + Ax)]

1e

x? eAx

1 - (1 + x-c)-

(5)

(6)

Plot of the hazard function for selected values of BIIIMW parameters are shown in Figure (2). The shape of the hazard function shows that BIIIMW distribution can accommodate both monotone and non monotone behavior such as monotonically decreasing, monotonically increasing, unimodal and bathtub shapes for different values of the parameters, which are more likely to be meet when dealing with survival and lifetime data.

Figure 2: Plot for hazard rate functions of BIIIMW distribution.

4. The Statistical Properties

In this section, some statistical properties of BIIIMW distribution such as quantile function and

order statistics are discussed.

4.1. Quantile Function

The quantile function of the BIIIMW distribution is obtained by solving the non-linear equation, that is FBIIIMW(x) = u,0 < u < 1.

1 - e-ax?eAx (1 - (1 + x-c)-k) = u, 1 - (1 + x-c)-k) = 1 - u.

k

Taking logarithm on both sides we have,

-c_<YPA*--

log (l - (1 + x-c) *) - axpeAx - log(1 - u)

0.

(7)

The quantile values of the BIIIMW distribution can be obtained by solving the equation (7) using numerical methods, where u denotes a uniformly distributed random variable on the interval [0,1]. Table 2 describe the quantile values of BIIIMW distribution for certain given parameter values.

Table 2: BIIIMW quantile values for selected parameters

(c,k,a,p,A,)

u (3,2,0.7,0.4,2) (3,2,1,2,2) (0.7,1,0.5,2,4) (1,0.5,3,4,3) (0.7,0.4,0.4,0.6,0.8)

0.1 0.00844 0.25199 0.04266 0.01010 0.000244

0.2 0.04566 0.33641 0.12545 0.04167 0.00274

0.3 0.10799 0.39880 0.22052 0.09871 0.01134

0.4 0.18256 0.45172 0.30147 0.18620 0.03183

0.5 0.26241 0.50027 0.36646 0.28989 0.07354

0.6 0.34600 0.54755 0.42171 0.37356 0.15240

0.7 0.43513 0.59646 0.47254 0.43637 0.29609

0.8 0.53557 0.65143 0.52404 0.49053 0.55078

0.9 0.66521 0.72343 0.58525 0.54787 0.99652

4.2. Order Statistics

Order statistics have a wide application in realibiity. Let X\, X2,..., Xn be a simple random sample from BIIIMW distribution with cdf and pdf given in (3) and (4), respectively. Let X(!) < X(2) < .... < X(n) denote the order statistics. Then the pdf of ith order statistc is given by,

fx,(Y)

n! f (y)

(i - 1)!(n - i)!

-[F(x)]i-1 [1 - F(Y)]'

i-1

fx (Y) = (i - 1)n(n - i)! e-<xxYex [ckx-c-1 (1 + Y-c )-k-1 + (1 - (1 + Y-c )-k) axp-1 eAx (p + Ay)]

* [1 - e-axPA (1 - (1 + Y-c)-k) The pdf of 1st order statistic X1 is given by,

e-axPA (1 - (1 + Y-c)-k

fx1 (y)

ne-axPeAx [ckx-c-1 (1 + Y-c)-k-1 + (1 - (1 + Y-c<yP-1 eAY(p + Ay)]

e-aYpeAx (1 - (1 + x-c)-k

n—1

The pdf of n order statistic Xn is given by,

c) k 1

fxn(x)

ne-axPeA% [ckx-c-1 (1 + x-c)-k-1 + (1 - (1 + x-c)-k) axp-1 eAx(p + Ax)]

1 _£-XYp eAx

k

n1

1 - (1 + Y-c) -5. Maximum Likelihood Estimation

In this section, we consider the estimation of unknown parameters of BIIIMW model using maximum lilelihood estimation technique. Let X1, X2,...., xn be a random sample from BIIIMW

*

*

distibution having parameters A=(c, k, a, ft, A) then the log likelihood function is given by,

l(x; c,k,a, ft, A) = ^ ln (f (xi, c,k, a, ft, A)). i=1

The log liklyhood of a single observation x of X is given by,

(8)

l (x; c, k, a, ft, A) = -axft eAx + ln

ckx

-c-1

(1 + x-c)-k-1 + (1 - (1 + x-c)-k) axft-1eAx(ft + Ax) . (9)

-k

The first order derivatives of the log-likelihood function with respect to the parameters A=(c, k, a, ft, A)T is given by,

da

1 - (1 + x- -c )-k (AxfteAx + ftxft-1 eAx)

ckx-c-1(1+ x -c )-k-1 + 1 - (1 + x-c)-k (aAxft eAx + aftx ft-1 eAx )

dl 1 - (1 + x-c)-k (aeAx (xft-1ln(x)ft + xft-1) + aAxfteAxln(x))

dft ckx-c-1 (1 + x-c )-k-1 + 1 - (1 + x-c)-k (aAxft eAx + aftxft-1 eAx)

dl 1 - (1 + x- -c )-k (axft (xAexA + exA) + aftxftexA)

dA ckx-c-1 (1 + x-c )-k- -1 + 1 - (1 + x-c)-k (aAxft eAx + aftx ft-1eAx)

- xft eAx. (10)

- axfteAx ln(x). (11)

- axft+1exA. (12)

dl k (-x-c-1ln(x)(1 + x-c)-k-1c + c(k + 1)x-2c-1ln(x)(1 + x-c)-k-2 + x-c-1 (1 + x-c)-k-1) - A dc =

ckx-c-1 (1 + x-c)-k-1 + 1 - (1 + x-c)-k (aAxfteAx + aftxft-1 eAx)

(13)

where,

A

k (aAxfteAx + aftxft-1 eAx) ln(x)(1 + x-c)-k-1

x .

dl _ cx-c-1 B + (aAxfteAx + aftxft-1 eAx) ln(1 + x-c) (1 + x-c) k dk ckx-c-1 (1 + x-c)-k-1 + [1 - (1 + x-c)-k 1 (aAxfteAx + aftxft-1 eAx)

where B = (1 + x-c)-k-1 [1 - ln (1 + x-c) k].

(14)

The total log-likelihood function of BIIIMW distribution based on a random sample of size n (xi,x2,...,xn) is given by l(A)=E"=i k(A) where k(A)(i=1, 2,..., n) is the log-likelyhood of ith observation. By setting the above partial derivatives to zero, the solution will yield the maximum likelyhood estimators of c, k, a, ft and A. These equations can be solved by using Newton-Raphson method. All the second order derivatives are exist for BIIIMW distribution. Then the observed information matrix is given by,

-E

(Vcc Vck Vca Vcft VcA\

Vkc Vkk Vka Vkft VkA

Vac Vak Vaa Vaft VaA

Vftc Vftk Vfta Vftft VftA

(15)

Wac VAk VAa VAft Vaa/

Here Vjj, j=c, k, a, ft, A denotes the second order derivatives of log-likelihood function with respect to the parameters and E(.) denotes the expexted value. The asymptotic variance and co-variances of these maximum likelihood estimators for c, k, a, ft and A can be obtained by solving this observed information matrix. From (15), the 100(1-£)% confidence intervals for the

parameters are approximately given as follows,

n

1

1

c ± Zs/2\fVcc, k ± Zg/2\fVk, a ± Zg/2\fVoM, P ± Zg/2\jVpp, X ± Zg

where Zg/2 is the upper gth percentile of the standard normal distribution.

6. Simulation and Data Analysis

In this section, we consider the simulation of BIIIMW distribution. Parameters are estimated by using optim CG method in R. We simulate 1000 samples for the true parameters values I : c = 8.5, k = 6, a = 2, p =5 , X = 0.8 and II: c=2, k=9, a=0.5, 0=9, X=1.5. Table 3 lists the means of the MLE estimates, bias and RMSE. From the table, it is obvious that as the sample size grows, estimates approach the true value of the parameter, whereas bias and RMSE decrease as expected.

Table 3: Simulation Results: Mean Estimates, Bias and MSE of BIIIMW distribution.

I II

Sample Size Parameter Mean Bias RMSE Mean Bias RMSE

c 7.8901 -0.6098 0.7106 2.0048 0.0048 0.0070

n=200 k a 7.6779 2.2418 1.6779 0.2418 1.8808 1.5350 9.0031 0.2244 0.0031 -0.2755 0.0036 0.2783

P 7.7369 2.7369 2.9288 9.9856 0.9856 1.0549

A 0.3435 0.2435 0.866 2.3495 0.8495 0.8946

c 7.9066 -0.5933 0.6594 2.0045 0.0045 0.0059

n=500 k 7.6698 1.6698 1.7860 9.0028 0.0028 0.0032

a 1.8931 -0.1068 1.094 0.2281 -0.2734 0.2780

P 7.6866 2.6866 2.786 9.9638 0.9667 0.9752

A 0.3193 0.2193 0.6565 2.3398 0.8398 0.8577

c 8.0117 -0.4882 0.568 2.0042 0.0042 0.0052

n=800 k 7.3521 1.3521 1.5066 9.0025 0.0025 0.0026

a 1.9233 -0.0766 0.3474 0.2286 -0.2712 0.2751

P 7.6809 2.6809 2.7521 9.9581 0.9581 0.9925

A 0.3089 0.2089 0.4682 2.3216 0.8276 0.8396

6.1. Data Analysis

We evaluate a data set in this part to demonstrate the importance and flexibility of the BIIIMW distribution and comparison is done with other well known distributions such as Burr III, modified Weibull, inverse Lomax, exponetiated Weibull and Rayleigh(BIII, MW, IL, EW, R). The data is obtained from Sylwia ([21]) on the lifetime of a certain device (30 device).

Data set : 0.0094, 0.05, 0.4064, 4.6307, 5.1741, 5.8808, 6.3348,7.1645, 7.2316, 8.2604, 9.2662, 9.3812, 9.5223, 9.8783, 9.9346, 10.0192, 10.4077,10.4791, 11.076, 11.325, 11.5284, 11.9226, 12.0294, 12.074, 12.1835,12.3549, 12.5381,12.8049,13.4615, 13.853.

The estimated values of the parameters, Akaike Information Criterion, Bayesian Information Criterion are presented. Also presented Kolmogorov-Smirnov(KS) (its p value), Cramer-Von Mises(W) and Anderson-Darling (A) statistic for hypothesis test. Goodness of fit is performed to test whether the proposed model fits better to the real data sets. In general, the distribution with smallest values of these statistics better fits for the data. Table 4 list the MLE's of model parameters of BIIIMW, BIII, MW, EW, IL, R and the statistics values of W, A and KS.

Table 4: Goodness of fit ofDataset

Model Estim Parameters ates Estimates W A KS P value AIC BIC

BIIIMW c k a ß A 2.693 1920.256 0.0202 0.2089 0.3164 0.076 0.585 0.124 0.692 156.0103 180.0223

BIII c k 0.6941 2.564 1.295 6.373 0.367 0.00038 231.1692 240.774

MW a A ß 3.797 0.021 0.0001 0.118 1.258 0.154 0.428 233.1692 247.5763

EW a ß 0.396 4.658 6.200 1.250 0.348 0.0009 213.734 223.339

IL a ß 0.9301 8.432 1.046 5.277 0.357 0.0006 215.5066 225.1113

R a 6.935 0.456 3.652 0.257 0.030 213.5066 218.3089

From the table 4, it is clear that BIIIMW distribution has the smallest values for these statistics, hence the proposed model is regarded as the better one. The variance- covariance matrix after substituting the unknown parameters of the MLE's in (15), we get,

i 0.00560 -0.012635 -0.000049 0.001199 0.00162 \

-0.01263 0.71116

-0.00004 0.00062

0.00119 -0.08420

V 0.00162 -0.00695

0.00062 -0.08420 -0.00695

0.00001 -0.00428 0.00018

-0.00428 1.56465 -0.14485

0.00018 -0.14485 0.02210

and the corresponding 95% confidence interval is given by c G (2.693 ± 1.96 * 0.0748), k G (1920.256 ± 1.96 * 0.8433), a G (0.0202 ± 1.96 * 0.0031), ft G (0.2089 ± 1.96 * 1.2508) and A G (0.3164 ± 1.96 * 0.1486). Plot of the fitted densities, histogram of the data and pp plot of the real data set is shown in figure (3).

0.0^-0.0

0.4 0.6

Observed

Figure 3: (a) Plot of the estimated pdfs over the histogram and (b) PP plot of the BIIIMW model for dataset.

x

Figure 4: Survival plot of BIIIMW distribution with the estimated parameter values.

Figure 5: (a) The expected parameter values for the BIIIMW distribution hazard rate function (b) The BIIIMW density function for the estimated parameter values.

Figure (4) describe the survival plot of BIIIMW distribution. Figure (5) gives a clear picture that the model satisfy both the necessary and sufficent condition of the distribution of lifetime of system ( both bathtub hazard rate and bimodal density function). The hazard rate function reaches minimum at time, x=0.05 (time period till which the system may undergo breakdown at the beginning of its use itself), followed by normal life till x=4.5 and thereafter exibit an increasing failure rate (wear-out process of the system). Also from figure (5b) it is evidant that the density function attains maximum at time x=0.05 and x=10.5. From the results, it is clear that the model can be used for modelling burn-in procedures.

7. Conclusion

In this paper, we introduce the Burr III-Modified Weibull (BIIIMW) distribution, a new five-parameter lifetime distribution. This distribution can be used to represent both monotone and non-monotone hazard rates in lifetime data. The statistical properties such as quantile fuction, hazard rate function and order statistics are presented. The performance of the new BIIIMW

distribution is performed by using simulation study. Finally, flexibility and applicability of this model is illustrated by a real data set. From the data analysis, it is obvious that the proposed model is highly desirable in the modeling of the lifetime of system. The model satisfy both the necessary and sufficent condition of the distribution of lifetime of system. The proposed method can be used to successfully plan a burn-in process and preventive maintenance of inoperable devices.

References

[1] Bebbington, Mark., Chin-Diew Lai. and Ricardas Zitikis. (2007). A flexible Weibull extension, Reliability Engineering and System Safety, 92, pp. 719-726.

[2] Burr, I. W. (1942). Cumulative frequency functions. The Annals of Mathematical Statistics, 13(2), pp. 215-232.

[3] Burr, I. W. and Cislak, P. J. (1968). On a general system of distributions: I. its curve shaped characteristics, II. the sample median. Journal of the American Statistical Association, 63, pp. 627-635.

[4] Carrasco, J. M. F., Ortega, E. M. M. and Cordeiro, G. M. (2008). A generalized modified Weibull distribution for lifetime modeling. Computational Statistics and Data Analysis. 53(2), pp. 450-462.

[5] Chandrasekaran, R. (1977). Optimal policies for burn-in procedures. Opsearch, 4, pp. 149-160.

[6] Dagum, C. A. (1977). New model of personal income distribution: specification and estimation. Economic Applique, 30(3), pp. 413-437.

[7] Gove, J. H., Ducey, M. J., Leak, W. B., and Zhang, L. (2008). Rotated sigmoid structures in managed uneven-aged northern hardwood stands: a look at the Burr Type III distribution. Forestry, 81(2), pp. 161-176.

[8] M^ltoft, J. (1983). Behind The Bathtub - Curve a New Model and Its Consequences, Microelectronics and Reliability. 23(3),pp 489-500 .

[9] Johnson, N. L., Kot, S. and Balakrishnan, N. (1995). Continuous Univariate Distributions, John Wiley and Sons, New York, USA, 2nd edition.

[10] Kleiber, C. and Kotz, S. (2003). Statistical size distribution in economics and actuarial sciences. New York, John Wiley and Sons.

[11] Lai, C. D., Xie, M., and Murthy, D. N. P. (2003). A modified Weibull distribution. IEEE Trans. Reliab, 52, pp. 33-37.

[12] Lawrence, M. J. (1966). An investigation of the burn-in problem. Technomtrics, 8(1), pp. 61-71.

[13] Meeker, W. Q. and Escobar, L. A. (1998). Statistical Methods for Reliability Data, New York, Wiley.

[14] Mielke, P. W. (1973). Another family of distributions for describing and analyzing precipitation data. J. Appl. Meterol, 12, pp. 275-280.

[15] Khan, M. N. (2015). The modified beta Weibull distribution, Journal of Mathematics and Statistics, 44 , pp. 1553-1568.

[16] Oluyede, B. O., Huang, S. and Yang, T. (2015). A new class of generalized modified weibull distribution with applications. Austrian Journal of Statistics, 44, pp. 45-68.

[17] Park, K. S. (1985). Effect of burn-in on mean residual life. IEEE Transactions on Reliability, 34(5), pp. 522-523.

[18] Pham, H. and Lai, C. D. (2007). On recent generalizations of the weibull distribution. IEEE Transaction on Reliability, 56, pp. 454-458.

[19] Sarhan, A. M. and Zaindin, M. (2009). Modified Weibull distribution. APPS Appl. Sci. 11, pp. 123-136.

[20] Silva, G. O., Ortega, E. M. M. and Cordeiro, G. M. (2010). The Beta ModifiedWeibull Distribution.Lfetime Data Anaysis, 16, pp. 409-430.

[21] Sylwia, K. B. (2007). Makeham's Generalised Distribution. Computational Methods In Science and Technology, 13, pp. 113-120.

[22] Way, K. (1984). Reliability Enhancement Through Optimal Burn-In. IEEE Trans. Reliab, 33(2).

[23] Yanez, S., Escobar, L. A. and Gonzalez N. (2014). Characteristics of two Competing Risks Models with Weibull Distributed Risks. Mathematicas, 38(148), pp. 298-311.

[24] Zwillinger, D. and Jeffrey, A. (2014). Table of Integrals, Series and Products. Academic Press-Elsevier, USA, 7th edition.