Statistical Inference for a Simple Step-Stress Model Based on Censored Data from the Kumaraswamy Weibull Distribution

The step-stress accelerated life tests allow increasing the st ress levels on test units at fixed time during the experiment. In this paper, accelerated life tests are considered when lifetime of a product follows a Kumaraswamy Weibull distribution. The shap e parameter is assumed to be a log linear function of the stress and a cumulative exposure m odel holds. Based on Type II and Type I censoring, the maximum likelihood estimates are obtained for the unknown parameters. The reliability and hazard rate functions are estimated at usual conditions of stress. In addition, confidence intervals of the estimators are constructed. Optimum test plans are obtained to minimize the generalized asymptotic variance of the maxim um likelihood estimators. Monte Carlo simulation is carried out to investigate the prec ision of the maximum likelihood estimates. An application using real data is used to indicate the prope rties of the maximum likelihood estimators.


Introduction
In order to obtain highly reliable products long life-spans, time consuming and expensive tests are often required to collect enough failure data. The standard life-testing methods are not appropriate in such situations and to overcome this difficulty accelerated life tests are applied; wherein the test units are run at higher stress levels (which includes temperature, voltage, pressure, vibration, cycling rate, etc.) to cause rapid failures. Accelerated life tests allow the experimenter to apply severe stresses to obtain information on the parameters of the lifetime distributions more quickly

Original Research Article
than under normal operating conditions. Such tests can reduce the testing time and save a lot of manpower, material sources and money. The stress can be applied in different ways: Commonly used methods are constant stress, progressive stress and step-stress (see [1,2,3]).
In step-stress accelerated life testing (ALT), the stress for survival units is generally changed to a higher stress level at a predetermined time. This model assumes that the remaining life of a unit depends only on the current cumulative fraction failed and current stress. Moreover, if it is held at the current stress, survivors will continue failing according to the cumulative distribution function (cdf) of that stress but starting at the age corresponding to the previous fraction failed. This model is called the cumulative exposure (CE) model. Some references in the field of the accelerated life testing include [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22]. [23] constructed a distribution with two shape parameters on (0, 1). Kumaraswamy (Kum) distribution is applicable to many natural phenomena whose outcomes have lower and upper bounds, such as heights of individuals, scores obtained in a test, atmospheric temperatures and hydrological data. Also, Kum distribution could be appropriate in situations where scientists use probability distributions which have infinite lower and or upper bounds to fit data, when in reality the bounds are finite (see [24]). A compound between Kum distribution and any distribution was constructed by [25].
It has three shape parameters. These parameters allow for a high degree for flexibility of the KumW distribution. Some special cases can be obtained from KumW distribution such as Kum exponential, Kum Rayleigh, exponentiated Weibull, exponentiated Rayleigh, exponentiated exponential, Weibull, Rayleigh and exponential distributions. It is wide applicable in reliability, engineering and in other areas of research. The KumW is discussed in details in [26].
The reliability function (rf) of KumW and the hazard rate function (hrf) corresponding to (2), can be written, respectively, as follows: and hሺt; θ, β, φ, λሻ = The rest of this paper is organized as follows: in Section 2, the k step-stress accelerated life testing is presented. The statistical inference for simple step-stress life testing based on Type II censoring is obtained in Section 3. In Section 4, the statistical inference for simple step-stress life testing based on Type I censoring is discussed.

The K Step-Stress Accelerated Life Testing
Assuming k step-stress accelerated life testing, the model of constant stress is considered in the first step. In this model, the lifetime of the unit is affected by a certain level of stress x 1 , where x 1 is larger than the usual stress x u .
In the consecutive steps, other stresses are considered as x 2 , x 3 , ..., x k , where x u <x 1 <x 2 <..<x k , then the cumulative exposure model reflects the effect of moving from one stress to another one. In the following subsection some basic assumptions are considered.
2-β, λ, φ are constants with respect to the stress x, and the shape parameter θ is affected by the stress x j , j=1, 2, ..., k, through the log linear model in the form where a and b are unknown parameters depending on the nature of the unit and the test method.
3-Suppose that, for a particular pattern of stress, units run at stress x j starting at time τ j-1 and reaching to time τ j , j=1, 2, ..., k, (τ 0 = 0). The behavior of such units is as follows: At Step 1, the population fraction F 1 (t) of units failing by time τ 1 under constant stress x 1 is If F(t) is the population cumulative distribution fraction of units failing under step-stress, then in the first step: where τ 1 is the time when the stress is raised from x 1 to x 2 .

