Package 'GammaFrailtySPC'

The Weibull Accelerated Failure Time Model with Gamma Frailty

Description

Density, distribution function, quantile function, and random generation for the Weibull Accelerated Failure Time (AFT) distribution integrated with gamma frailty as defined in Asadzadeh (2022).

Usage

dgamma_aft(
  x,
  beta0,
  beta1,
  covariates = 0,
  shape_kappa,
  frailty_lambda,
  log = FALSE
)

pgamma_aft(
  q,
  beta0,
  beta1,
  covariates = 0,
  shape_kappa,
  frailty_lambda,
  lower.tail = TRUE,
  log.p = FALSE
)

qgamma_aft(
  p,
  beta0,
  beta1,
  covariates = 0,
  shape_kappa,
  frailty_lambda,
  lower.tail = TRUE,
  log.p = FALSE
)

rgamma_aft(n, beta0, beta1, covariates = 0, shape_kappa, frailty_lambda)

Arguments

x	Vector of non-negative quantiles or random generation sample size (rgamma_aft).
beta0	Intercept parameter of the AFT log-scale model ($\beta_0$).
beta1	Vector or scalar of regression coefficients ($\beta_1$).
covariates	Matrix or vector of observed covariate values matching beta1.
shape_kappa	Positive shape parameter of the Weibull baseline distribution ($\kappa > 0$).
frailty_lambda	Positive parameter of the gamma frailty distribution ($\lambda > 0$).
log, log.p	Logical; if TRUE, probabilities/densities $p$ are given as $\log(p)$.
q	Vector of quantiles.
lower.tail	Logical; if TRUE (default), probabilities are $P[Y \le y]$, otherwise $P[Y > y]$.
p	Vector of probabilities.
n	Number of observations for random generation.

Value

dgamma_aft gives the density, pgamma_aft gives the distribution function, qgamma_aft gives the quantile function, and rgamma_aft generates random deviates.

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

# Density calculation
dgamma_aft(x = 10, beta0 = 4.12, beta1 = -0.15, covariates = 3,
           shape_kappa = 1.72, frailty_lambda = 1.17)

# Cumulative distribution calculation
pgamma_aft(q = 10, beta0 = 4.12, beta1 = -0.15, covariates = 3,
           shape_kappa = 1.72, frailty_lambda = 1.17)

# Quantile calculation
qgamma_aft(p = 0.5, beta0 = 4.12, beta1 = -0.15, covariates = 3,
           shape_kappa = 1.72, frailty_lambda = 1.17)

# Random sample generation
set.seed(123)
rgamma_aft(n = 5, beta0 = 4.12, beta1 = -0.15, covariates = rep(3, 5),
           shape_kappa = 1.72, frailty_lambda = 1.17)

---

Average Run Length (ARL) Simulation Study

Description

