Multiple-Steps Step-Stress Accelerated Degradation Modeling

Introduction

The MultiStepSSAD package provides generalized statistical tools for modeling Multiple-Steps Step-Stress Accelerated Degradation Testing (SSADT) data, implementing the methods developed by Pan and Balakrishnan (2010).

In step-stress accelerated degradation experiments, products are subjected to elevated stress levels (\(S_1 < S_2 < \dots < S_K\)). Unlike traditional step-stress tests where stress levels are elevated at pre-determined fixed time points for all units, Pan and Balakrishnan (2010) proposed elevating stress levels when an individual product’s degradation crosses pre-specified degradation threshold values (\(\omega_1, \dots, \omega_{K-1}\)). Consequently, the stress transition times (\(\tau_{i,k}\)) vary randomly from product to product.

This package supports: 1. Wiener Process Degradation Models: Drift parameter changes with stress level; diffusion parameter \(\sigma\) remains constant. 2. Gamma Process Degradation Models: Monotone degradation process with independent Gamma-distributed increments. Stress transition times are modeled using the Birnbaum-Saunders distribution approximation. 3. Acceleration Models: Arrhenius model \(A(S) = \exp(a + b / (273 + S))\) and Power model \(A(S) = a \cdot S^b\). 4. Estimation Methods: Maximum Likelihood Estimation (MLE) and Bayesian Markov Chain Monte Carlo (MCMC). 5. Diagnostics & Lifetime Prediction: Geweke’s MCMC convergence diagnostic, Mean Time To Failure (MTTF), and reliability \(R(t)\) evaluation under normal use stress \(S_0\).


1. Package Datasets

The package includes the exact simulated datasets from Tables 2–5 of Pan and Balakrishnan (2010):

library(MultiStepSSAD)

# Load Wiener-Arrhenius dataset (Table 2)
data(wiener_arrhenius)
head(wiener_arrhenius)
#>   time   unit1   unit2   unit3   unit4   unit5
#> 1   72  4.3195  5.2906  5.5159  3.9891  5.4399
#> 2  144  8.9338 11.7680 10.4270  8.5626  9.8537
#> 3  216 14.0880 15.9770 13.9230 13.7280 15.0460
#> 4  288 18.5820 22.2880 18.0950 19.7180 18.9540
#> 5  360 23.3680 26.5740 23.9320 24.4780 24.4950
#> 6  432 28.4320 32.3390 30.2910 29.8310 29.0020

# Load Gamma-Arrhenius dataset (Table 4)
data(gamma_arrhenius)
head(gamma_arrhenius)
#>   time   unit1   unit2   unit3   unit4   unit5
#> 1   72  4.2125  4.3152  2.6938  6.3355  4.2028
#> 2  144  9.3001  8.3011  6.4768 10.6550  9.7827
#> 3  216 13.6430 13.3310 10.9490 14.5920 15.0170
#> 4  288 18.5890 20.5590 16.2550 19.7330 19.2050
#> 5  360 24.8850 24.4380 21.7880 25.3580 24.2830
#> 6  432 30.3840 27.0450 27.3590 28.7250 32.2570

2. Fitting Step-Stress Accelerated Degradation Models

2.1. Wiener Process Model with Arrhenius Acceleration

We fit a 3-step step-stress Wiener-Arrhenius model to wiener_arrhenius data using Maximum Likelihood Estimation (MLE):

fit_w_mle <- ssad_fit(
  data = wiener_arrhenius,
  process = "wiener",
  model = "arrhenius",
  stress_levels = c(45, 65, 85),
  thresholds = c(90, 160),
  method = "mle"
)

summary(fit_w_mle)
#> 
#> =======================================================
#>   Detailed Summary: SSADT Model Fit
#> =======================================================
#> Process       : wiener 
#> Acceleration  : arrhenius 
#> Method        : mle 
#> 
#> Model Selection Criteria:
#>   Log-Likelihood : -224.643 
#>   AIC            : 455.287 
#>   BIC            : 464.319 
#> 
#> Maximum Likelihood Estimates:
#>            Estimate    Std.Error
#> a      5.257242e+00  0.161984667
#> b     -2.511549e+03 54.228973360
#> sigma  9.596252e-02  0.005448289
#> =======================================================