The Cumulative Exposure Model for the Remaining Steps
When Step 2, starts, units have equivalent age u 1 , which have produced the same fraction failed seen at the end of Step 1. In other words the survivors at time τ 1 will be switched to the stress x 2 beginning at the point u 1 , which can be determined as the solution of where ∆ 0 = τ 1 − τ 0 , u 0 = 0 and ∆ j-2 = τ j-1 − τ j-2 , j=2, 3, ..., k, by solving (9), one obtains , by taking the logarithm for two sides, it follows that the cumulative exposure model for j steps can be written as follows: Substituting u j-1 in (11), it is seen that F(t), for a step-stress pattern which consists of segments of the cdf, F 1 , F 2 , ..., F k , can be written in the form: and the associated pdf, f(t), has the following form The maximum likelihood (ML) method is applied to the step-stress model. The pdf for each test is shown by (13), which is the time derivative of the cdf given by (12). The likelihood function is the product of such pdf's evaluated at failure times if complete sampling is used or of observed pdf's of such survival functions evaluated at censoring times when censoring is applied. It is shown that F(t), differs for units with different step-stress patterns. The likelihood function is used to obtain maximum likelihood estimators (MLEs) of the parameters a, b and β.

Inference and Optimal Simple Step-Stress Accelerated Life Tests Based on Type II Censoring
The experiment based on Type II censoring step-stress has the following assumptions: 1-There are k levels of stress x 1 , x 2 , ..., x k , where x 1 <x 2 < ... <x k , are applied such that each unit is initially put under stress x 1 . 2-The experiment begins with n units. The stress x 1 is applied at the first step and the result is n 1 failure times t i1 , i=1, 2, ..., n 1 of test units are observed. When stress x 2 is applied at the second step, n 2 failure times t i2 , i=1, 2, ..., n 2 are observed. Finally, at Step j, stress x j is applied, n j failure times t ij , are observed. 3-The test begins at stress level x 1 if the unit doesn't fail till the predetermined n 1 failures, the stress is raised to x 2 and held until n 2 failures. If it doesn't fail, stress is raised to x 3 . In general, if the unit doesn't fail until the occurrence of n j-1 at stress x j-1 , then the stress is raised to x j at τ j-1 , j=2, 3, ..., k, and held n j failure in case of censored samples, then the test is continued until the occurrence of a predetermined number of failures ∑ n ୨ ୩ ୨ୀଵ . Then there are n c units still survived, at the Step k, the data would be the failure times of (n-n c ) failed units arranged in order and the units which survived beyond t c ( τ nk). 4-Then, it is shown from the previous points that: The stress x j-1 is raised to x j at τ j-1 , j=2, 3, ..., k, when exactly n j-1 failures are observed. It is assumed that the test is continued until all units fail when exactly n k failures are observed. Then, the failure n j-1 , j=2, 3, ..., k is predetermined but τ j-1 , j=2, 3, ..., k and t c are random variables.
The failure time distribution is assumed to be KumW distribution and the shape parameter is shown as a function of the stress through the log linear model. The likelihood function of the experiment is assumed to have the following form: It is shown from (14) that the likelihood function consists of three parts. The first one represents the likelihood of the first step which is the same as the case of constant stress. The second part shows the likelihood function of the (k-1) other stresses. The third part shows the likelihood function of the survived units by time t c . Considering the cumulative exposure model to relate cdf under step-stress to the cdf under constant stress and using the previous assumptions, it is clear that: The failure time distribution at j-th step where As a special case, let k=2, τ ୨ିଵ = τ ଵ and τ ୬ଶ = t ୡ , it is shown that: then hence, the population cumulative fraction of specimens failing in Step 2, by time t is given by: When k=2 there are two steps only with two levels of stress x 1 and x 2 . In this case, the likelihood function has the following form: Suppose λ and φ are known, the logarithm of the likelihood function in (20), denoted by ℓ 2 is given by: where The first derivatives of the logarithm of the likelihood function (21), with respect to a, b and β are obtained.
Therefore, the MLEs can be obtained by equating the first derivatives of ℓ 2 to zero. As shown they are nonlinear equations, the estimates â 2 , ܾ ଶ and ߚ መ ଶ are numerically obtained using Newton Raphson method. Depending on the invariance property of the MLEs, the MLE of the shape parameter, θ u , of the KumW distribution at usual stress x u , can be estimated using the following equation.
also, the MLE of the rf under the same usual conditions ܴ ଶ௨ ሺ‫ݐ‬ ሻ, can be given by and the MLE of the hrf under the same usual conditions ℎ ଶ௨ ሺ‫ݐ‬ ሻ, is given as follows where t 0 is a mission time.
The asymptotic Fisher information matrix can be written as follows: where ߰ ଵ = ܽ, ߰ ଶ = ܾ, ߰ ଷ = ߚ and the elements of the information matrix (25) are derived . percentile. The two sided approximate 100(1-α) % confidence intervals for a, b and β will be respectively, as follows:

