Reliability Engineering and System Safety 112 (2013) 38–47
Contents lists available at SciVerse ScienceDirect
Reliability Engineering and System Safety journal homepage: www.elsevier.com/locate/ress
A Bayesian reliability evaluation method with integrated accelerated degradation testing and field information Lizhi Wang a, Rong Pan b,n, Xiaoyang Li a, Tongmin Jiang a a b
School of Reliability and System Engineering, Beihang University, Beijing, PR China Ira A. Fulton School of Engineering, Arizona State University, Tempe, AZ, USA
a r t i c l e i n f o
a b s t r a c t
Article history: Received 30 March 2012 Received in revised form 17 September 2012 Accepted 25 September 2012 Available online 29 November 2012
Accelerated degradation testing (ADT) is a common approach in reliability prediction, especially for products with high reliability. However, oftentimes the laboratory condition of ADT is different from the field condition; thus, to predict field failure, one need to calibrate the prediction made by using ADT data. In this paper a Bayesian evaluation method is proposed to integrate the ADT data from laboratory with the failure data from field. Calibration factors are introduced to calibrate the difference between the lab and the field conditions so as to predict a product’s actual field reliability more accurately. The information fusion and statistical inference procedure are carried out through a Bayesian approach and Markov chain Monte Carlo methods. The proposed method is demonstrated by two examples and the sensitivity analysis to prior distribution assumption. & 2012 Elsevier Ltd. All rights reserved.
Keywords: Degradation analysis Information fusion model Bayesian inference Sensitivity analysis
1. Introduction A product’s lifetime or degradation data is the major data source for inferring its failure time distribution model. However, it is very time-consuming to obtain enough failure data from products operated under the field-use condition. To address this problem, accelerated life tests are widely used to obtain product failure information in a short time frame. But for highly reliable products, even accelerated life testing (ALT) may be inadequate due to a large amount of censoring. Therefore, accelerated degradation testing (ADT), in which product degradation data are collected and analyzed for reliability extrapolation, becomes a common approach to the reliability prediction for those highly reliable products. However, the extrapolated reliability measurement using ADT only represents the reliability under the laboratory testing condition, which could be quite different from the product’s field condition. For example, in ADT the testing tool cannot reproduce all the stresses and stress variations that a product will experience in the field. Therefore, one should pay close attention to the model-based extrapolation that has to be carried out for inferring the field reliability. The field failure data are usually sparse, but the laboratory testing data could be abundant, thus it is meaningful to establish the correlation between the laboratory and field data and
n
Corresponding author. Tel.: þ1 480 965 4259. E-mail address:
[email protected] (R. Pan).
0951-8320/$ - see front matter & 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ress.2012.09.015
to predict the field reliability using the information from both sources. This paper aims to develop a reliability prediction model with data from multiple data sources. A brief literature review on ADT and ADT data analysis is given in the next section. Section 3 describes the degradation model, lifetime distribution and acceleration model considered in this paper. Section 4 provides the Bayesian inference approach to model parameter estimation. A synthetic example is used to demonstrate the proposed method in Section 5, along with the sensitivity analysis of this method. Section 6 describes a real product test and data analysis. Finally, the paper is concluded in Section 7.
2. Literature review Accelerated degradation testing is an effective approach to the product reliability prediction, especially for a product with high reliability and long lifetime, where even accelerated life testing is difficult to obtain sufficient failure data within the limited time and budget. Same as ALT, ADT analysis is a method based on the assumption of identical failure mechanism; that is, the life characteristic of the product under the use stress level can be extrapolated from the degradation data obtained from the accelerated stress level using a life-stress acceleration model. In this kind of test, the stress is applied in the same way as in ALT, but it is not necessary to observe the occurrence of failure event; instead, the degradation measurement of some life characteristic of interest is monitored.
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
In literature, Nelson [1] described a degradation model and the statistical method for analyzing age-degradation data on product performance. Meeker and Hamada [2] presented a conceptual degradation-based reliability model to describe the role of, and need for, integration of reliability data sources. Meeker and Escobar [3] applied the Arrhenius acceleration model on the ADT analysis. Bae et al. [4] investigated the link between a practitioner’s selected degradation model and the resulting lifetime model. Lu [5] assumed that the degradation process follows a Wiener process, so the first passage time distribution to a barrier is distributed as an inverse Gaussian distribution. Pettit and Young [6] extended the analysis in Lu by using a fully Bayesian approach to model estimation and prediction. Barker and Newby [7] addressed the problem of determining inspection and maintenance strategy for a system whose state is described by a multivariate Wiener process. Carey and Koenig [8], Tseng and Yu [9], Tseng et al. [10] and Onar and Padgett [11,12] presented their relevant methods of ADT based on film resistance, MOS, logic chip, LED, carbolic fiber, respectively. Unlike the field-use condition, the stresses applied on the test units in a laboratory testing condition can be precisely controlled, so the extrapolated result from an ADT may not accurately represent the real situation in the field, where most stress variables are in fact random. To solve this problem, Liao and Elsayed [13] extended an ADT model obtained from constantstress ADT experiments to predict field reliability by considering the field stress variations. Meeker et al. [14] developed a use rate model to relate ALT data and field data. Monroe and Pan [15] made an effort of using global climate data to adjust the acceleration factor for any specific market place of electronic products. However, it is difficult to predict all stresses and stress variations that the product will experience under its field condition. Integrating reliability information from various sources is crucial to the solution and Bayesian data analysis provides a feasible means for information fusion. In a series of paper by Hamada and his colleagues [16–18] the methods of using hierarchical Bayes model to integrate reliability information were investigated. Touw [19] studied the Bayesian estimation of mixed Weibull distribution and others. Pan [20] introduced a calibration factor to lifetime model to make a more accurate prediction of field failure time by considering the data from both ALT and field. This calibration factor statistically calibrates the ALT model’s parameter and prediction result without explicitly considering specific stress or influencing factors. In this paper we formulate calibration factors for the ADT model and use Bayesian methods to integrate ADT and field data for reliability prediction.
3. Models 3.1. Degradation model We assume that the degradation follows a Wiener process, and then the first passage time of this process to a threshold follows the inverse Gaussian distribution. Chhikara and Folks [21] discussed the use of this distribution as a lifetime distribution model. Whitmore and Schenkelberg [22] and Lu [5] had considered it for degradation data analysis. The degradation model of Wiener process is YðtÞ ¼ sBðtÞ þ dðsÞt þy0 ,
ð1Þ
where Y(t) is the performance degradation process of product, B(t) is the standard Brownian motion with the mean value of 0 and the variance of t, denoted as B(t) N(0, t), s is the diffusion coefficient, which is a constant and does not change with stress or time, y0 is a known initial value of product performance, and d(s) is the drift
39
coefficient, which represents the degradation rate of product and is a function of stress. The acceleration model of degradation rate is given by dðsÞ ¼ exp ½b0 þ bjðsÞ,
ð2Þ
where j(s) is a function of stress factor s. Known from the property of the Wiener process, the degradation increment Dy during the unit time Dt is subject to a normal distribution with the mean of d(s)Dt and the variance of s2Dt, i.e.,
Dy N ðdðsÞDt, s2 DtÞ
ð3Þ
Let l be the failure threshold and the product fails when y(t) lo0, then the first passage time t to the threshold l has an inverse Gaussian distribution. The probability density function of an inverse Gaussian distribution is rffiffiffiffiffiffiffiffiffiffiffi 2 l l exp 2 xm , ð4Þ f x; m, l ¼ 3 2px 2m x and the cumulative distribution function is rffiffiffi rffiffiffi ! ! l x 2l l x F ðxÞ ¼ PðX r xÞ ¼ F 1 þ exp F þ1 : x m m x m ð5Þ For a Weiner process, its first passage time t follows an inverse Gaussian distribution with m ¼((l y0)/(d(s))) and l ¼(l y0/s)2, thus its reliability function is given by
2dðsÞ ly0 ly0 dðsÞt ly0 þ dðsÞt pffiffi p ffiffi Rðt Þ ¼ F F exp : ð6Þ s2 s t s t From the reliability function above, we know that the key parameters for reliability prediction are d(s) and s2. After these parameters are obtained, the reliability prediction can be made. 3.2. Calibration factor The purpose of ADT is to predict the lifetime and the reliability of a product deployed in the field. However, the extrapolated reliability value from ADT data only represents the reliability in the laboratory environment, which may be quite different from the field environment. For example, some stresses and variability of stresses, experienced by the product in its field, cannot be simulated in a testing laboratory and it is difficult to replicate the environmental noise in a lab too. As d(s) and s2 are the key parameters in the degradation model, we introduce two calibration factors, k1 and k2, to update the estimates of these two parameters so as to solve the problem above. Calibration factor k1 is introduced to take into the consideration of additional environmental stresses and their effects on the degradation path, d(s). For instance, there are multiple active stress variables, s1, s2y sm, under the field condition, but only the stress s1 can be simulated in ADT, while other stresses are ignored. In this paper, assuming these stress variables are independent, the field acceleration model is df ðSÞ ¼ exp½b0 þ b1 j1 ðs1 Þ þ b2 j2 ðs2 Þ þ UUU, þ bm jm ðsm Þ ¼ exp½b0 þ b1 j1 ðs1 ÞUexp½b2 j2 ðs2 Þ þ UUU, þ bm jm ðsm Þ:
ð7Þ
But the ADT acceleration model is dADT ðs1 Þ ¼ exp½b0 þ b1 j1 ðs1 Þ: Let the effect of other stresses be represented by the calibration factor k1, then k1 is defined as k1 ¼ exp½b2 j2 ðs2 Þ þUUU, þ bm jm ðsm Þ,
ð8Þ
and we have df ðSÞ ¼ k1 UdADT ðs1 Þ:
ð9Þ
40
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
The exponential acceleration function used in (7) (or the loglinear relationship between stress and accelerated parameter) is very common in reliability analysis, but it is possible to use other type of acceleration model depending on specific failure mechanism and the definition of calibration factor will be adjusted accordingly. Although the parameter s2 is not directly related to stress variables, it governs the variability of lifetime distribution. Considering that in the field environment the deployed products can be different from those tested in laboratory, due to, for examples, manufacturing variations, changed assembling conditions, etc., we introduce another factor to calibrate this parameter. The degradation increments in ADT, DyADT, and under the field condition, Dyf, are, respectively,
DyADT ¼ dADT ðsÞUDt þ eADT ,
ð10Þ
and
Dyf ¼ df ðSÞUDt þ ef :
ð11Þ 2 N 0, sADT , the calibration factor, k2, is additive to s2,
As eADT i.e., ef N 0, s2ADT þk2
ð12Þ
or
s2f ¼ s2ADT þ k2 :
ð13Þ
3.3. Extension to gamma process Although we assume that the degradation process follows a Wiener process above, our calibration method of ADT inference can be easily extended to other stochastic processes as long as the distribution function of degradation increment and the degradation model are well-defined. Gamma process is a stochastic process that is always positive and strictly increasing [23]. It can be used to model a degradation in which the system deterioration is measured and the deterioration cannot be reversed. A Gamma process {Y(t),tZ0} has independent, non-negative increments Dy¼ Y(tþ t)Y(t) that follow a gamma distribution as [24,25]
Dy Gðb, aðt þ tÞaðtÞÞ,
ð14Þ
where b (b 40) is a constant scale parameter and a(t) is a timevarying shape parameter with a(0)¼0. Let T be the first passage time to the failure threshold l and y0 be the initial value of the process, then the reliability function is given by P ðT 4 t Þ ¼ P ðY ðt Þ olÞ ¼ P Y ðt Þy0 o ly0 Z ly0 1 y yaðtÞ1 exp dy: ð15Þ ¼ b 0 GðaðtÞÞbaðtÞ The shape parameter a(t) is assumed to be a function of stress variable in Park and Padgett [23] and Pan and Balakrishnan [24] as
aðtÞ ¼ dðsÞUt,
ð16Þ
where d(s) is the drift coefficient, representing the degradation rate as defined by the acceleration model in Eq. (2). So, the degradation increment Dy during the unit time Dt is subject to a gamma distribution with the shape parameter of d(s)Dt and the scale parameter of b, i.e.,
Dy Gðb,dðsÞDtÞ:
ð17Þ
When we take into consideration of additional environmental stresses and their effects on the degradation path, we can introduce a calibration factor k1 into the field drift coefficient, df(s), as in Eq. (9). The parameter b of gamma distribution is a scale parameter and when b is small the distribution is more concentrated. We can
introduce another calibration factor k2 to the effect of the field environmental noise on this parameter, then
bf ¼ k2 UbADT
ð18Þ
Now we have DyADT G(bADT,d(s)Dt) and Dyf G(k2bADT,k1d(s)Dt) and the lab testing data and field data can be combined together to infer the model parameters of interest.
4. Bayesian inference To understand the relationship between different environmental conditions and to estimate the reliability accurately, it is important to evaluate the parameters, such as b0, b1 and s2, and the calibration factors k1 and k2. It is also a difficult task when the model is complex, with a high dimensional parameter space. Bayesian approach is a feasible method to integrate all available information to infer unknown parameters. In Bayesian approach parameters are treated as random variables and their probabilistic models are obtained by posterior distributions. Markov Chain Monte Carlo (MCMC) is a simulation approximation method that is widely used in Bayesian inference. By this method sampling of an unknown variable is implemented through establishing an appropriate Markov chain, as the stationary distribution of the Markov chain converges to the variable’s posterior distribution. A special case of MCMC is the Gibbs sampler, which involves simulating from the conditional posterior distribution of each parameter given the data and all of the other parameters. It is effective in handling high dimensional problems and has been applied on various reliability engineering applications. 4.1. Degradation data from both ADT and field When it is possible to collect degradation data, yfield, from field, we may combine them with the ADT data, yADT, and model them by
DyA NðdA ðsÞUDtA , s2A UDtA Þ
ð19Þ
and
Dyf N df ðs UDtf , s2f UDtf Þ:
ð20Þ
Applying the two calibration factors, it becomes Dyf Nðk1 UdA ðsÞUDtf , s2A þ k2 UDt f Þ: In general, the degradation data is a vector, including Dyi, si and ci, where Dyi is degradation increment, si is the stress level, ci is the condition indicator variable (if it is the data from ADT, ci ¼0; if it is from field, ci ¼1), so a combined degradation model is as Dyi N m0 , s20 , ð21Þ where
m0 ¼ dA ðsi ÞU Dt A þci k1 UDtf DtA Þ, and
s20 ¼ s2A UDt A þci k2 UDtf þ s2A UDtf s2A UDt A : Then the posterior distribution of unknown parameters is p Y9D ¼ p b0 , b1 , s2A ,k1 ,k2 9Dy,s,c !1=2 ð22Þ 2 ! n Y Dyi m0 1 p exp Up b0 , b1 ,k1 ,k2 , s2A 2 2 s0 2Us0 i¼1
where n is the number of degradation observations, pðb0 , b1 ,k1 ,k2 , s2A Þ is the joint prior distribution of b0, b1, k1, k2 and s2A .
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
4.2. Degradation data from ADT and failure time data from field Sometimes continuous product monitoring cannot be implemented in field, so field degradation data are not available; instead, we may have a few of field failure time data. Based on the Weiner process model, these field failure time data, X ¼ x1 ,x2 ,. . .,xn1 , are inverse Gaussian distributed, so the probability density function is !1=2 ( ) 2 ldf ðs xÞ2 l hðxÞ ¼ exp : ð23Þ 2px3 sf 2 2xsf 2 Thus, !1=2
2
hðxÞ ¼
l 2px3 sA 2 þk2
( exp
2
ðlk1 UdA ðsÞxÞ 2x sA 2 þ k2
) :
ð24Þ
From (3), the degradation increment data DY ¼ Dy1 , Dy2 ,. . ., Dyn2 Þ from ADT has a normal distribution with the following probability density function, ! 1 ðDydA ðsÞDtÞ2 f ðDyÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi exp : ð25Þ 2sA 2 Dt 2psA 2 Dt So the total likelihood of both ADT and field data is LðdðsÞ, s2 ,k1 ,k2 Þ ¼
n1 Y
5. A synthetic example and sensitivity analysis 5.1. Prior distribution selection As the nature of prior distribution defined in Bayesian inference is subjective, it is important to select and compare feasible prior distributions. For b0, b1and s2 are parameters in a normal regression model, according to Ntzoufras [26], the simplest approach is to assume that these parameters have independent priors with p Y f bj f s2 , f b, s2 ¼ j¼0
bj N mbj , e2j ,
s2 IGaða,bÞ,
j¼1
where n1 is the number of field failure time observations and n2 is the number of ADT degradation increment observations. To facilitate the MCMC computation, a zeros-ones trick [26] of defining an arbitrary likelihood function in WINBUGS is applied. Using this method, either the Bernoulli or the Poisson distribution can be utilized to indirectly specify any arbitrary likelihood function. Assume there are n independent observations and each of them is from a Poisson distribution with the Poisson parameter being the negative value of log-likelihood of a single failure time or degradation increment, then the likelihood function can be written as 0 n n n Y Y Y eðli Þ ðli Þ L z0 9y ¼ ¼ eli ¼ f P ð0; li Þ 0! i¼1 i¼1 i¼1
As one can see from (22) and (28), they are very complex model and would require the MCMC and Gibbs algorithms to obtain the posterior distribution of the parameters of interest with WinBUGS. The reliability inference can be obtained by (9), (13) and (6) after the model parameter estimation is conducted.
and
n2 Y hðxi Þ f Dyj :
i¼1
41
ð26Þ
Hence, the model likelihood can be written as the product of the densities of new pseudo-random variables (all observations are zeros) that follow the Poisson distribution with mean equal to the negative log-likelihood value. Let lf i ¼ logh xi 9b0 , b1 , s2A ,k1 ,k2 ,
that is, bj follows the normal distribution and s2 follows the inverse gamma distribution. For the calibration factor k1 and k2, either normal or gamma distribution can be a feasible solution. We test four sets of candidate prior distributions with the equal mean values. These prior distributions are shown in Table 1. The deviance information criterion (DIC) is typically used as a measure in model comparison, but DIC must not be used when the posterior distributions are highly skewed or bimodal (Ntzoufras [26]). We also find that it is not a suitable method in this application. To select the best prior distributions of k1 and k2, we choose a simulation study and analyze the accuracy of the estimated parameter values with their assumed values used in the model. First, the data generated from an assumed model are used to estimate the parameters, b0, b1, k1, k2 and s2, by the method proposed in this paper with different k1 and k2 prior distributions. Then, we evaluate the accuracy of these estimates. Finally, the best k1 and k2 prior distributions are selected. The analysis would be processed under the situation of degradation data from both ADT and field, and the situation of degradation data from ADT and failure time data from field, respectively.
and lAj ¼ logf Dyj 9b0 , b1 , s2A : Let ci be the testing condition indicator variable, then li ¼ ci Ulf i þ ð1ci ÞUlAi : The data include zi, si and ci, where zi represents the failure time data ti or the degradation increment data Dyi. Then the model likelihood is þ n2 n1Y þ n2 ðli Þ 0 n1Y þ n2 n1Y e ðli Þ L z,s,c9b0 , b1 , s2A ,k1 ,k2 ¼ eli ¼ f P ð0; li Þ: ¼ 0! i¼1 i¼1 i¼1
The posterior distribution is p Y9D ¼ p b0 , b1 , s2A ,k1 ,k2 9z,s,c pL z,s,c9b0 , b1 , s2A ,k1 ,k2 Up b0 , b1 ,k1 ,k2 , s2A ,
ð27Þ
5.1.1. Degradation data from both ADT and field Suppose there is a temperature step-stress ADT applied to a product, where there are 4 test units and 4 temperature levels (60 1C, 80 1C, 100 1C, 120 1C). The time of each level is 1250, 750, 500 and 500 h, and the inspection time interval is 5 h, so there are 600 (250þ150þ100þ100¼ 600) observations for each test unit and 2400 observations in total. There is 1 field sample which is tested at 25 1C over 5000 h and the inspection time interval is also 5 h, so there are 1000 observations. Suppose that the degradation of this product obeys the Weiner process. The initial performance value y0 is 100 and the failure Table 1 Candidate prior distribution.
ð28Þ
where pðb0 , b1 ,k1 ,k2 , s2A Þ is the joint prior distribution of b0, b1, k1, k2 and s2A .
k1 k2
1
2
3
4
Normal Normal
Normal Inverse-gamma
Gamma Inverse-gamma
Gamma Normal
42
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
threshold l is 50. As the stress factor is temperature, the Arrhenius acceleration function is applicable, so j(s)¼1/T, where T is the absolute temperature, and dðTÞ ¼ exp½b0 þ b1 =T: In the simulation, the true values of b0, b1, k1, k2 and s are given in Table 2. To compare the evaluation results without other interference, all prior distributions are defined as having their mean values as shown in Table 3, and their variances are equal to 105. A total of 50 sets of simulation data are generated, and for every set of prior distributions, the posterior estimations of these parameters are computed. When the MCMC simulation starts, checking the distribution convergence from the sampled values of the parameter posterior distribution is necessary. The Gelman–Rubin ratio Table 2 True value of model parameters.
(Gelman and Rubin [27]) is taken to monitor the convergence, it is the value of width of the central 80% interval of the pooled runs over the average width of the 80% intervals within individual runs. We start with two Markov chains with different initial values. When the Markov chain simulation sample size increases, if the process converges, the Gelman–Rubin ratio converges to 1. We run 40,000 iterations for each Markov chain in the simulation. Fig. 1 shows the Gelman–Rubin ratio of b0, b1, k1, k2 and s from one set of simulation data. One can see that the Gelman–Rubin statistics of all parameters become stable and converge to 1 after 15,000 iterations, so the remained 25,000 iterations of each chain (50,000 in total) are selected as the samples from each parameter’ posterior distribution. Fig. 2 provides the histograms of posterior distributions. The average relative errors of these 5 parameters are shown in Table 4. As b0 and b1 are parameters in d(s), the accuracy of d(s) is more meaningful than these parameters. The relative errors of d(s) are included in Table 4, where s is set to 25 1C. The relative error is relative error ¼ 9evaluation valuetrue value9=true value
Parameter
b0
b1
s
k1
k2
True value
9.2103
4800
0.025
1.5
6e-4
ð29Þ
From Table 4, one can see that the gamma distribution and inverse gamma distribution for k1 and k2, respectively, is the best Table 4 Relative error.
Table 3 Mean value of prior distributions. Parameter
b0
b1
s
k1
k2
True value
9.2103
4800
0.025
1
0
1 2 3 4
b0
b1
d (25 1C)
k1
k2
s
0.126 0.141 0.138 0.137
0.0099 0.0100 0.0105 0.0106
0.044 0.059 0.040 0.041
0.284 0.282 0.219 0.221
0.129 0.061 0.061 0.130
0.037 0.032 0.032 0.038
Fig. 1. Monitoring MC convergence by using Gelman–Robin statistic (red line). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)
Fig. 2. Posterior distributions of b0, b1, k1, k2 and s.
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
In this case, there is not a single prior candidate that is the best at every parameter estimation. We compute the reliability values at the product’s third year in field use for the 4 sets of prior distribution using the 50 sets of simulation data and compare their relative errors. The average relative errors are shown in Table 6. Here we find that normal distribution as the prior distribution for k1 and k2 is the best. Let k1 N mk1 , s2k1 k2 N mk2 , s2k2 ,
candidate. Let the prior distribution of k1 and k2 be k1 G ak1 ,bk1 k2 IGa ak2 ,bk2 , then (22) can be written as p Y9D ¼ p b0 , b1 , s2A ,k1 ,k2 9Dy,s,c !1=2 2 ! n Y Dyi m0 1 p exp s20 2Us20 i¼1 b m b m bak1 0 1 k b0 b1 a þ1 Uf U 1 ðk1 Þ k1 Uf c0 c1 G ak1
ak2
bk 1 ebk1 k1 U 2 G ak2 k2
ak
2
þ1
a
b 1 ebk2 k2 U GðaÞ s2A
!a þ 1 eb=s
43
2 A
5.1.2. Degradation data from ADT and failure time data from field We follow the previous example, but let the field information be 20 failure times. The true values of model parameters b0, b1, k1, k2 and s used in the simulation are the same as in Table 2. The average relative errors of the 5 parameters, along with the relative errors of d(s), are shown in Table 5.
then (28) can be written as p Y9D ¼ p b0 , b1 , s2A ,k1 ,k2 9z,s,m pL z9b0 , b1 , s2A ,k1 ,k2 ! ! b0 mb0 b1 mb1 k1 mk1 k2 mk2 Uf Uf Uf Uf c0 c1 sk1 sk2 a
U
b 1 GðaÞ s2A
!a þ 1 2
eb=sA
Table 5 Relative error.
1 2 3 4
b0
b1
d (25 1C)
k1
k2
s
0.119 0.1179 0.1129 0.1169
0.0094 0.0099 0.0091 0.0095
0.038 0.042 0.039 0.041
0.0409 0.0411 0.0384 0.0424
0.0036 0.00053 0.00052 0.00058
0.0084 0.0508 0.0514 0.0512
Table 6 Relative error of reliability (3 years).
Relative error
1
2
3
4
0.0263
0.0524
0.0495
0.0489
Fig. 4. Sensitivity plot of reliability for different values of b0 prior mean.
Fig. 3. Sensitivity plot of b0 posterior mean for different values of b0 prior mean. (For interpretation of the references to color in this figure, the reader is referred to the web version of this article.)
44
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
Fig. 5. Sensitivity plots for different values of b1 prior mean.
Fig. 6. Sensitivity plots for different values of k1 prior mean.
Fig. 7. Sensitivity plots for different values of k2 prior mean.
5.2. Sensitivity analysis In this section, sensitivity analysis is implemented with different values of the prior mean to study the robustness of the Bayesian method. There are 5 parameter prior means, thus we carry out an orthogonally designed experiment to study their influence on posterior evaluation. Again, this analysis would be processed under the situation of degradation data from both ADT and field, and the situation of degradation data from ADT and failure time data from field. When the degradation data are from ADT and field, followed the example in Section 5.1.1, we vary the prior mean value of the parameter mb0 and observe changes in its posterior mean, which is
shown in Fig. 3 (the blue line is posterior mean, the green line is the posterior quantile value at 97.5%, the red line is the posterior quantile value at 2.5%). One can see that the posterior mean is quite robust with values from 9.1 to 9.4, when the prior values are from 8.6 to 9.6. Then we study the sensitivity of the reliability evaluation result (3 years) over different values of the prior parameter mb0 , the sensitivity is shown in Fig. 4, one can see that the reliability evaluation value is also robust with values from 0.92 to 0.93, when the prior values vary from 8.6 to 9.6. Same to mb0 , the sensitivity of the posterior mean and reliability over different values of b1, k1, k2 and s’s prior means are shown in Figs. 5–8. We can see that the posterior mean of b1 is from 4770 to 4820, reliability evaluation value is from 0.93 to 0.94, when the prior
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
45
Fig. 8. Sensitivity plots for different values of s prior mean.
Table 7 L16(4^5) orthogonal array of sensitivity analysis for the example in Section 5.1.1.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 Ij IIj IIIj IVj Rj
b0
b1
k1
k2
s2
1
2
3
4
5
120% 120% 120% 120% 80% 80% 80% 80% 140% 140% 140% 140% 60% 60% 60% 60% 0.0914 0.146 0.1306 0.0824 0.0636
120% 80% 140% 60% 120% 80% 140% 60% 120% 80% 140% 60% 120% 80% 140% 60% 0.1215 0.0919 0.0944 0.1426 0.0507
120% 80% 140% 60% 80% 120% 60% 140% 140% 60% 120% 80% 60% 140% 80% 120% 0.0809 0.0365 0.2390 0.0940 0.2025
120% 80% 140% 60% 140% 60% 120% 80% 60% 140% 80% 120% 80% 120% 60% 140% 0.0838 0.1396 0.0947 0.1323 0.0558
120% 80% 140% 60% 60% 140% 80% 120% 80% 120% 60% 140% 140% 60% 120% 80% 0.14 0.1266 0.1045 0.0793 0.0607
Table 8 L16(4^5) orthogonal array of sensitivity analysis for the example in 5.1.2.
REOR (3 years)
0.0168 0.0089 0.0434 0.0223 0.0069 0.0263 0.0236 0.0892 0.0749 0.0252 0.0186 0.0119 0.0229 0.0315 0.0088 0.0192 T¼ 0.4504
values are from 2800 to 6800; the posterior mean of k1 is from 1.3 to 1.6, reliability evaluation value is from 0.93 to 0.94, when the prior values are from 0.9 to 2.1; the posterior mean of k2 is from 5e-4 to 7e4, reliability evaluation value is from 0.915 to 0.925, when the prior values are from 1e-4 to 1.2e-3; the posterior mean of s is from 0.25 to 0.26, reliability evaluation value is from 0.925 to 0.935, when the prior values are from 0.012 to 0.035. From here, we conclude that the reliability prediction is robust to a certain extent, but the degrees of influence of the five parameters are different. To which parameter that the reliability prediction is more sensitive is a question we should discuss next. To study the influence of the priors and the sensitive difference among these five parameters, a contrast test is applied, where the relative error of reliability (REOR) is computed with Eq. (29). The mean values of b0, b1, k1, k2 and s are the factors, and we let the 120%, 80%, 140% and 60% of their true values be the 1st, 2nd, 3rd and 4th level of these factors. If we pair every factor to every level, there will be 4^5 ¼1024 combinations. A L16(4^5) orthogonal array [28] is selected so as to reduce the number of experiments to 16 and the sensitivity analysis of different parameters can still be carried out by the analysis of range method. Using the example in Section 5.1.1, the results are shown in Table 7.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 I II III IV R
b0
b1
k1
k2
s2
Relative error of reliability (3 years)
120% 120% 120% 120% 80% 80% 80% 80% 140% 140% 140% 140% 60% 60% 60% 60% 0.0560 0.0481 0.0420 0.0554 0.0140
120% 80% 140% 60% 120% 80% 140% 60% 120% 80% 140% 60% 120% 80% 140% 60% 0.0578 0.0374 0.0586 0.0477 0.0212
120% 80% 140% 60% 80% 120% 60% 140% 140% 60% 120% 80% 60% 140% 80% 120% 0.0403 0.0607 0.0510 0.0495 0.0204
120% 80% 140% 60% 140% 60% 120% 80% 60% 140% 80% 120% 80% 120% 60% 140% 0.0380 0.0585 0.0468 0.0582 0.0205
120% 80% 140% 60% 60% 140% 80% 120% 80% 120% 60% 140% 140% 60% 120% 80% 0.0618 0.0169 0.0923 0.0305 0.0754
0.0131 0.0070 0.0258 0.0101 0.0106 0.0187 0.0024 0.0164 0.0061 0.0090 0.0071 0.0198 0.0280 0.0027 0.0233 0.0014 T¼ 0.2015
In Table 7, Ij, IIj, IIIj and IVj are the sum of 1st, 2nd, 3rd and 4th level’s REOR values in parameter j(j¼1,y,5), and T is the sum of the 16 REOR values. Take I1 for example, which is the sum of 1st level’ REOR in the column of b0, so the corresponding REOR values are 0.0168, 0.0089, 0.0434 and 0.0223, and I1 ¼ 0:0168 þ 0:0089 þ 0:0434þ 0:0223 ¼ 0:0914 Rj is the range of parameter j, i.e., Rj ¼ maxðIj ,IIj ,IIIj ,IVj Þmin Ij ,IIj ,IIIj ,IVj The value of Rj indicates the influence of the changed levels to the target (REOR) for this factor, so it stands for the evaluation of the sensitivities of REOR to different parameters. From Table 6 one can see that the ranking of sensitivities is as k1 4 b0 4 s 4k2 4 b1; therefore, when one should exercise more cautions when to assign the prior distributions for those parameters that the reliability prediction is sensitive to. The same type of analysis is performed on the example in Section 5.1.2, when the degradation data from ADT and the failure data from field condition are available. The results are shown in Table 8. In this case, the ranking of the sensitivities is as s 4 b1 4k2 4k1 4 b0.
46
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
209.45 209.4
output optical power
209.35 209.3 209.25 209.2 209.15 209.1 209.05 209 208.95
Fig. 9. Structure of SLD.
0
20
40
60
80
100
120
140
160
180
200
Time/day Fig. 12. Degradation data of SLD under its working environment.
Fig. 10. Stress profile of step-stress ADT.
18 16
output optical power
14 12 10 8 6 4 2
0
500
1000
1500
2000
2500
3000
3500
4000
Time/hour Fig. 11. Degradation data of SLD.
6. A real example Super luminescent diode (SLD) is a semiconductor product that has a wide range of applications in optical fiber sensors, optical fiber communication systems, clinical systems, and so on. This product enjoys long lifetime and high reliability, as well as other design and manufacturing advantages. The structure of SLD is shown in Fig. 9. As failure time of SLD is difficult to obtain, degradation of SLD can be monitored to predict its reliability. The performance degradation measurement of SLD is its output optical power, which is sensitive to temperature. This type of SLD had been preliminary used for many years, but there was no systematical study of the effect of temperature on its lifetime. We conducted a step-stress ADT experiment with
4 progressive stress levels on 4 SLDs. The test plan was followed as 60 1C of 2400 h, 80 1C of 800 h, 100 1C of 300 h, and finally 110 1C of 300 h. The stress profile of ADT is shown in Fig. 10. The degradation data were recorded every 10 min with power tester automatically controlled with a computer. These ADT data are plotted in Fig. 11. SLD’s typical work condition is at around 25 1C. But in the field, it is difficult to measure its output optical power as frequent as in the experiment. We have measurements of one SLD every 6 h for 200 days and they are plotted in Fig. 12. One can see that they have significant larger variability than the ADT data. A hypothesis testing (p¼0.05) was applied on the degradation increments at different temperatures and the result showed that they followed normal distributions, so we elected to use the Wiener process to describe the degradation process. We also know the initial value of the degradation increment’s mean and variance from the hypothesis testing, it would be used to define the parameters’ prior values. The failure threshold value of SLD is the 60% of its initial output optical power value. Using the Weiner process model and the ADT data, we can predict that the 5 years’ reliability of SLD was 0.99. But, known from experience, the field reliability of SLD in 5 years is much lower than 0.99, it is roughly between 0.87 and 0.97. The field reliability estimated from the field degradation data is 0.93 and the 95% confidence interval is [0.89, 0.97], which is closer to engineer’s knowledge. However, using only one field sample cannot provide us an accurate result, so we combined the 200 days’ field data with the ADT data via the method in this paper, the prior of parameter bj of drift coefficient follows the normal distribution and diffusion coefficient s2 follows the inverse gamma distribution, and the parameter values of these priors are defined with the hypothesis testing results above and experience. Gamma distribution and inverse gamma distribution are selected as the prior distribution for k1 and k2, with their means to be 1 and 0, respectively. Their prior parameters are defined such that the gamma and inverse gamma distributions have such mean values and variance as 105 The calibration factor k1 was calculated to be 1.3 and k2 was 0.001, and the 5-year reliability of SLD became 0.935, with the 95% confidence interval as [0.93, 0.94]. This result is more precise than that with field data only and it is much closer to field experience than using ADT data only. In addition, with the estimated degradation model, we can study the effect of temperature on the SLD degradation rate.
7. Conclusions In this paper a Bayesian approach is proposed to integrate the product’s reliability information from both the ADT and field-use conditions. The degradation model and calibration factors are introduced to bridge the reliability evaluation from the accelerated
L. Wang et al. / Reliability Engineering and System Safety 112 (2013) 38–47
degradation testing in a testing laboratory to the product’s actual performance in the field. These calibration factors model the uncertainty of stress fluctuations and the complexity of failure behaviors of the product in its field-use environment, so our method can solve the dilemma that the abundant ADT information cannot represent a product’s reliability accurately when the modelbased extrapolation has to be conducted for reliability prediction. We provide two examples and a discussion on prior distribution selection. The sensitivity analysis gives the guideline for prior distribution assignment.
References [1] Nelson W. Analysis of performance-degradation data from accelerated tests. IEEE Transactions on Reliability 1981;30(2):149–54. [2] Meeker WQ, Hamada M. Statistical tools for the rapid development and evaluation of high-reliability products. IEEE Transactions on Reliability 1995;44(2):391–400. [3] Meeker WQ, Escobar LA. Statistical Methods for Reliability Data. 1st ed. New York: John Wiley and Sons; 1998. [4] Bae SJ, Kuo W, Kvam PH. Degradation models and implied lifetime distributions. Reliability Engineering and System Safety 2007;92:601–8. [5] Lu J, Degradation processes and related reliability models, Ph.D. thesis, McGill University; 1995. [6] Pettit LI, Young KDS. Bayesian analysis for inverse Gaussian lifetime data with measures of degradation. Journal of Statistical Computation and Simulation 1999;63(3):217–34. [7] Barker CT, Newby MJ. Optimal non-periodic inspection for a multivariate degradation model. Reliability Engineering and System Safety 2009;94:33–43. [8] Carey MB, Koenig RH. Reliability assessment based on accelerated degradation a case study. IEEE Transactions on Reliability 1991;40(5):499–506. [9] Tseng ST, Yu HF. A termination rule for degradation experiments. IEEE Transactions on Reliability 1997;46(1):130–3. [10] Tseng ST, Hamada M, Chiao CH. Using degradation data from a factorial experiments to improve fluorescent lamp reliability. Journal of Quality Technology 1995;27(4):363–9. [11] Onar A, Padgett WJ. Accelerated test models with the inverse Gaussian distribution. Journal of Statistical Planning and Inference 2000;89:119–33.
47
[12] Onar A, Padgett WJ. Inverse Gaussian accelerated test models based on cumulative damage. Journal of Statistical Computation and Simulation 2000;66(3):233–47. [13] Liao H, Elsayed EA. Reliability inference for field conditions from accelerated degradation testing. Naval Research Logistics 2006;53(6):576–87. [14] Meeker WO, Escobar LA, Hong Y. Using accelerated life tests results to predict product field reliability. Technometrics 2009;51(2):146–61. [15] Monroe E, Pan R. Knowledge-based reliability assessment for time-varying climates. Quality and Reliability Engineering International 2009;25:111–24. [16] Hamada M, Martz HF, Reese CS, Graves T, Johnson V, Wilson AG. A fully Bayesian approach for combining multilevel failure information in fault tree quantification and optimal follow-on resource allocation. Reliability Engineering and System Safety 2004;86:297–305. [17] Graves TL, Hamada MS, Klamann R, Koehler A, Martz HF. A fully Bayesian approach for combining multi-level information in multi-state fault tree quantification. Reliability Engineering and System Safety 2007;92:1476–83. [18] Graves TL, Hamada MS, Klamann R, Koehler A, Martz HF. Using simultaneous higher-level and partial lower-level data in reliability assessments. Reliability Engineering and System Safety 2008;93:1273–9. [19] Touw AE. Bayesian estimation of mixed Weibull distributions. Reliability Engineering and System Safety 2009;94:463–73. [20] Pan RA. Bayes approach to reliability prediction utilizing data from accelerated life tests and field failure observations. Quality and Reliability Engineering International 2009;25(2):229–40. [21] Chhikara RS, Folks JL. The inverse Gaussian distribution as a lifetime model. Technometrics 1977;19(4):461–8. [22] Whitmore GA, Schenkelberg F. Modelling accelerated degradation data using Wiener diffusion with a time scale transformation. Lifetime Data Analysis 1997;3(1):27–45. [23] Park C, Padgett WJ. Accelerated degradation models for failure based on geometric brownian motion and gamma processes. Lifetime Data Analysis 2005;11:511–27. [24] Pan Z, Balakrishnan N. Reliability modeling of degradation of products with multiple performance characteristics based on gamma processes. Reliability Engineering and System Safety 2011;96:949–57. [25] Guida M, Postiglione F, Pulcini G. A time-discrete extend gamma process for time-dependent degradation phenomena. Reliability Engineering and System Safety 2012;105:73–9. [26] Ntzoufras I. Bayesian Modeling using WinBUGS. 1st ed. New York: John Wiley and Sons; 2009. [27] Gelman A, Rubin DB. Inference from iterative simulation using multiple sequences. Statistical Science 1992;7:457–72. [28] The University of York [Internet]. Available from: /http://www.york.ac.uk/ depts/maths/tables/l16.htmS.