Evaluates the detectability and performance of the PrL, EWMA, and CUSUM control schemes under various mean shift sizes ($\upsilon'$) and censoring rates via Monte Carlo simulation as detailed in Section 4 of Asadzadeh (2022).

Usage

spc_arl_sim(
  n_sim = 100,
  max_rl = 1000,
  beta0 = 4.12,
  beta1 = -0.15,
  shape_kappa = 1.72,
  frailty_lambda = 1.17,
  shifts = c(1, 0.975, 0.95, 0.9, 0.8, 0.7),
  censor_rate = 0,
  chart_type = c("cusum", "ewma", "prl"),
  lcl = -2.5,
  omega = 0.05,
  alpha = 0.005,
  seed = NULL
)

Arguments

n_sim	Number of simulation replicates (default 100 for fast checks, 10000 for paper-exact precision).
max_rl	Maximum run length cap to prevent infinite loops (default 1000).
beta0	Intercept parameter $\beta_0$ (default 4.12).
beta1	Regression coefficient $\beta_1$ (default -0.15).
shape_kappa	Weibull shape parameter $\kappa$ (default 1.72).
frailty_lambda	Gamma frailty parameter $\lambda$ (default 1.17).
shifts	Vector of scale parameter shift multipliers $\upsilon'$ (default c(1.0, 0.975, 0.95, 0.90, 0.80, 0.70)).
censor_rate	Target right-censoring rate (0 for uncensored, 0.2, 0.5, 0.8).
chart_type	Control chart type to evaluate: "cusum", "ewma", or "prl".
lcl	Lower control limit for the specified chart.
omega	Smoothing parameter for EWMA chart (default 0.05).
alpha	False alarm probability for PrL chart (default 0.005).
seed	Random seed for reproducibility.

Value

A data frame containing shift sizes ('upsilon_prime'), calculated 'ARL', and standard errors ('SE').

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

# Fast demonstration with 10 simulation replicates
res_cusum <- spc_arl_sim(n_sim = 10, shifts = c(1.0, 0.95),
                         chart_type = "cusum", lcl = -2.5, seed = 123)
print(res_cusum)

---

CUSUM Control Chart for Heterogeneous Reliability Data

Description

Computes the likelihood-ratio based one-sided Cumulative Sum (CUSUM) control chart statistics for monitoring downward mean shifts under the Weibull AFT gamma frailty model as derived in Asadzadeh (2022). Also supports competing CUSUM schemes: "iuc" (ignoring unobserved frailty) and "iuoc" (ignoring both unobserved frailty and observed covariates).

Usage

spc_cusum(
  y,
  x,
  delta = NULL,
  fit = NULL,
  beta0 = NULL,
  beta1 = NULL,
  shape_kappa = NULL,
  frailty_lambda = NULL,
  upsilon = 0.95,
  lcl = -2.5,
  type = c("proposed", "iuc", "iuoc")
)

Arguments

y	Vector of observed lifetimes (or censoring thresholds).
x	Matrix or vector of observed covariate values.
delta	Vector of event indicators (1 for observed event/uncensored, 0 for right-censored). If NULL, all observations are assumed uncensored.
fit	Optional spc_fit object from spc_gamma_frailty_fit.
beta0	Intercept parameter $\beta_0$. Required if fit is NULL.
beta1	Vector or scalar of regression coefficients $\beta_1$. Required if fit is NULL.
shape_kappa	Weibull shape parameter $\kappa > 0$. Required if fit is NULL.
frailty_lambda	Gamma frailty parameter $\lambda > 0$. Required if fit is NULL.
upsilon	Shift magnitude factor $\upsilon < 1$ (e.g. 0.95 for a 5 percent downward shift).
lcl	Lower control limit threshold (default -2.5).
type	Type of CUSUM chart to compute: "proposed" (AFT-Frailty CUSUM), "iuc" (ignores unobserved covariate), or "iuoc" (ignores unobserved and observed covariates).

Value

An S3 object of class 'spc_cusum' containing: - 't_stat': CUSUM control statistics $t_i$. - 'scores': Computed score sequence $w_i$. - 'lcl': Lower control limit threshold. - 'signals': Logical vector indicating out-of-control signals ('t_stat < lcl'). - 'signal_indices': Indices of out-of-control observations. - 'upsilon': Shift magnitude factor $\upsilon$. - 'type': Type of CUSUM chart computed.

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

set.seed(123)
x_val <- rnorm(30, mean = 3, sd = 1)
y_val <- rgamma_aft(n = 30, beta0 = 4.12, beta1 = -0.15,
                    covariates = x_val, shape_kappa = 1.72, frailty_lambda = 1.17)
cusum_res <- spc_cusum(y = y_val, x = x_val, beta0 = 4.12, beta1 = -0.15,
                       shape_kappa = 1.72, frailty_lambda = 1.17, upsilon = 0.95, lcl = -2.5)
print(cusum_res)

---

EWMA Control Chart with Conditional Expected Values (CEV)

Description

Computes exponentially weighted moving average (EWMA) control statistics with Conditional Expected Values (CEV) for right-censored observations under the Weibull AFT gamma frailty model as proposed by Asadzadeh (2022).

Usage

spc_ewma(
  y,
  x,
  delta = NULL,
  fit = NULL,
  beta0 = NULL,
  beta1 = NULL,
  shape_kappa = NULL,
  frailty_lambda = NULL,
  omega = 0.05,
  lcl = -3
)

Arguments

y	Vector of observed lifetimes (or censoring thresholds).
x	Matrix or vector of observed covariate values.
delta	Vector of event indicators (1 for observed event/uncensored, 0 for right-censored). If NULL, all observations are assumed uncensored.
fit	Optional spc_fit object from spc_gamma_frailty_fit.
beta0	Intercept parameter $\beta_0$. Required if fit is NULL.
beta1	Vector or scalar of regression coefficients $\beta_1$. Required if fit is NULL.
shape_kappa	Weibull shape parameter $\kappa > 0$. Required if fit is NULL.
frailty_lambda	Gamma frailty parameter $\lambda > 0$. Required if fit is NULL.
omega	Smoothing parameter $\omega \in (0, 1]$ (default 0.05).
lcl	Lower control limit threshold (default -3.0).

Value

An S3 object of class 'spc_ewma' containing: - 'q_stat': EWMA statistics $Q_i$. - 'cev': Conditional expected values for observations. - 'y_std': Standardized variables $y^*_i$. - 'lcl': Lower control limit threshold. - 'signals': Logical vector indicating out-of-control signals ('q_stat < lcl'). - 'signal_indices': Indices of out-of-control observations. - 'omega': Smoothing parameter $\omega$.

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

set.seed(123)
x_val <- rnorm(30, mean = 3, sd = 1)
y_val <- rgamma_aft(n = 30, beta0 = 4.12, beta1 = -0.15,
                    covariates = x_val, shape_kappa = 1.72, frailty_lambda = 1.17)
ewma_res <- spc_ewma(y = y_val, x = x_val, beta0 = 4.12, beta1 = -0.15,
                     shape_kappa = 1.72, frailty_lambda = 1.17, omega = 0.05, lcl = -2.5)
print(ewma_res)

---

Comprehensive Statistical Process Control for Heterogeneous Reliability Data

Description

Master function for performing statistical process control (SPC) analysis on heterogeneous reliability observations using Weibull AFT Gamma Frailty models. Fits Phase I model parameters (or accepts user-provided parameters) and constructs Phase II control charts: Probability Limits (PrL), EWMA with CEV, and likelihood-based CUSUM.

Usage

spc_gamma_frailty(
  y,
  x,
  delta = NULL,
  phase1_y = NULL,
  phase1_x = NULL,
  phase1_delta = NULL,
  fit = NULL,
  beta0 = NULL,
  beta1 = NULL,
  shape_kappa = NULL,
  frailty_lambda = NULL,
  alpha = 0.005,
  omega = 0.05,
  ewma_lcl = -3,
  upsilon = 0.95,
  cusum_lcl = -2.5
)

Arguments

y	Vector of Phase II response lifetimes.
x	Matrix or vector of Phase II observed covariates.
delta	Vector of Phase II event indicators (1 for observed event, 0 for right-censored). If NULL, all Phase II observations are assumed uncensored.
phase1_y	Optional Phase I historical lifetime vector for parameter estimation.
phase1_x	Optional Phase I covariate matrix/vector for parameter estimation.
phase1_delta	Optional Phase I event indicator vector for parameter estimation.
fit	Optional pre-fitted spc_fit object from spc_gamma_frailty_fit.
beta0	Intercept parameter $\beta_0$ if not using fit or Phase I data.
beta1	Regression coefficients $\beta_1$ if not using fit or Phase I data.
shape_kappa	Weibull shape parameter $\kappa$ if not using fit or Phase I data.
frailty_lambda	Gamma frailty parameter $\lambda$ if not using fit or Phase I data.
alpha	False alarm rate for PrL chart (default 0.005).
omega	Smoothing parameter for EWMA chart (default 0.05).
ewma_lcl	Lower control limit for EWMA chart (default -3.0).
upsilon	Shift magnitude factor for CUSUM chart (default 0.95).
cusum_lcl	Lower control limit for CUSUM chart (default -2.5).

Value

An S3 object of class 'spc_gf' containing: - 'fit': Phase I model fit object ('spc_fit'). - 'prl': PrL control chart results ('spc_prl'). - 'ewma': EWMA control chart results ('spc_ewma'). - 'cusum': CUSUM control chart results ('spc_cusum'). - 'summary_stats': Data frame of per-observation statistics and signals across all charts. - 'call': The matched function call.

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

# Generate synthetic Phase I and Phase II datasets
set.seed(456)
p1_x <- rnorm(100, mean = 3, sd = 1)
p1_y <- rgamma_aft(n = 100, beta0 = 4.12, beta1 = -0.15,
                   covariates = p1_x, shape_kappa = 1.72, frailty_lambda = 1.17)
p1_d <- ifelse(p1_y <= 60, 1, 0)
p1_y <- pmin(p1_y, 60)

p2_x <- rnorm(40, mean = 3, sd = 1)
p2_y <- rgamma_aft(n = 40, beta0 = 4.12, beta1 = -0.15,
                   covariates = p2_x, shape_kappa = 1.72, frailty_lambda = 1.17)
# Introduce a 10% downward shift after observation 20
p2_y[21:40] <- p2_y[21:40] * 0.90

res <- spc_gamma_frailty(y = p2_y, x = p2_x, phase1_y = p1_y,
                         phase1_x = p1_x, phase1_delta = p1_d)
print(res)
summary(res)
plot(res)

---

Fit Weibull AFT Gamma Frailty Model (Phase I Parameter Estimation)

Description

Estimates the parameters of the Weibull Accelerated Failure Time (AFT) model integrated with continuous gamma frailty for complete or right-censored Phase I dataset using Maximum Likelihood Estimation (MLE) as described in Asadzadeh (2022).

Usage

spc_gamma_frailty_fit(y, x, delta = NULL, init = NULL)

Arguments

y	Vector of observed lifetimes (or censoring times).
x	Matrix, data frame, or vector of observed covariates.
delta	Vector of event indicators (1 for observed event/uncensored, 0 for right-censored). If NULL, all observations are assumed uncensored.
init	Initial parameter values for optimization (beta0, beta1, shape_kappa, frailty_lambda). If NULL, smart starting values are derived automatically.

Value

An S3 object of class 'spc_fit' containing: - 'coefficients': Estimated parameters ($\beta_0$, $\beta_1$, shape_kappa, frailty_lambda). - 'se': Standard errors of parameter estimates. - 'vcov': Variance-covariance matrix of parameter estimates. - 'loglik': Maximum log-likelihood value. - 'aic': Akaike Information Criterion. - 'bic': Bayesian Information Criterion. - 'n_obs': Total number of observations. - 'n_events': Number of uncensored events. - 'call': The matched function call.

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

set.seed(123)
x_val <- rnorm(100, mean = 3, sd = 1)
y_val <- rgamma_aft(n = 100, beta0 = 4.12, beta1 = -0.15,
                    covariates = x_val, shape_kappa = 1.72, frailty_lambda = 1.17)
delta_val <- ifelse(y_val <= 60, 1, 0)
y_obs <- pmin(y_val, 60)

fit <- spc_gamma_frailty_fit(y = y_obs, x = x_val, delta = delta_val)
print(fit)
summary(fit)

---

Print and Summary Methods for Statistical Process Control Objects

Description

Standard S3 print, summary, and plot methods for fitted models and control chart objects generated by GammaFrailtySPC.

Usage

## S3 method for class 'spc_fit'
print(x, ...)

## S3 method for class 'spc_fit'
summary(object, ...)

## S3 method for class 'spc_prl'
print(x, ...)

## S3 method for class 'spc_ewma'
print(x, ...)

## S3 method for class 'spc_cusum'
print(x, ...)

## S3 method for class 'spc_gf'
print(x, ...)

## S3 method for class 'spc_gf'
summary(object, ...)

## S3 method for class 'spc_gf'
plot(x, ...)

Arguments

x	An object of class spc_fit, spc_prl, spc_ewma, spc_cusum, or spc_gf.
...	Additional arguments passed to low-level functions.
object	An object of class spc_fit or spc_gf.

Value

Invisibly returns the input object x or object.

---

Probability Limits Control Chart (PrL)

Description

Computes the probability-limits-based control chart (PrL) statistics and lower control limits (LCL) for monitoring reliability data under the Weibull AFT gamma frailty model as derived in Asadzadeh (2022).

Usage

spc_prl(
  y,
  x,
  fit = NULL,
  beta0 = NULL,
  beta1 = NULL,
  shape_kappa = NULL,
  frailty_lambda = NULL,
  alpha = 0.005
)

Arguments

y	Vector of Phase II response quality observations (lifetimes).
x	Matrix or vector of observed covariate values.
fit	Optional spc_fit object from spc_gamma_frailty_fit.
beta0	Intercept parameter $\beta_0$. Required if fit is NULL.
beta1	Vector or scalar of regression coefficients $\beta_1$. Required if fit is NULL.
shape_kappa	Weibull shape parameter $\kappa > 0$. Required if fit is NULL.
frailty_lambda	Gamma frailty parameter $\lambda > 0$. Required if fit is NULL.
alpha	Desired false alarm probability (default 0.005).

Value

An S3 object of class 'spc_prl' containing: - 'y': Observed lifetime values. - 'lcl': Computed lower control limit values $LCL_{x,i}$. - 'signals': Logical vector indicating out-of-control signals ('y < lcl'). - 'signal_indices': Indices of out-of-control observations. - 'alpha': Specified false alarm probability.

References

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

set.seed(123)
x_val <- rnorm(30, mean = 3, sd = 1)
y_val <- rgamma_aft(n = 30, beta0 = 4.12, beta1 = -0.15,
                    covariates = x_val, shape_kappa = 1.72, frailty_lambda = 1.17)
prl_res <- spc_prl(y = y_val, x = x_val, beta0 = 4.12, beta1 = -0.15,
                   shape_kappa = 1.72, frailty_lambda = 1.17, alpha = 0.005)
print(prl_res)

---

Textile Industry Reliability and Tensile Strength Dataset

Description

Data from a textile company in Kashan, Iran, analyzing thread tensile strength (response lifetime $Y$, in cN) as a function of the proportion of polyester to wool ratio (observed covariate $X$) under unobserved fiber maturity (gamma frailty), as presented in Section 5 of Asadzadeh (2022).

Usage

textile_data

Format

A data frame with 150 rows and 4 variables:

y

Thread tensile strength endurance lifetime (in cN).

x

Proportion of polyester to wool ratio (observed covariate).

delta

Censoring indicator (1 = uncensored event, 0 = right-censored observation).

phase

Phase indicator ("PhaseI" for the first 100 observations used for parameter estimation, and "PhaseII" for the subsequent 50 observations monitored for process shifts).

Source

Asadzadeh, S. (2022). Statistical process control based on gamma-frailty models for heterogeneous reliability observations. *Journal of Statistical Computation and Simulation*, 92(2), 337-351. doi:10.1080/00949655.2021.1959582

Examples

data(textile_data)
head(textile_data)
summary(textile_data)

---