The Confidence Intervals Based on Type II Censoring
where ߪ ௪ ෝ is the standard deviation and in this study ‫ݓ‬ ෝ is â, ܾ or ߚ መ , respectively (see [27]).

Optimum Test Plans Based on Type II Censoring
The generalized asymptotic variance (GAV) of the MLE of the model parameters is the reciprocal of the determinant of the asymptotic Fisher information matrix Ĩ 2 (see [11]).
That is Thus, minimization of GAV is equivalent to maximization of the determinant of Ĩ 2 . Newton Raphson method is applied to determine numerically the best choice of the sub sample proportion allocated to each level of stress which minimizes GAV as defined previously. Accordingly, the corresponding optimal numbers of items allocated to each step of stress can be obtained by getting the first partial derivatives of | Ĩ 2 | with respect to G 1 and G 2 .
Then setting (28) equal to zero, where G j are the sub sample proportion which can be optimally determined by solving them simultaneously and applying Newton Raphson method. The determinant can be obtained as follows

Remark
When r=n all the results obtained for Type II censoring, results reduce to those of the complete sample.

Numerical Results
This section aims to investigate the precision of the theoretical results of both estimation and optimal design plans on basis of simulated and real data.

Simulation algorithm
• Several data sets are generated from KumW distribution for a combination of the initial parameter values of a, b and β, and for sample sizes 20, 30, 60 and 100 using 1000 replications for each sample size. • The transformation between uniform distribution and KumW distribution in step j=l is • While the transformation between uniform distribution and KumW distribution in step j=2 is where • The whole sample size n is with initial values of the parameters a=0.5, b=1.5 and β =1.2, given n 1 = 0.4n, n 2 = 0.5n and n c = 0.ln. • It is assumed that there are only two different levels of stress (k=2), x 1 =l and x 2 =1.5, which are higher than the stress at usual condition, x u =0.5. • Number of test units is allocated to each level of stress where G 1 =0.4, G 2 =0.5, r j =90%(n j ), j=l, 2. • The initial parameter values of a, b and β are used in this simulation study to generate t ij , j = 1, 2 and i=l, 2, ..., r j . • Computer program is used depending on MathCad 14 using Newton Raphson method to solve the derived nonlinear logarithmic likelihood equations simultaneously. • Once the values of â 2 , ܾ 2 and ߚ መ 2 are obtained, the estimates are used to obtain, depending on (22) and the design stress, x u =0.5, the shape parameter under this stress, θ u , is estimated as ߠ ଶ௨ = exp ሺâ ଶ + bƹ ଶ x̂ଶ ௨ ሻ. Also, the rf, the hrf and their relative absolute bias are estimated at different values of mission times under usual conditions using (23) and (24). • The performance of the â 2 , ܾ 2 and ߚ መ 2 has been evaluated through some measurements of accuracy. In order to study the precision and variation of MLEs (E 2 ), it is convenient to use, the relative absolute bias (RAB 2 ), the mean square error (ER 2 ) and the relative error (RE 2 ). • The two sided approximate 100(1-α) % confidence intervals for a, b and β will be obtained using (26). The different sample sizes of n=20, 30, 60, 100 are considered.

