E-Bayesian and H-Bayesian Inferences for a Simple Step-Stress Model with Competing Failure Model under Progressively Type-II Censoring

In this paper, we discuss the statistical analysis of a simple step-stress accelerated competing failure model under progressively Type-II censoring. It is assumed that there is more than one cause of failure, and the lifetime of the experimental units at each stress level follows exponential distribution. The distribution functions under different stress levels are connected through the cumulative exposure model. The maximum likelihood, Bayesian, Expected Bayesian, and Hierarchical Bayesian estimations of the model parameters are derived based on the different loss function. Based on Monte Carlo Simulations. We also get the average length and the coverage probability of the 95% confidence intervals and highest posterior density credible intervals of the parameters. From the numerical studies, it can be seen that the proposed Expected Bayesian estimations and Hierarchical Bayesian estimations have better performance in terms of the average estimates and mean squared errors, respectively. Finally, the methods of statistical inference discussed here are illustrated with a numerical example.


Introduction
In recent years, it has been challenging to obtain sufficient failure data during the general service condition with increasing reliability of products. Accelerated life testing (ALT) is conducted to overcome such difficulties. In ALT, the test units are subjected to higher levels of stress condition for shortening the total testing time, and sufficient failure data can be obtained for reliability assessment. Some authentic books in the field of the ALT include Nelson et al. [1], Nelson [2], Meeker et al. [3], and Bagdonavicius et al. [4]. The constant-stress accelerated life testing (C-SALT) and the step-stress accelerated life testing (S-SALT) are two special types of ALT. S-SALT an advantage of yielding more failure data in a limited testing time and changes the stress level at a prefixed time or a prefixed number of failures during the testing. To analyze the data from S-SALT, one requires a model that relates to the failure lifetimes under different stress levels. The cumulative exposure model (CEM) is the most studied in the literature, and the model was first introduced by Sedyakin [5]. The S-SALT under the assumption of CEM has attracted great attention, such as Balakrishnan et al. [6,7], Lee et al. [8], and Tang [9]. Furthermore, it is common that the products' failure may be caused by more than one cause. Therefore, each of these causes would compete with each other to result in the final failure. We call those possible failure causes competing failures. The competing failure model plays an important role and the competing failure data are analyzed by Cox [10] and David et al. [11]. In S-SALT, some researchers have discussed the competing failure model, such as Balakrishnan et al. [12], Beltrami [13,14], Shi and Liu [15], Srivastava et al. [16], Xu et al. [17], Zhang et al. [18,19], Ganguly et al. [20], Han et al. [21][22][23], Varhgese et al. [24], Liu et al. [25], Abu-Zinadah et al. [26] and Aljohaniet al. [27]. The maximum likelihood estimation (MLE) and the Bayesian estimation (BE) based on the different loss function (LF) are the common inference for analyzing statistical data. Lindley [28] introduced the Hierarchical Bayesian estimation (H-BE) primarily, and it was examined by Han [29]. The Expected Bayesian estimation (E-BE) is the expectation of BE, and it was introduced by Han [30]. Han [31,32] derived the E-BE and H-BE of the reliability parameter under testing data. However, there are few works concerning the E-Bayesian and H-Bayesian inference for the step-stress accelerated competing failure model.
We will discuss the EB and HB inference for a simple step-stress accelerated competing failure model which has only two stress levels, and the layout of the paper is organized as follows. Section 2 gives the model and basic assumptions. Asymptotic confidence intervals (CIs) and Bootstrap confidence intervals (BCIs) are constructed in Section 3. Bayesian estimations (BEs) are derived based on the different LFs in Section 4. The E-BE is derived based on the different LFs in Section 5. The H-BEs is derived based on the different LFs in Section 6. The Simulation data analysis is provided in Section 7. Section 8 illustrates the proposed methods by a numerical example. Finally, we provide concluding remarks and future research proposes.

Model Description and MLEs
In this section, we discuss a simple-stress life testing and provide the MLEs of the unknown parameters based the observed data.

