Inverse Gaussian Process Degradation Models with Frailty"
In IGPFrailty: Inverse Gaussian Process Degradation Models with Frailty

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 7,
  fig.height = 5
)
library(IGPFrailty)

Introduction

In reliability engineering, components designed for high-reliability systems rarely fail during traditional life testing. Monitoring performance degradation over time provides rich statistical information to assess component reliability and predict remaining useful life without waiting for catastrophic failures.

The IGPFrailty package implements classical Inverse Gaussian Process (IGP) degradation models as well as advanced frailty extensions (Gamma frailty and Inverse Gaussian frailty) that capture unobserved unit-to-unit heterogeneity.

The methodology implemented in this package is based on:

Morita, L. H. M., Tomazella, V. L. D., Balakrishnan, N., Ramos, P. L., Ferreira, P. H., & Louzada, F. (2021). Inverse Gaussian process model with frailty term in reliability analysis. Quality and Reliability Engineering International, 37(2), 763–784. .

Application 1: GaAs Laser Degradation Data

The laser dataset consists of 15 Gallium Arsenide (GaAs) laser devices tested at 80 degrees Celsius over 4,000 hours with 16 equidistant inspection intervals. Failure is defined as a 10% increase in operating current.

data(laser)
head(laser)

Model Fitting and Selection

We fit the Classical IGP, IGP-Gamma frailty, and IGP-IG frailty models:

fit_none  <- igp_fit(laser, time_col = "t", deg_col = "increase", unit_col = "unit", frailty = "none")
fit_gamma <- igp_fit(laser, time_col = "t", deg_col = "increase", unit_col = "unit", frailty = "gamma")
fit_ig    <- igp_fit(laser, time_col = "t", deg_col = "increase", unit_col = "unit", frailty = "ig")

# Model comparison table
model_comp <- data.frame(
  Model = c("Classical IGP", "IGP-Gamma Frailty", "IGP-IG Frailty"),
  theta = c(coef(fit_none)["theta"], coef(fit_gamma)["theta"], coef(fit_ig)["theta"]),
  eta = c(coef(fit_none)["eta"], coef(fit_gamma)["eta"], coef(fit_ig)["eta"]),
  xi = c(NA, coef(fit_gamma)["xi"], coef(fit_ig)["xi"]),
  logLik = c(logLik(fit_none), logLik(fit_gamma), logLik(fit_ig)),
  AIC = c(AIC(fit_none), AIC(fit_gamma), AIC(fit_ig)),
  BIC = c(BIC(fit_none), BIC(fit_gamma), BIC(fit_ig))
)
model_comp[, -1] <- round(model_comp[, -1], 4)
knitr::kable(model_comp)

Posterior Individual Frailty Estimation

frail_gamma <- individual_frailty(fit_gamma)
print(frail_gamma)

Implied Lifetime Distribution and Quantiles

lt_gamma <- lifetime_dist(fit_gamma, threshold = 10, probs = c(0.01, 0.05, 0.1, 0.5, 0.8))
print(lt_gamma)

Diagnostic Visualizations

plot(fit_gamma, type = "all", threshold = 10)

Application 2: Fatigue Crack Size Growth Data

The crack dataset contains crack growth measurements for 21 specimens tested up to 0.12 million cycles with initial crack length 0.90 inches and failure threshold 1.60 inches (rho = 0.5754 under transformed metric log(D(t)/0.9)).

data(crack)
fit_crack_gam <- igp_fit(crack, time_col = "t", deg_col = "deg", unit_col = "specimen", frailty = "gamma")
summary(fit_crack_gam)

lt_crack <- lifetime_dist(fit_crack_gam, threshold = 0.5754, probs = c(0.01, 0.05, 0.1, 0.5, 0.8))
print(lt_crack)

Simulation Study

You can easily simulate degradation paths from any of the three models using sim_igp():

set.seed(42)
sim_paths <- sim_igp(n = 10, times = seq(0, 4, length.out = 11), theta = 2, eta = 15, xi = 0.2, frailty = "gamma")
head(sim_paths)

References

Morita, L. H. M., et al. (2021). Inverse Gaussian process model with frailty term in reliability analysis. Quality and Reliability Engineering International, 37(2), 763–784. .
Meeker, W. Q., & Escobar, L. A. (1998). Statistical Methods for Reliability Data. John Wiley & Sons. .
Lu, C. J., & Meeker, W. Q. (1993). Using degradation measures to estimate a time-to-failure distribution. Technometrics, 35(2), 161–174. .