•
The results are displayed in Tables 1-4.

Concluding remarks
• It is clear from Table l that the MLEs (E 2 ) are very close to the initial values of the parameters as the sample size increases. Also, as shown in the numerical results the RAB 2 , ER 2 and RE 2 are decreasing when the sample size is increasing. For all sample sizes we noted that: • ߚ መ 2 performs better than other estimates. • ܾ 2 performs better than â 2 .
• Table 2 indicates that the reliability decreases when the mission time t 0 increases. The results get better in the sense that the aim of an accelerated life testing experiments is to get large number of failures (reduce the reliability) of the device with high reliability. As t 0 increases the RAB R2 increases and when the sample size increases, the rf increases. Also, the RAB R2 for the rf decreases when the sample size increases. The hrf increases when the mission time t 0 increases and when t 0 increases the RAB h2 decreases. • The two-sided 95% central asymptotic confidence intervals for the parameters of KumW are displayed in Table 3. This table contains the standard error (SE 2 ), lower bound (L 2 ), upper bound (U 2 ) and the length of the intervals. The interval estimate of the parameters becomes narrower as the sample size increases. For all sample sizes, it is clear that: • The length of the interval for β is shorter than the other lengths.
• The length of the interval for b is shorter than the length of the interval for a.
• Optimum test plans are developed numerically; it can be observed from the numerical results presented in Table 4, that the optimum test plans do not allocate the same number of the test units to each stress. Also, Table 4, includes the expected number of items that must be allocated to each level of stress represented by ݊ ଵ * , ݊ ଶ * which minimizes the GAV. As indicated from the results, the optimal GAV of the MLE of the model parameters decreased as the sample size n increased.

Application
The main aim of this subsection is to demonstrate how the proposed method can be used in practice. [25] used Kolmogorov-Smirnov goodness of fit test and data points representing failure time. The data were taken from [28]. The data were 30 items (n=30) tested with test stopped after 20 th failure (r=20). It is assumed that k=2, i.e. there are only two different levels of stresses x 1 =0.6 and x 2 =1, which are higher than the stress at usual conditions, x u =0.5. The failure times in the first step are [0.0014, 0.0623, 1.3826, 2.0130, 2.5274, 2.8221, 3.1544, 4.9835, 5.5462, 5.8196, 5.8714, 7.4710] and the failure times in the second step are [7.5080, 7.6667, 8.6122, 9.0442, 9.1153, 9.6477, 10.1547, 10.7582].
The initial parameter values of a, b and β used in this application are a=0.5, b =1.5, β =2, λ=2 and φ=2. Once the estimate values of a, b and β are obtained, the estimators are used to estimate θ u , as ߠ ଶ௨ = ݁‫ݔ‬൫â ଶ + ܾ ଶ ‫ݔ‬ ଶ௨ ൯. Letting the design stress, x u = 0.5. Also, the reliability function is estimated at different values of mission times under usual conditions depending on (22). Moreover, the precision and variation of MLEs (E 2 ) are studied through some convenient measures such as the RAB 2 , ER 2 and RE 2 . These measures are computed for each parameter in Table 5.
The estimated shape parameter, rf, RAB R2 , the hrf and RAB h2 under usual condition are shown in Table 6. Table 7 and Table 8 indicate confidence intervals of the parameters at confidence levels 95% and 99%. These tables contain the standard error (SE 2 ), lower bound (L 2 ), upper bound (U 2 ) and the length of the intervals.
The relationship between the stress and the shape parameter is tested through testing the significance of the coefficient b. Hypothesis test is obtained when α=0.05 and with one degree of freedom, assuming the null hypothesis is b=0. It is rejected and the relationship between the level of the stress and the shape parameter exist.