Basic Assumption
To describe the simple S-SALT clearly, some assumptions are made as follows: (1) The unit fails only due to one of the two independent competing failure causes with failure times T 1 and T 2 , respectively. A failure time is recorded as joint random variable (T, δ). Let δ = 1, the failure is caused by the first cause 2, the failure is caused by the sec ond cause , that denote the indicator variable for the cause of failure time T = min(T 1 , T 2 ). (2) The lifetime follows an exponential distribution with scale parameter λ ij . Let 1/λ ij be the mean time-to-failure of a test unit at the stress level S i by the failure cause j for i, j = 1, 2. The cumulative distribution function (CDF) and probability density function (PDF) are given as follows, respectively where t ≥ 0, λ ij > 0 and i, j = 1, 2. (3) The scale parameter λ ij agrees with a log-linear function of stress where a j and b j are unknown coefficient parameters, ϕ(S i ) = 1/S i that is chosen as the Arrhenius model [2] in this paper, i, j = 1, 2, but this not used explicitly in this paper. (4) Lifetime distribution at different stress levels is related by CEM. The CEM assumes that the remaining lifetime of a unit depends only on the cumulative exposure accumulated at the current stress level, regardless of how the exposure is actually accumulated. The failure probability of product working time t 1 under stress S 1 is equivalent to the failure probability of product working time t 2 under stress S 2 . At the time τ when the stress level increases from S 1 to S 2 , the CDF of the lifetime of the test unit failed due to cause j for j = 1, 2 can be written as follows: where F 1j and F 2j are the CDFs of the lifetime of the test unit failed under stress S 1 and S 2 , respectively, and a j is such that it satisfies F 1j (τ) = F 2j (τ − a j ) and Under these assumptions, the CDF of the lifetime of the test unit failed due to cause j is and the corresponding PDF is given by where j = 1, 2.

Maximum Likelihood Estimates
Based on the assumptions (1)-(4) and the lifetime T = min(T 1 , T 2 ) of the test unit, the CDF and PDF of T can be obtained as follows: Then the joint distribution of (T, δ) is given by With the life-testing scheme described above, the following ordered failure time will be observed: t 1:n < t 2:n < · · · < t N 1 :n < t N 1 +1:n < · · · < t N 1 +N 2 :n n 1j is the number of units that fail under stress S 1 due to the failure cause j; n 2j is the number of units that fail under stress S 2 due to the failure cause j; for j = 1, 2 such that N 1 = n 11 + n 12 and N 2 = n 21 + n 22 . The observed failure times t = (t 1:n , t 2:n , · · · , t N 1 :n , t N 1 +1:n , · · · , t N 1 +N 2 :n ) and the Θ = (λ 11 , λ 12 , λ 21 , λ 22 ). Then, the likelihood function can be written as where By using Equation (10), the log-likelihood function can be written as: Taking the first partial derivative of Equation (11) with respect to λ ij , and let it equal zero Using simple algebra calculations, we obtain

Interval Estimations
In this section, we propose the asymptotic confidence intervals (ACIs) and the bootstrap confidence intervals (BCIs) for the parameters λ ij .

Asymptotic Confidence Intervals (ACIs)
According to the asymptotic likelihood theory, we get the information matrix of Θ, the elements of which are Then, the observed Fisher information matrix of Θ is that The approximate asymptotic variance-covariance matrix can be given byÎ(Θ) where Z γ/2 is the γ/2th upper percentile of the standard normal distribution.

Bootstrap-t Method
Hall (1988) developed Bootstrap-t method, and the Bootstrap-t samples are generated as follows: Step 1. Same as Bootstrap-p method step 1.
Step 2. Same as Bootstrap-p method step 2.

