Optimal Constant-Stress Accelerated Life Test Plans for One-Shot Devices with Components Having Exponential Lifetimes under Gamma Frailty Models

: Optimal designs of constant-stress accelerated life test plans is one of the important topics in reliability studies. Many devices produced have very high reliability under normal operating conditions. The question then arises of how to make the optimal decisions on life test plans to collect sufﬁcient information about the corresponding lifetime distributions. Accelerated life testing has be-come a popular approach to tackling this problem in reliability studies, which attempts to extrapolate from the information obtained from accelerated testing conditions to normal operating conditions. In this paper, we develop a general framework to obtain optimal constant-stress accelerated life test plans for one-shot devices with dependent components, subject to time and budget constraints. The optimal accelerated test plan considers an economical approach to determine the inspection time and the sample size of each accelerating testing condition so that the asymptotic variance of the maximum likelihood estimator for the mean lifetime under normal operating conditions is minimized. This study also investigates the impact of the dependence between components on the optimal designs and provides practical recommendations on constant-stress accelerated life test plans for one-shot devices with dependent components.


Introduction
One-shot device test data analysis has received increased attention in reliability studies. One-shot devices, such as vehicle airbags, electro-explosive devices, and missiles, are immediately destroyed after testing. The single-use devices result in left-censoring and right-censoring. The exact lifetime of such a device cannot be observed, no matter whether the life test is conducted successfully or not. Furthermore, one-shot device test data can be obtained in destructive tests, wherein each specimen/item is inspected only once, for instance, Li/SO 2 storage batteries [1], rolling ball bearing [2], and grease-based magnetorheological fluids [3]. In addition, many devices/products have lifetimes of years or decades under normal operating conditions. When practitioners are dealing with one-shot device test data, they are usually facing great challenges with a limited time and resources, in addition to a lack of information inherently presented in the data collection. In this regard, accelerated life testing has been widely adopted in reliability studies to collect sufficient lifetime information by subjecting specimens/devices to higher-than-normal operating conditions to induce rapid failures. A recent book [4] provides an overview of accelerated life testing of one-shot devices and presents several inferential methods for analyzing one-shot devices. Balakrishnan and his collaborators have made substantial efforts in the past few years to further advance the methodologies for analysis and data collection of one-shot devices [5][6][7][8][9][10][11][12]. time and budget constraints, based on V-optimality, for one-shot devices with dependent components having exponential lifetime distributions under frailty models.
The rest of this paper is organized as follows. Section 2 presents one-shot device test data with multiple failure modes collected from a CSALT and introduces the gamma frailty model with exponential component lifetime distributions as well as the mean lifetime under normal operating conditions. Section 3 presents the Fisher information matrix and the asymptotic variance of the maximum likelihood estimator for the mean lifetime. These are required to develop a general framework to obtain the optimal CSALT plan with time and budget constraints for one-shot devices with dependent components. In Section 4, a simulation study is carried out to evaluate the performance of the developed framework as well as the impact on the optimal design of CSALT under various degrees of dependence. Based on the simulation results, some practical recommendations on CSALT plans will be provided. A real example from a diabetic retinopathy study is presented to demonstrate the proposed framework in Section 5. Finally, some concluding remarks are given in Section 6.