Inference and Optimal Simple Step-Stress Accelerated Life Tests Based on Type I Censoring
In Type I censoring step-stress, the stress x j-1 is raised to x j at τ j-1 , j=2, 3, ..., k. It is assumed that the test is continued until all units fail or until time t c . The difference between time step-stress and failure step-stress, is that in failure step-stress, n j-1 , j=2, 3, ..., k+1, are predetermined but τ j-1 , j=2, 3, ..., k and t c are random variables. On the other hand, τ j-1 , and t c are predetermined in time stepstress and n j-1 , j=2, 3, ..., k+1, are random variables.

The Maximum Likelihood Estimation Based on Type I Censoring when there are 2 Steps of Stress as a Special Case
As a special case, let k=2, τ 1 is the time at which the stress changes from x 1 to x 2 and t c is the time at which the experiment is terminated (censoring time). The likelihood function of the experiment is considered to have the same form as (21) but τ j-1 and t c are predetermined in time step-stress and n j-1 , j=2, 3, ..., k+1 are random variables. Then the maximum likelihood estimates are obtained for the unknown parameters. The reliability and the hazard rate functions are estimated.
In addition, confidence intervals of the estimators are constructed. Optimum test plans are obtained to minimize the generalized asymptotic variance of the maximum likelihood estimators.

Numerical Results
This subsection aims to illustrate the precision of the theoretical results of both estimation and optimal design problems on basis of simulated data.

Simulation algorithm
The same steps of the algorithm in Subsection (3.3) will be considered in this algorithm with the following data: • The values τ 0 = 0, τ 1 =2 and t c = 5.5 are given. Once the values of â 1 , ܾ 1 and ߚ መ 1 are obtained, the estimates are used to obtain, depending on (22) and the design stress, x u =0.5, the shape parameter under this stress, θ u , is estimated as ߠ ଵ௨ =exp (ܽ ො 1 +ܾ 1 x 1u ). Also, the reliability function, the hazard rate function and their relative absolute bias are estimated at different values of mission times under usual conditions using (23) and (24). • The performance of the estimates â 1 , ܾ 1 and ߚ መ 1 has been evaluated through some measurements of accuracy. In order to study the precision and variation of MLEs, it is convenient to use the relative absolute bias (RAB 1 ), the mean square error (ER 1 ) and the relative error (RE 1 ). Depending on the same procedure as in Section 3, the numerical results of the experiment are displayed in Tables 9-12.

Concluding remarks
• It is clear from Table 9 that the MLEs (E 1 ) are very close to the initial values of the parameters as the sample size increases. Also, as shown in the numerical results the RAB 1 , ER 1 and RE 1 are decreasing when the sample size is increasing. For all sample sizes we noted that: • ߚ መ 1 performs better than other estimates. • ܾ 1 performs better than â 1 .
• Table 10, indicates that the reliability decreases when the mission time t 0 increases. The results get better in the sense that the aim of an accelerated life testing experiments is to get large number of failures (reduce the reliability) of the device with high reliability. As t 0 increases the RAB R1 increases and when sample size increases, the rf increases. Also, the RAB R1 for the rf decreases when the sample size increases. The hrf increases when the mission time t 0 increases and when t 0 increases the RAB h1 decreases. • The two-sided 95% central asymptotic confidence intervals for the parameters of KumW are displayed in Table 11. This table contains the standard error (SE 1 ), lower bound (L 1 ), upper bound (U 1 ) and the length of the intervals. The interval estimate of the parameters becomes narrower as the sample size increases.
• As shown in Section 3, by setting the డ|ூ ሚ| డఛ భ = 0, ߬ ଵ and t c can be optimally determined by solving them simultaneously. For all sample sizes, it is clear that: • The length of the interval for β is shorter than the other lengths.
• The length of the interval for b is shorter than the length of the interval for a.
• optimum test plans are developed numerically. The expected time, ߬ ଵ * , at which the stress changes from x 1 to x 2 and the expected time, ‫ݐ‬ * , at which the experiment is terminated are displayed in Table 12. As indicated from the results, the optimal GAV of the MLE of the model parameters is decreasing as the sample size n is increasing.

Remarks
The results obtained in this paper can be modified to obtain results for sub-models of KumW distribution under Type I and Type II censored samples such as • The Kum exponential distribution if φ = 1.