Bayesian Analysis
In this section, we consider the BEs of λ ij based on the squared error loss function (SELF), entropy loss function (ELF), and LINEX (linear-exponential) loss function (LLF). The loss functions expressions are in the Appendix A.
As the conjugate prior, an independent gamma was chosen as prior distribution Based on Equations (10) and (19), the joint posterior density function of λ ij can be written as: where ξ= (n ij , t i , R i , τ, i = 1, · · · , N; j = 1, 2 .

Bayesian Estimation of λ ij under SELF
The BE of λ ij under the SELF is the expectation of the posterior distribution. Therefore, λ ij(BS) can be derived aŝ Thus, the BEs of λ 1j and λ 2j are obtained under Equations (21) and (22) aŝ

Expected Bayesian Analysis
Since the values of the hyper-parameters are not easy to determine, it has certain randomness. The Expected Bayesian estimation is the average of Bayes estimates of θ with hyper-parameters a and b in the domain Θ. In this section, we will obtain the E-BE of based on the SELF, ELF, and LLF. The definition of E-BE was originally addressed by Han [29] as follows.
is called the E-BE of θ, which is assumed to be finite, where Θ is the domain of a and b,θ B (a, b) is the BE of θ with hyper-parameters a and b, and π(a, b) is the joint density function of a and b over Θ. By Definition 1, the expected Bayesian is the expectation ofθ B (a, b) with hyper-parameters a and b.
According to the prior information, the λ ij is large with a low probability, and the λ ij is small with a high probability, so the hyper-parameters α ij and β ij should be chosen to guarantee that π(λ ij α ij , β ij ) is a decreasing function of λ ij . For more details, see Han [29]. The derivative of π(λ ij α ij , β ij ) with respect to λ ij is, Given 0 < α ij < 1, then the larger the value of β ij , the thinner the tail of the density function. Berger [33] showed that the hyper parameter β ij should be chosen under the restriction 0 < β ij < c (c is a constant). Suppose that α ij and β ij are independent, the joint density of α ij and β ij is given by π α ij , β ij = π α ij π β ij . Depending on different distribution of the parameters α ij and β ij , E-BE estimation of λ ij is obtained. Several authors have applied the E-BE method to analyze data, such as Abdul-Sathar et al. [34] and Shahram [35].
E-BE of λ ij is obtained relying on different distributions of α ij and β ij . These distributions are used to describe the effect of the different prior distributions on E-BE of λ ij . In this paper, the π α ij , β ij may be used:

E-Bayesian Estimation of λ ij under SELF
The E-BEs of λ ij are obtained under Equations (23) and (30a) aŝ The E-BEs of λ ij are obtained under Equations (23) and (30b) aŝ

E-Bayesian Estimation of λ ij under ELF
The E-BEs of λ ij are obtained under Equations (25) and (30a) aŝ The E-BEs of λ ij are obtained under Equations (25) and (30b) aŝ

E-Bayesian Estimation of λ ij under LLF
The E-BEs of λ ij are obtained under Equations (27) and (30a) aŝ The E-BEs of λ ij are obtained under Equations (27) and (30b) aŝ

Hierarchical Bayesian Estimation
In this section, we obtain the H-BEs of λ ij based on the SELF, ELF, and LLF. The definition of H-BE was originally addressed by Lindley and Smith [28] as follows.
Definition 2. If a and b are hyper-parameters in the parameter of θ , the prior density function of θ is π(θ|a, b) , and the prior density function of the hyper-parameters of a and b is π(a, b) , then the hierarchical prior (H-P) density function of θ is defined as follows: where Θ is the domain of a and b.
The H-P density functions of the parameters λ ij are obtained under Equations (19) and (30a) as: We get the posterior density function of λ ij under Equations (20) and (38) as: In the same way, The H-P density function of the parameters λ ij are obtained under Equations (19) and (30b) as: We get the posterior density function of λ ij under Equations (20) and (40) as:

H-Bayesian Estimation of λ ij under SELF
The H-BEs of λ ij are obtained under Equations (22) and (39) as: The H-BEs of λ ij are obtained under Equations (22) and (41) as:

H-Bayesian Estimation of λ ij under LLF
The H-BEs of λ ij are obtained under Equations (26) and (39) as: The H-BEs of λ ij are obtained under Equations (26) and (41) as:

Highest Posterior Density (HPD) Credible Intervals (CRIs)
Here, we propose the following algorithm to compute the associated CRI and HPD CRI (reference [36]).
Step 2. Arrange the sample in ascending order to obtain the CRIs of the parameters ij , i, j = 1, 2.

Simulation Study
In order to investigate the proposed methods, we use Monte Carlo simulations to compare different methods for different sample size and different progressive censoring schemes, which are shown in Table 1. The values of the parameters are chosen to be λ 11 = 2.0, λ 12 = 1.0, λ 21 = 4.0, λ 22 = 2.0, N * = 1000, c = 1/5 and k = 4. For different censoring schemes, we compute the Average Estimates (AEs) and the mean square errors (MSEs) of the MLEs, BEs, E-BEs, and H-BEs, respectively. AEs and MSEs of the estimator of λ ij can be calculated as ij is the kth estimator of the parameter λ ij (i, j = 1, 2). Simulation study has been done according to the following steps: Step 1. For given λ 11 = 2.0, λ 12 = 1.0, λ 21 = 4.0, generate a PT-IIC samples based on different sample sizes.
Step 5. Repeat step 1 to step 4 N * times and calculate the AEs and MSEs for each estimate. The results are presented in Tables 2-13.
Step 6. Compute the average lengths (Als) and the coverage probabilities (CPs) of the 95% asymptotic Cis of the MLEs, the Bootstrap-p method, the Bootstrap-t method, the Bayesian method, the EB method, and the HB method. The results are presented in Tables 14-17.
Step 7. Compute the Als and the CPs of the HPD CRIs using the Bayesian method, the EB method, and the HB method, and the results are presented in Tables 18-21.

Results Analysis
From Tables 2-21, some conclusions are summarized as follows: (1) From Tables 2-13, we observe that the AEs of λ ij (i, j = 1, 2) are close to the true values and the MSEs of λ ij (i, j = 1, 2) decrease as n increasing for all estimates. This indicates that the number of failure values of test units affect the estimation accuracy of parameters. (2) From Tables 2-13, the Bayesian performances are better than that of MLE, and the E-Bayesian or H-Bayesian performances are better than that of the Bayesian for fixed n, N 1 , N 2 , and censoring scheme. The results show that the Bayesian method improves the estimation accuracy of model parameters due to combining the prior information.

An Illustrative Example
In this section, we simulate a PT-IIC sample from a simple step-stress accelerated competing failure model. The dataset is generated with the following choices of the parameters: λ 11 = 1.0, λ 12 = 1.5, λ 21 = 2.0 and λ 22 = 3.0, and n= 30, N 1 = 10, R 1 = 4, N 2 = 12, and R 2 = 4. The data are given in Table 22. From this dataset, we have n 11 = 4, n 12 = 6, n 21 = 5 and n 22 = 7, and the Average Estimates (AEs) of MLEs, BEs, EBs, and HBs of the parameters are derived based on the SELF. The results are presented in Table 23. From Table 23, it is clearly observed that the E-Bayesian or H-Bayesian performances are better than the MLEs.
We constructed the ALs of the 95% asymptotic CIs and the results are presented in Table 24. We also consider the HPDCRIs, and the results are presented in Table 25. From  Tables 24 and 25, we observe that the H-Bayesian CIs are always narrower than the other CIs. The HPDCRIs are always narrower than the other CIs.

Conclusions
This paper has proposed a simple S-SALT based on the CEM and under different PT-IIC schemes. It has been assumed that the lifetimes have an exponential distribution with different scale parameters. We have derived the maximum likelihood estimations, Bayesian estimations, Expected Bayesian estimations, and Hierarchical Bayesian estimations of the scale parameters based on the different LF. Based on Monte Carlo Simulations, we also obtained the ALs and the CPs of the 95% ACIs, BPCIs, BTCIs, and HPDCRIs for all the unknown parameters. Results show that the MSEs of the scale parameters decrease as the experimental units N increasing for all estimators, the H-BEs are the best in all cases of the larger samples sizes, and the E-BEs are the best in all cases of the smaller samples sizes. From the perspective of the CIs and CRIs, it has been observed that the H-Bayesian CIs are always narrower than the other Cis, and the HPDCRIs are always narrower than the other CIs under the same loss function. Note that in this paper we assume that the two competing failures are independent. In future work, the simple step-stress accelerated dependent competing failure model will be considered, such as Copula [37,38], which is one of the popular models for releasing the restriction.  where k determines the shape of the loss function and k = 0.