Model Description
In this section, we describe one-shot device test data with multiple failure modes collected from CSALT. Interested readers may refer to [33] in connection to the likelihood inference for one-shot devices in detail. Suppose that N one-shot devices consisting of M components are partitioned into I higher-than-normal operating conditions. For i = 1, 2, . . . , I, N i devices are placed at a higher-than-normal stress level s i and inspected at K i equally-spaced time points. There are ∑ I i=1 K i test groups in total. In particular, N i,k devices are selected and inspected at a specific time τ i,k , for k = 1, 2, . . . , K i . As each one-shot device can perform its intended function only once at a specified inspection time, a practitioner can observe either a success or a failure. A successful test indicates that all the components within the device have not malfunctioned at the inspection time. On the other hand, if the test fails, it is assumed that the practitioner can determine which components of the failed device have malfunctioned after a careful investigation or autopsy of the device. Here, let P (Ω) be the power set of Ω, where Ω = {1, 2, . . . , M}. The power set with 2 M elements represent the set of all possible combinations of malfunctioned components. Let n i,k,X denote the number of devices with malfunctioned components X ∈ P (Ω) at inspection time τ i,k . Finally, we use z = {s i , τ i,k , n i,k,X , i = 1, 2, . . . , I, k = 1, 2, . . . , K i , X ∈ P (Ω)} to denote the observed data, with ∑ X∈P (Ω) n i,k,X = N i,k .
In the i-th operating condition, where devices are exposed to an elevated stress level s i , let T m denote the failure time of the m-th component, for m = 1, 2, . . . , M. For the j-th device, it is assumed that, conditioned on a latent (unobserved) frailty, γ i,j > 0, T m follows an exponential distribution with rate parameter γ i,j λ i,m > 0. The conditional probability density functions (pdf) are then given by The corresponding conditional reliability function is It is further assumed that the failure rate of the m-th component is related to the elevated operating condition s i by using a log-linear function [34], namely, λ i,m = exp(a m0 + a m1 s i ).
Furthermore, the frailty γ i,j is assumed to follow a gamma distribution with scale parameter β > 0 and shape parameter 1/β > 0. The pdf of γ is then Note that the mean of the frailty is 1 and the variance is β so that the model is identifiable [17]. Moreover, Ling et al. [33] found that some constraints on β are imposed to ensure that the mean (β < 1) and the variance (β < 0.5) of component lifetimes are finite.
Furthermore, conditioned on γ i,j , the failure times of those M components (T 1 , T 2 , . . . , T M ) are assumed to be independent, namely, By integrating R T Ω (t 1 , t 2 , . . . , t M |y) f γ (y) with respect to y over (0, ∞), the (unconditional) joint reliability function is then given by ( In addition, from (2), we readily find the (unconditional) reliability function for the m-th component, namely, It is worth noting that T m has a Lomax (or Pareto Type II) distribution with scale parameter (βλ i,m ) −1 and shape parameter 1/β (see [35,36]). Its Kendall's tau is β/(β + 2), which is commonly used to measure the correlation between variables. It is therefore evident that an increase in β results in a stronger dependence between these M components. In addition, these components are independent when β tends to zero.
If the device can perform its intended function when all the M components work, from the results in [33], the mean lifetime of such a device at the normal operating condition s 0 is given by where β < 1 and λ 0,m = exp(a m0 + a m1 s 0 ). From (4), it is worth noting that an increase in the dependence between the components in the device, i.e., β, results in an increase in the mean lifetime of the device.

Optimal CSALT Plans
In this section, we present a general framework to determine the optimal ALT plan with time and budget constraints, which minimizes the asymptotic variance of the maximum likelihood estimator of µ D in (4) by selecting the inspection time (τ i,K i ) and the sample size (N i,k ) of each accelerating testing condition. The required optimization problem needs the Fisher information matrix as well as the first-order derivatives of the mean lifetime with respect to the model parameters.

Asymptotic Variance
For a one-shot device with M components, let θ = (a 10 , a 11 , . . . , a M0 , a M1 , β) denote the vector of 2M + 1 model parameters. The log-likelihood function, based on the observed data z = {s i , τ i,k , n i,k,X , i = 1, 2, . . . , I, k = 1, 2, . . . , K i , X ∈ P (Ω)}, is given by where X c is the complement of X, and c is a constant. We can further obtain the Fisher information matrix, which is the negative of the expectation of the second-order derivatives of the log-likelihood function in (5) with respect to the model parameters, namely, where with I(m ∈ {Y, X c }) denoting an indicator function that takes the value 1 if m ∈ {Y, X c } and takes the value 0 if m / ∈ {Y, X c }. Interested readers may refer to [33] for the derivation of the first-order derivatives.
The asymptotic variance-covariance matrix of the maximum likelihood estimators for the model parameters is the inverse of the Fisher information matrix, that is, V(θ) = (I(θ)) −1 . We further apply the delta method in [33] to obtain the asymptotic variance of the maximum likelihood estimator for µ D , which requires the first-order derivatives of µ D with respect to the model parameters θ = (a 10 , a 11 , . . . , a M0 , a M1 , β), namely, Consequently, the asymptotic variance by using the delta method is given by However, the MLE of µ D is not normally distributed unless the sample size is sufficiently large. Ling et al. [33] revealed that the logarithm of the MLE has an approximately normal distribution, and thus the asymptotic variance of the MLE after log-transformation is considered, that is It is obvious that the optimal CSALT minimizes both V ln(µ D ) and V µ D simultaneously.
On the other hand, recycling and waste management are important environmental issues worldwide. For example, iron, copper, zinc, and rare metals, such as platinum and palladium, are used in the manufacture of automobiles. These valuable metals can be collected from deregistered automobiles and traded as a valuable secondary resource. In addition, components such as batteries, engines, tires, and bumpers are recyclables and valuable materials for secondary use and, thus, very often re-sold or recycled as secondary products [37]. In this regard, it is of great interest to consider the recycling of components in failed devices and estimate the revenue from the sale of recycled components.
In testing of one-shot devices with multiple components, suppose that a one-shot device has a set of X ∈ Ω malfunctioned components. The rest of the components, X c , have not malfunctioned and can be re-sold as valuable secondary resources. We let R(X c ) denote the revenue from selling the recycled components in the device. The total revenue for the recycling is then given by For i = 1, 2, . . . , I and k = 1, 2, . . . , K i , n i,k,X has a multinominal distribution with the sample size N i,k and the corresponding probabilities P i = (P i (∅, τ i ), · · · , P i (Ω, τ i )) . Therefore, the expected revenue is and the corresponding variance is obtained as where P R = (R(Ω), · · · , R(∅)) is a vector of the corresponding revenues.

Procedure of Obtaining the Optimal CSALT Plan
Here, we assume that the cost of conducting a CSALT consists of two terms: the cost linked to the purchase of the devices and the cost associated with operation at each accelerating operating condition. In this study, it is also assumed that the operation cost increases as a result of an increase in the elevated operating condition s i , as well as the duration of the CSALT, because more resources will normally be consumed to increase the elevated operating condition and conduct the CSALT for a longer duration. Given a CSALT plan ξ = {K i , τ i,k , N i,k , i = 1, 2, . . . , I, k = 0, 1 . . . , K i }, the total cost of conducting the CSALT is, therefore, given by where τ i,K i and C i,oper are the duration and the operation cost, respectively, for the i-th operating condition, and C item is the cost of each device. Moreover, time and budget constraints are usually imposed to CSALT in practice, and therefore, the total cost TC(ξ) cannot exceed a specified budget C budget and the duration τ i,K i cannot exceed a specified termination time τ ter . In short, the optimal CSALT plan can be determined by solving the following objective function subject to time and budget constraints.  Here, · is a truncated integer, f represents the inspection frequency, with τ i,k = f k, and K * i represents the maximum number of inspections in the i-th operating condition. It is worth noting that when the number of test groups is large, the optimization function "optim()" provides reasonable but not optimal sample sizes due to a high dimensional objective function involved in the optimization problem. Therefore, some inspection times τ i,k are discarded if the sample sizes are very small, e.g., N i,k ≤ N min , so that more resources will be concentrated on the remaining test groups and the results obtained by the optimization function will be more stable due to a lower-dimensional objective function. At the same time, the duration for some operating conditions may be shortened due to the discard of inspection times, resulting in more devices being purchased for testing. It is worth noting that the median or quantile can also be used to replace the mean lifetime in the objective function. One should follow the same process in Section 3.1 to determine the corresponding asymptotic variance of the MLE.
The optimal CSALT plans for different settings are presented in Table 1. Here, TC represents the total cost of the optimal CSALT plan, V = V ln(µ D ) is the theoretical variance of the MLE for ln(µ D ) in (7) at the normal operating condition s 0 = 25, MSE is the mean square error (based on 1000 simulated samples) for the MLE of ln(µ D ).
First, the procedure guarantees that the total costs of the determined CSALT plans are under the corresponding budgets. Second, the theoretical and empirical variances of the MLE (V and MSE) are similar in all the considered cases. In general, the inspection times of the optimal CSALT plans are robust. When more devices are tested under CSALTs with a larger budget, it is not surprising that the asymptotic variance of the MLE becomes smaller. More importantly, the optimal CSALT plans require the lowest and highest stress levels. It is also observed that a long-term experiment reduces the asymptotic variance of the MLE significantly. It is recommended that more devices be placed at the lowest stress level. The dependence influences the proportions of devices being placed at the lowest and highest stress levels. When the components in devices are nearly independent, the optimal CSALT plan with a high proportion of devices placed at the lowest stress level is obtained. Table 1. Optimal CSALTs with various budgets (C budget ), termination times (τ ter ), and degrees of dependence (β), along with the corresponding asymptotic variance (V) and MSE for ln(µ D ) at normal operating conditions s 0 = 25.

Setting
Optimal CSALT (ξ) In addition, we consider the income from the recycling of components under the optimal CSALT plans. It is assumed that the amounts of income generated from the recycling of components C1, C2, C3, and C4 are $100, $50, $30, and $10, respectively. For example, if failure modes X = {1, 3} are identified from a device, components C2 and C4 can then be re-sold, and the corresponding revenue from the recycling of the device is $60. However, for X = ∅, those devices are destroyed after successful tests, and the corresponding revenue is zero. Table 2 reports the theoretical mean and standard derivation (µ R , σ R ) and the empirical mean and standard derivation (x R , s R ) of revenue from recycling the components under various settings. It is realized that the theoretical and empirical results are similar. Moreover, more revenue would be generated from recycling the malfunctioned devices with independent components. Table 2. Optimal CSALTs with various budgets (C budget ), termination times (τ ter ), and degrees of dependence (β), along with the corresponding theoretical mean and standard derivation (µ R , σ R ) and the empirical mean and standard derivation (x R , s R ) of revenue from recycling the components.

Setting
Optimal

Eye Data from Diabetic Retinopathy Study
In this section, we consider the diabetic retinopathy study in [38]. This study examined the effectiveness of laser photocoagulation treatment in delaying the onset of blindness in patients with diabetic retinopathy. These data contain age, and the patients are classified into juveniles (aged below 20) and adults (aged 20 or above). For each patient, one eye was randomly selected for treatment and the other eye did not receive any treatment. Therefore, M = 2 and we let T 1 and T 2 represent the time to blindness for treated and untreated eyes, respectively, and s is the age covariate that influences the time to blindness. The study showed that there is a moderate association between the treated and untreated eye for each patient. As β is bounded between 0 and 0.5 for finite mean and variance, it is therefore reasonable to set β = 0.4. Table 3 presents the sample mean times to blindness for treated and untreated eyes of juvenile and adult patients as well as the sample means of age for juvenile and adult patients. With (4), it is reasonable to set (a 10 , a 11 ) = (−3, −0.006) for treated eyes and (a 20 , a 21 ) = (−3, 0.003) for untreated eyes, so that the mean times to vision loss of treated and untreated eyes for juvenile patients (s 1 = 10) are 35.55 and 32.49 months, respectively, while the mean times for adults patients (s 2 = 35) are 41.30 and 30.14 months, respectively. It is evident that the mean times to vision loss based on these settings are close to the sample means reported in Table 3. Suppose that eye doctors check the vision of both eyes for each patient only once after receiving the treatment. There are four possible outcomes for each patient, namely, no vision loss in both eyes, vision loss in the treated eye alone, vision loss in the untreated eye alone, and vision loss in both eyes. Suppose that there is no operation cost for juvenile and adult patients (C 1,oper = C 2,oper = $0), the cost of eye doctors visiting each patient is C item = $100, the budget is C budget = $20,000, and the study allows 5 years of follow-up (τ ter = 60). The proposed procedure can then be used to schedule the appointments for the patients so that the asymptotic variance of the MLE of mean time to blindness for patients of age 25 (s 0 = 25) is minimized. The optimal appointment schedule is presented in Table 4, which suggests that juvenile and adult patients are to be visited after receiving the treatment for 5 years and the ratio of juvenile patients to adult patients is 2:3. The asymptotic variance V ln(µ D ) = 0.025.

Concluding Remarks
In this paper, we have considered a one-shot device with dependent components having exponential lifetime distribution under gamma frailty and developed a procedure for obtaining the optimal CSALT plan with budget and time constraints. The optimal CSALT plan minimizes the asymptotic variance of the MLE for the mean lifetime under normal operating conditions. The simulation results showed that the procedure is reliable to obtain CSALT plans under specified budgets and experimental duration. Under the exponential distribution, the optimal CSALT plan requires only two stress levels. More devices are recommended to be placed at the lowest stress level. In addition, a short-term CSALT is not recommended because it results in a large variation of the MLE. The R codes can be available from the author upon request.
On the other hand, we investigate the impact of dependence between components in a device on the optimal CSALT plans. It is observed that the inspection times are robust, but more devices are placed at the highest stress level when the components are highly dependent. The dependence also influences the revenue from the recycling of the malfunctioned devices. It is found that more revenue would be generated when the components are nearly independent. Furthermore, as CSALT plans with two stress levels are the optimal design, it is natural to determine the explicit forms of the proportions of devices placed at the lowest and highest stress levels.
In this study, it is assumed that a practitioner can determine which components of the failed device have malfunctioned after a careful investigation or autopsy of the device. In some cases, some malfunctioned components may not be identified or some failure modes of interest are hindered by other failure modes, which lead to masked data or competing risk data. In this regard, it will be of practical interest to consider the analysis of one-shot devices based on masked data or competing risk data.
In addition, it will be of great interest to study optimal CSALT plans under more flexible lifetime distributions for components such as Weibull and gamma. Furthermore, the replacement of conditional independence given frailty with a conditional copula function described in [16] and the consideration of a log-linear model for the frailty variance as in [39] can be treated as extensions of this current work. We are currently working on these problems and hope to report the findings in a future paper.

Conflicts of Interest:
The author declares no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Abbreviations
The following abbreviations are used in this manuscript: ALT accelerated life test CSALT constant-stress ALT MLEs maximum likelihood estimates pdf probability density function MSE mean square error