2.2. Gamma Process Model with Arrhenius Acceleration

Next, we fit a 3-step step-stress Gamma-Arrhenius model to gamma_arrhenius data using Bayesian MCMC:

fit_g_mcmc <- ssad_fit(
  data = gamma_arrhenius,
  process = "gamma",
  model = "arrhenius",
  stress_levels = c(45, 65, 85),
  thresholds = c(90, 160),
  method = "mcmc",
  n_iter = 1000,
  burnin = 200
)

summary(fit_g_mcmc)
#> 
#> =======================================================
#>   Detailed Summary: SSADT Model Fit
#> =======================================================
#> Process       : gamma 
#> Acceleration  : arrhenius 
#> Method        : mcmc 
#> 
#> Model Selection Criteria:
#>   Log-Likelihood : -297.185 
#>   AIC            : 600.37 
#>   BIC            : 609.401 
#> 
#> Bayesian MCMC Posterior Statistics:
#>          Mean    Std.Dev    MC.Error        X2.5.         X50.       X97.5.
#> a     6.90579  0.2127810 0.007522943     6.535857     6.924744     7.352241
#> b -2562.00038 61.0113208 2.157075933 -2703.025851 -2548.418621 -2498.378078
#> u     4.47113  0.3561228 0.012590841     3.862860     4.530658     4.934119
#> 
#> MCMC Acceptance Rate: 3.2 %
#> 
#> Geweke's MCMC Convergence Diagnostics:
#>   parameter   mean_first    mean_last    z_score      p_value converged
#> 1         a     6.934648     7.001323 -0.8457088 3.977152e-01      TRUE
#> 2         b -2599.741597 -2572.848200 -1.1057978 2.688141e-01      TRUE
#> 3         u     4.092808     4.738739 -6.0642863 1.325409e-09     FALSE
#> =======================================================

3. MCMC Convergence Diagnostics

The package includes Geweke’s MCMC convergence diagnostic (geweke_diag()):

if (!is.null(fit_g_mcmc$chain)) {
  g_res <- geweke_diag(fit_g_mcmc$chain)
  print(g_res)
}
#>   parameter   mean_first    mean_last    z_score      p_value converged
#> 1         a     6.934648     7.001323 -0.8457088 3.977152e-01      TRUE
#> 2         b -2599.741597 -2572.848200 -1.1057978 2.688141e-01      TRUE
#> 3         u     4.092808     4.738739 -6.0642863 1.325409e-09     FALSE

4. Lifetime Prediction and Reliability Analysis

Using the fitted parameters, we can predict the reliability function \(R(t)\) and Mean Time To Failure (MTTF) under normal operating stress conditions (\(S_0 = 25^\circ\text{C}\)) and failure threshold \(\omega_F = 200\):

rel_pred <- predict(
  fit_w_mle,
  t = seq(100, 10000, by = 200),
  S0 = 25,
  omega_F = 200
)

cat("Predicted MTTF under use stress S0 = 25:", round(rel_pred$mttf, 2), "hours\n")
#> Predicted MTTF under use stress S0 = 25: 4765.18 hours

5. Simulating Step-Stress Degradation Data

We can generate simulated step-stress degradation datasets for power calculation and experimental planning using ssad_simulate():

set.seed(123)
sim_df <- ssad_simulate(
  n_units = 5,
  times = seq(72, 2160, by = 72),
  stress_levels = c(45, 65, 85),
  thresholds = c(90, 160),
  process = "wiener",
  model = "power",
  params = c(a = 7.39e-4, b = 1.2, sigma = 0.1)
)

head(sim_df)
#>   time     unit1     unit2     unit3     unit4     unit5
#> 1   72  4.651050  5.488497  5.448765  5.969646  5.226456
#> 2  144  9.582368 10.364750 10.149159 11.561606  9.549127
#> 3  216 16.031606 16.250919 14.993053 16.890806 14.259505
#> 4  288 21.218064 22.122670 19.255392 21.484640 19.168833
#> 5  360 26.454398 27.946434 23.472577 27.765821 25.860032
#> 6  432 33.036308 33.657394 28.856759 32.383114 30.433464

References