Package {MultiFrailty}

Title: Shared Frailty Regression Models with Inverse Gaussian, Generalized Lindley, and Gamma Frailty Distributions
Version: 0.1.0
Description: Implements shared frailty regression models for survival data under eight censoring mechanisms: exact, right censoring (Kalbfleisch and Prentice, 2002), left censoring, interval censoring (Sun, 2006), progressive Type I censoring, and progressive Type II censoring (Balakrishnan and Aggarwala, 2000 <doi:10.1007/978-1-4612-1334-5>). Combines four frailty distributions – Gamma (Clayton, 1978), Inverse Gaussian (Hougaard, 1984), and two variants of the Generalized Lindley (GL) distribution: GL Type 1, a two-component gamma mixture with distribution-specific scale/shape linkage (Pandey, Hanagal, and Tyagi, 2022), and GL Type 2, a two-component gamma mixture with a common rate parameter (Pandey and Tyagi, 2021 <doi:10.1134/S1995080222010140>) – with two baseline hazard distributions: the two-parameter Weibull distribution (Weibull, 1951) and the three-parameter Generalized (Exponentiated) Weibull distribution (Mudholkar and Srivastava, 1993 <doi:10.1109/24.229504>). A no-frailty baseline-only model is also supported for nested model comparison. Maximum likelihood estimation is conducted using Newton-Raphson and Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithms via the 'maxLik' package (Henningsen and Toomet, 2011 <doi:10.1007/s00180-010-0217-1>). Provides standard errors, confidence intervals, hypothesis tests, Akaike Information Criterion (AIC, Akaike, 1974 <doi:10.1109/TAC.1974.1100705>), Bayesian Information Criterion (BIC, Schwarz, 1978 <doi:10.1214/aos/1176344136>), corrected Akaike Information Criterion (AICc, Hurvich and Tsai, 1989), Hannan-Quinn Information Criterion (HQIC, Hannan and Quinn, 1979), a bootstrap approximation of the Widely Applicable Information Criterion (WAIC, Watanabe, 2010), k-fold cross-validation, frailty variance estimation, survival, hazard, median, risk, and marginal predictions, Cox-Snell (Cox and Snell, 1968), martingale (Barlow and Prentice, 1988), and deviance residuals with a Kolmogorov-Smirnov goodness-of-fit test, influence diagnostics (leverage, Cook's distance, difference in fits (DFFITS), difference in betas (DFBETAS); Belsley, Kuh, and Welsch, 1980), random data generation under all eight censoring mechanisms, a Monte Carlo simulation-study function, and a diagnostic and survival plotting suite.
License: GPL-3
Depends: R (≥ 4.0.0)
Imports: survival, maxLik, numDeriv, stats, graphics, grDevices, utils
Suggests: testthat (≥ 3.0.0), knitr, rmarkdown
VignetteBuilder: knitr
Encoding: UTF-8
Language: en-US
RoxygenNote: 7.3.3
Config/testthat/edition: 3
NeedsCompilation: no
Packaged: 2026-07-30 23:00:46 UTC; shikhar tyagi
Author: Shikhar Tyagi ORCID iD [aut, cre], Arvind Pandey [aut], Bhupendra Singh [aut], Vrijesh Tripathi [aut]
Maintainer: Shikhar Tyagi <shikhar1093tyagi@gmail.com>
Repository: CRAN
Date/Publication: 2026-08-07 19:40:06 UTC

MultiFrailty: Shared Frailty Regression Models with Inverse Gaussian, Generalized Lindley, and Gamma Frailty Distributions

Description

Implements shared frailty regression models for survival data under eight censoring mechanisms: exact, right censoring (Kalbfleisch and Prentice, 2002), left censoring, interval censoring (Sun, 2006), progressive Type I censoring, and progressive Type II censoring (Balakrishnan and Aggarwala, 2000 doi:10.1007/978-1-4612-1334-5). Combines four frailty distributions – Gamma (Clayton, 1978), Inverse Gaussian (Hougaard, 1984), and two variants of the Generalized Lindley (GL) distribution: GL Type 1, a two-component gamma mixture with distribution-specific scale/shape linkage (Pandey, Hanagal, and Tyagi, 2022), and GL Type 2, a two-component gamma mixture with a common rate parameter (Pandey and Tyagi, 2021 doi:10.1134/S1995080222010140) – with two baseline hazard distributions: the two-parameter Weibull distribution (Weibull, 1951) and the three-parameter Generalized (Exponentiated) Weibull distribution (Mudholkar and Srivastava, 1993 doi:10.1109/24.229504). A no-frailty baseline-only model is also supported for nested model comparison. Maximum likelihood estimation is conducted using Newton-Raphson and Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithms via the 'maxLik' package (Henningsen and Toomet, 2011 doi:10.1007/s00180-010-0217-1). Provides standard errors, confidence intervals, hypothesis tests, Akaike Information Criterion (AIC, Akaike, 1974 doi:10.1109/TAC.1974.1100705), Bayesian Information Criterion (BIC, Schwarz, 1978 doi:10.1214/aos/1176344136), corrected Akaike Information Criterion (AICc, Hurvich and Tsai, 1989), Hannan-Quinn Information Criterion (HQIC, Hannan and Quinn, 1979), a bootstrap approximation of the Widely Applicable Information Criterion (WAIC, Watanabe, 2010), k-fold cross-validation, frailty variance estimation, survival, hazard, median, risk, and marginal predictions, Cox-Snell (Cox and Snell, 1968), martingale (Barlow and Prentice, 1988), and deviance residuals with a Kolmogorov-Smirnov goodness-of-fit test, influence diagnostics (leverage, Cook's distance, difference in fits (DFFITS), difference in betas (DFBETAS); Belsley, Kuh, and Welsch, 1980), random data generation under all eight censoring mechanisms, a Monte Carlo simulation-study function, and a diagnostic and survival plotting suite.

Author(s)

Maintainer: Shikhar Tyagi shikhar1093tyagi@gmail.com (ORCID)

Authors:

Baseline Hazard and Cumulative Hazard Evaluation

Description

Computes baseline hazard phi0(t) and cumulative baseline hazard Phi0(t) for Weibull and Generalized Weibull (GW) distributions.

Usage

baseline_hazard(t, baseline = c("weibull", "gw"), par)

Arguments

t

Positive numeric vector of time points.

baseline

Character string specifying baseline distribution: "weibull" or "gw".

par

Numeric vector of parameters: c(lambda, gamma) for Weibull, or c(delta, zeta, xi) for GW.

Value

A named list with components:

H0

Cumulative baseline hazard Phi0(t).

h0

Baseline hazard phi0(t).

References

Mudholkar, G. S., & Srivastava, D. K. (1993). Exponentiated Weibull family for analyzing bathtub failure-rate data. IEEE Transactions on Reliability, 42(2), 299-302.

Examples

bh_weib <- baseline_hazard(1:5, baseline = "weibull", par = c(2, 1.5))
bh_gw <- baseline_hazard(1:5, baseline = "gw", par = c(0.5, 1.2, 1.1))

Bootstrapped WAIC Computation

Description

Bootstrapped WAIC Computation

Usage

bootstrap_waic(fit, B = 200)

Arguments

fit

Fitted multifrailty_fit object.

B

Number of bootstrap replications. Default is 200.

Value

Numeric value of bootstrapped WAIC.

Model Comparison, Bootstrapped WAIC, and Cross-Validation for MultiFrailty Models

Description

Compares multiple fitted shared frailty models using AIC, BIC, AICc, HQIC, optional bootstrapped WAIC, and k-fold cross-validation.

Usage

compare_models(
  ...,
  criteria = c("AIC", "BIC", "AICc", "HQIC"),
  compute_waic = FALSE,
  waic_B = 100L
)

Arguments

...

Fitted "multifrailty_fit" objects passed individually or as a list.

criteria

Character vector of selection criteria to include.

compute_waic

Logical; if TRUE computes bootstrapped WAIC. Default is FALSE.

waic_B

Number of bootstrap samples for WAIC. Default is 100.

Value

A data frame summarizing model comparisons ranked by AIC.

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
dat <- r_frailty(n = 60, baseline = "weibull", bpar = c(2, 1.5), frailty = "gamma", fpar = c(0.8))
fit1 <- fit_frailty(time = dat$time, status = dat$status, baseline = "weibull", frailty = "none")
fit2 <- fit_frailty(time = dat$time, status = dat$status, baseline = "weibull", frailty = "gamma")
cmp <- compare_models(fit1, fit2)
print(cmp)

K-Fold Cross Validation for MultiFrailty Models

Description

K-Fold Cross Validation for MultiFrailty Models

Usage

cv_frailty(
  time,
  status,
  x = matrix(nrow = length(time), ncol = 0),
  baseline = "weibull",
  frailty = "gamma",
  k = 5L,
  time2 = NULL
)

Arguments

time

Time vector.

status

Status vector.

x

Covariate matrix.

baseline

Baseline distribution.

frailty

Frailty distribution.

k

Number of folds. Default is 5L.

time2

Optional upper interval bound vector.

Value

Out-of-sample total log-likelihood score.

Diagnostic Summary Table

Description

Diagnostic Summary Table

Usage

diagnostics_table(fit)

Arguments

fit

Fitted multifrailty_fit object.

Value

Data frame of influence diagnostics.

Low-Level Maximum Likelihood Estimator for MultiFrailty Models

Description

Fits shared frailty regression models using Maximum Likelihood Estimation (MLE) across all 10 baseline-frailty combinations with robust optimizer fallbacks.

Usage

fit_frailty(
  time,
  status,
  x = matrix(nrow = length(time), ncol = 0),
  baseline = c("weibull", "gw"),
  frailty = c("none", "gamma", "ig", "gl1", "gl2"),
  time2 = NULL,
  prog_cen = NULL,
  init = NULL,
  method = "NR",
  ...
)

Arguments

time

Primary survival/censoring time vector.

status

Event indicator vector (0 = right-censored, 1 = event, 2 = left-censored, 3 = interval-censored).

x

Design matrix of covariates (n x p). Default is a 0-column matrix.

baseline

Character string for baseline hazard: "weibull" or "gw".

frailty

Character string for frailty family: "none", "gamma", "ig", "gl1", or "gl2".

time2

Vector of upper interval bounds when status == 3. Default is NULL.

prog_cen

Vector of progressive censoring counts. Default is NULL.

init

Vector of initial values on the estimation scale. Default is NULL (automatic).

method

Optimization method passed to maxLik: "NR" (Newton-Raphson) or "BFGS".

...

Additional arguments passed to optimization algorithms.

Value

An object of class "multifrailty_fit" containing parameter estimates, standard errors, information criteria, variance-covariance matrix, and diagnostic statistics.

References

Hougaard, P. (1984). Life table methods for heterogeneous populations: distributions of frailties. Biometrika, 71(1), 75-83.

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
dat <- r_frailty(n = 100, baseline = "weibull", bpar = c(2, 1.5),
                 frailty = "gamma", fpar = c(0.8),
                 x = matrix(rnorm(100), ncol = 1), beta = 0.5)
fit <- fit_frailty(time = dat$time, status = dat$status, x = as.matrix(dat[, "X1", drop=FALSE]),
                   baseline = "weibull", frailty = "gamma")
print(fit)

Future Survival Forecast

Description

Future Survival Forecast

Usage

forecast_frailty(fit, horizon, n_grid = 200, newdata = NULL)

Arguments

fit

Fitted multifrailty_fit object.

horizon

Forecast horizon time.

n_grid

Number of evaluation points. Default 200.

newdata

Optional new data frame.

Value

Matrix of forecasted survival probabilities.

Unified Survival, Density, Hazard, and Cumulative Hazard Combiner

Description

Evaluates unconditional survival function S(t), probability density function f(t), hazard function h(t), and cumulative hazard function H(t) for any specified combination of frailty distribution and baseline hazard.

Usage

frailty_functions(
  t,
  eta = 1,
  frailty = c("none", "gamma", "ig", "gl1", "gl2"),
  fpar = numeric(0),
  baseline = c("weibull", "gw"),
  bpar
)

Arguments

t

Positive numeric vector of time points.

eta

Linear predictor risk scores rho = exp(X %*% beta) (numeric vector of same length as t or scalar).

frailty

Character string specifying frailty family: "none", "gamma", "ig", "gl1", or "gl2".

fpar

Numeric vector of parameters for frailty distribution.

baseline

Character string specifying baseline hazard: "weibull" or "gw".

bpar

Numeric vector of parameters for baseline hazard.

Value

A named list containing:

S

Unconditional survival probabilities S(t).

f

Unconditional density values f(t).

h

Unconditional hazard values h(t).

H

Unconditional cumulative hazard values H(t).

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

gf <- frailty_functions(1:5, eta = 1, frailty = "gl1", fpar = c(2.13, 1.53),
                        baseline = "gw", bpar = c(0.71, 2.13, 1.54))

Laplace Transform and Derivatives for Frailty Distributions

Description

Computes Laplace transform L(s), first derivative L'(s), and frailty variance Var(W) for five frailty families: "none", "gamma", "ig", "gl1", and "gl2".

Usage

frailty_laplace(
  s,
  frailty = c("none", "gamma", "ig", "gl1", "gl2"),
  par = numeric(0)
)

Arguments

s

Non-negative numeric vector.

frailty

Character string specifying frailty family: "none", "gamma", "ig", "gl1", or "gl2".

par

Numeric vector of parameters per family: numeric(0) for "none", c(theta) for "gamma", c(eta) for "ig", c(eta, epsilon) for "gl1", or c(theta, mu) for "gl2".

Value

A named list with components:

L

Laplace transform L(s).

L1

First derivative L'(s).

Var

Variance of frailty distribution Var(W), evaluated as L”(0) - 1.

References

Hougaard, P. (1984). Life table methods for heterogeneous populations: distributions of frailties. Biometrika, 71(1), 75-83.

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

fl_gamma <- frailty_laplace(s = 0.5, frailty = "gamma", par = c(0.8))
fl_gl1   <- frailty_laplace(s = 0.5, frailty = "gl1", par = c(1.2, 0.5))
fl_gl2   <- frailty_laplace(s = 0.5, frailty = "gl2", par = c(1.5, 0.8))

Influence Diagnostics for MultiFrailty Models

Description

Computes leverage (hat values), Cook's distance, DFBETAS, and DFFITS for regression parameters.

Usage

influence_frailty(fit)

Arguments

fit

A fitted object of class "multifrailty_fit".

Value

A list containing leverage, Cook's distance, DFBETAS matrix, DFFITS vector, and a summary diagnostic table.

References

Belsley, D. A., Kuh, E., & Welsch, R. E. (1980). Regression diagnostics: Identifying influential data and sources of collinearity. John Wiley & Sons.

Examples

set.seed(123)
dat <- r_frailty(n = 50, baseline = "weibull", bpar = c(2, 1.5), frailty = "gamma", fpar = c(0.8))
fit <- fit_frailty(time = dat$time, status = dat$status, baseline = "weibull", frailty = "gamma")
inf <- influence_frailty(fit)
head(inf$diagnostics_table)

Numerically stable log(1 - exp(-x))

Description

Accurately evaluates log(1 - exp(-x)) for any x >= 0 avoiding catastrophic loss of precision for small or large x.

Usage

log1mexp(x)

Arguments

x

Numeric vector of non-negative values.

Value

Numeric vector of log(1 - exp(-x)).

Numerically stable log(1 + exp(x))

Description

Accurately evaluates log(1 + exp(x)) for all real x.

Usage

log1pexp(x)

Arguments

x

Numeric vector.

Value

Numeric vector of log(1 + exp(x)).

Log-sum-exp helper for two log-scale components

Description

Computes log(exp(a) + exp(b)) in log-space.

Usage

log_sum_exp(a, b)

Arguments

a

Numeric vector of log-scale component 1.

b

Numeric vector of log-scale component 2.

Value

Numeric vector of log(exp(a) + exp(b)).

Log-Likelihood Function for MultiFrailty Models

Description

Computes the log-likelihood for shared frailty models across all 10 baseline-frailty combinations with support for right, exact, left, and interval censoring, plus optional progressive censoring.

Usage

loglik_frailty(
  par_all,
  time,
  status,
  x = matrix(nrow = length(time), ncol = 0),
  baseline = c("weibull", "gw"),
  frailty = c("none", "gamma", "ig", "gl1", "gl2"),
  time2 = NULL,
  prog_cen = NULL
)

Arguments

par_all

Vector of all model parameters on estimation scale (log/logit transformed).

time

Primary event/censoring time vector.

status

Event status vector (0 = right-censored, 1 = exact event, 2 = left-censored, 3 = interval-censored).

x

Matrix of covariates (n x p). Default is 0-column matrix.

baseline

Baseline hazard distribution ("weibull" or "gw").

frailty

Frailty distribution ("none", "gamma", "ig", "gl1", or "gl2").

time2

Vector of upper interval bounds when status == 3. Default is NULL.

prog_cen

Vector of progressive censoring counts R_i per observation. Default is NULL.

Value

Scalar log-likelihood value. Returns -1e12 sentinel on numerical invalidity.

Examples

par_all <- c(log(2), log(1.5), log(0.8), 0.1) # Weibull + Gamma + 1 beta
time <- c(1, 2, 3, 4)
status <- c(1, 0, 1, 0)
x <- matrix(c(0.5, -0.2, 0.1, 0.8), ncol = 1)
ll <- loglik_frailty(par_all, time, status, x, baseline = "weibull", frailty = "gamma")

Formula Interface for Shared Frailty Regression Models

Description

Fits shared frailty regression models using a formula specifying survival response (e.g., survival::Surv(time, status) ~ x1 + x2) and data frame.

Usage

multifrailty(
  formula,
  data,
  baseline = c("weibull", "gw"),
  frailty = c("none", "gamma", "ig", "gl1", "gl2"),
  method = "NR",
  ...
)

Arguments

formula

Model formula of the form survival::Surv(time, status) ~ x1 + x2.

data

Data frame containing variables in formula.

baseline

Character string specifying baseline hazard: "weibull" or "gw".

frailty

Character string specifying frailty family: "none", "gamma", "ig", "gl1", or "gl2".

method

Optimization method for maxLik: "NR" (Newton-Raphson) or "BFGS".

...

Additional arguments passed to fit_frailty.

Value

Object of class c("multifrailty", "multifrailty_fit").

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
df <- data.frame(
  time = stats::rexp(50, rate = 0.1),
  status = sample(c(0, 1), 50, replace = TRUE),
  age = stats::rnorm(50, mean = 50, sd = 10),
  sex = sample(c(0, 1), 50, replace = TRUE)
)
fit <- multifrailty(survival::Surv(time, status) ~ age + sex, data = df,
                    baseline = "weibull", frailty = "gamma")
print(fit)

Full Diagnostic Plot Suite

Description

Displays comprehensive multi-panel diagnostic plot suite.

Usage

plot_all(fit, ask = grDevices::dev.interactive(), save_to_file = NULL)

Arguments

fit

Fitted model.

ask

Logical; if TRUE prompts user.

save_to_file

Optional file path to save plot grid.

Value

No return value, called for side effects (displays a 2x3 grid of diagnostic plots or saves them to a file).

Baseline Hazard Plot

Description

Baseline Hazard Plot

Usage

plot_baseline(
  baseline = c("weibull", "gw"),
  par,
  t_range = c(0.01, 3),
  n_grid = 300
)

Arguments

baseline

Baseline distribution.

par

Baseline parameters.

t_range

Time range vector.

n_grid

Grid resolution.

Value

No return value, called for side effects (displays baseline hazard and cumulative baseline hazard plots).

Coefficient Forest Plot

Description

Coefficient Forest Plot

Usage

plot_coef_forest(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays a forest plot of parameter estimates with 95% confidence intervals).

DFBETAS Plot

Description

DFBETAS Plot

Usage

plot_dfbetas(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays boxplots of DFBETAS values across parameters).

Frailty Density Overlay Plot across Families

Description

Frailty Density Overlay Plot across Families

Usage

plot_frailty_density(t_max = 3)

Arguments

t_max

Maximum grid value.

Value

No return value, called for side effects (displays frailty density overlay curves across distributions).

Leverage Index Plot

Description

Leverage Index Plot

Usage

plot_leverage(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays an index plot of observation leverage values).

Q-Q Plot of Cox-Snell Residuals vs Exp(1)

Description

Q-Q Plot of Cox-Snell Residuals vs Exp(1)

Usage

plot_qq_residuals(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays a Q-Q plot of Cox-Snell residuals versus standard exponential distribution).

Residuals vs Fitted Plot

Description

Residuals vs Fitted Plot

Usage

plot_residuals_fitted(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays a plot of martingale residuals against fitted survival probabilities).

Residuals vs Leverage Plot

Description

Residuals vs Leverage Plot

Usage

plot_residuals_leverage(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays a plot of standardized residuals against leverage values).

Scale-Location Plot

Description

Scale-Location Plot

Usage

plot_scale_location(fit)

Arguments

fit

Fitted model.

Value

No return value, called for side effects (displays a scale-location plot of square-root absolute standardized residuals against fitted survival probabilities).

Diagnostic and Visualization Suite for MultiFrailty Models

Description

Provides specialized plotting capabilities for fitted shared frailty models, including frailty density overlays across all five families, baseline hazard curves, Kaplan-Meier survival curves, residual diagnostic plots, coefficient forest plots, and influence diagnostic plots.

Usage

plot_survival_km(fit)

Arguments

fit

A fitted object of class "multifrailty_fit".

Value

No return value, called for side effects (displays a Kaplan-Meier survival curve with fitted MultiFrailty model overlay).

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
dat <- r_frailty(n = 60, baseline = "weibull", bpar = c(2, 1.5), frailty = "gamma", fpar = c(0.8))
fit <- fit_frailty(time = dat$time, status = dat$status, baseline = "weibull", frailty = "gamma")
plot_survival_km(fit)

Predictions for MultiFrailty Regression Models

Description

Computes survival probabilities, hazard rates, median survival times, expected survival times, risk scores, marginal survival curves, or future survival forecasts for a fitted multifrailty model.

Usage

predict_frailty(
  fit,
  newdata = NULL,
  newtime = NULL,
  type = c("survival", "hazard", "median", "expected", "risk", "marginal", "forecast"),
  window = NULL
)

Arguments

fit

A fitted object of class "multifrailty_fit".

newdata

Optional data frame of new covariate values. If NULL, uses training data.

newtime

Optional vector of evaluation time points. If NULL, uses default grid.

type

Type of prediction: "survival", "hazard", "median", "expected", "risk", "marginal", or "forecast".

window

Optional forecast window or horizon parameter.

Value

Vector or matrix of predictions depending on type.

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
dat <- r_frailty(n = 60, baseline = "weibull", bpar = c(2, 1.5), frailty = "gamma", fpar = c(0.8))
fit <- fit_frailty(time = dat$time, status = dat$status, baseline = "weibull", frailty = "gamma")
pred_surv <- predict_frailty(fit, type = "survival", newtime = c(1, 2, 3))

S3 Methods for MultiFrailty Fit Objects

Description

Provides standard print, summary, coef, vcov, logLik, AIC, BIC, and confint S3 methods for objects of class "multifrailty_fit".

Usage

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

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

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

## S3 method for class 'multifrailty_fit'
coef(object, ...)

## S3 method for class 'multifrailty_fit'
vcov(object, ...)

## S3 method for class 'multifrailty_fit'
logLik(object, ...)

## S3 method for class 'multifrailty_fit'
AIC(object, ...)

## S3 method for class 'multifrailty_fit'
BIC(object, ...)

## S3 method for class 'multifrailty_fit'
confint(object, parm, level = 0.95, ...)

Arguments

x

An object of class "multifrailty_fit".

...

Additional arguments.

object

An object of class "multifrailty_fit".

parm

Optional parameter vector or names.

level

Confidence level for interval. Default 0.95.

Value

Formatted console output or numeric objects depending on method.

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Random Data Generation under Eight Censoring Schemes

Description

Generates survival times and censoring indicators for shared frailty regression models across all 10 baseline-frailty combinations under eight distinct censoring mechanisms.

Usage

r_frailty(
  n,
  baseline = c("weibull", "gw"),
  bpar,
  frailty = c("none", "gamma", "ig", "gl1", "gl2"),
  fpar = numeric(0),
  x = matrix(nrow = n, ncol = 0),
  beta = numeric(0),
  cen_type = c("none", "right", "left", "interval", "type1", "type2", "progressive",
    "progressive_type1"),
  cen_rate = 0.2,
  left_threshold = NULL,
  int_width = NULL,
  cen_time = NULL,
  r_failures = NULL,
  prog_scheme = NULL,
  prog_times = NULL
)

Arguments

n

Number of observations to generate.

baseline

Baseline hazard distribution ("weibull" or "gw").

bpar

Baseline parameter vector.

frailty

Frailty distribution ("none", "gamma", "ig", "gl1", or "gl2").

fpar

Frailty parameter vector.

x

Matrix of covariates (n x p). Default is 0-column matrix.

beta

Regression coefficient vector matching columns of x.

cen_type

Censoring mechanism: "none", "right", "left", "interval", "type1", "type2", "progressive", or "progressive_type1".

cen_rate

Exponential rate for right-censoring time generation. Default is 0.2.

left_threshold

Threshold for left censoring. Default is 20th percentile.

int_width

Width of censoring window for interval censoring. Default is 20% of mean time.

cen_time

Fixed cutoff time for Type-I censoring. Default is 70th percentile.

r_failures

Target number of failures for Type-II censoring. Default is floor(0.7 * n).

prog_scheme

Vector of progressive removal counts for progressive censoring.

prog_times

Inspection time points for progressive Type-I censoring.

Value

A data frame containing generated time, time2 (for interval), status, and covariates.

References

Hougaard, P. (1984). Life table methods for heterogeneous populations: distributions of frailties. Biometrika, 71(1), 75-83.

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
dat <- r_frailty(n = 100, baseline = "weibull", bpar = c(2, 1.5),
                 frailty = "gl1", fpar = c(1.2, 0.5),
                 cen_type = "right", cen_rate = 0.1)
head(dat)

Random Generation for Generalized Lindley Type 1 (GL1) Frailty Distribution

Description

Generates random variates from GL Type 1 frailty distribution (two-component Gamma mixture).

Usage

r_gl1(n, eta, epsilon)

Arguments

n

Number of observations to generate.

eta

Parameter eta > 0.

epsilon

Parameter epsilon > 0.

Value

Numeric vector of length n.

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Examples

set.seed(123)
w_gl1 <- r_gl1(100, eta = 1.2, epsilon = 0.5)

Random Generation for Generalized Lindley Type 2 (GL2) Frailty Distribution

Description

Generates random variates from GL Type 2 frailty distribution (two-component Gamma mixture with common rate).

Usage

r_gl2(n, theta, mu)

Arguments

n

Number of observations to generate.

theta

Parameter theta > 0.

mu

Parameter mu in (0, 1 + theta).

Value

Numeric vector of length n.

References

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
w_gl2 <- r_gl2(100, theta = 1.5, mu = 0.8)

Random Generation for Generalized Weibull (GW) Baseline Distribution

Description

Generates random variates from the 3-parameter Generalized Weibull baseline distribution.

Usage

r_gw(n, delta, zeta, xi)

Arguments

n

Number of observations to generate.

delta

Scale parameter (delta > 0).

zeta

Shape parameter (zeta > 0).

xi

Shape parameter (xi > 0).

Value

Numeric vector of length n containing random samples.

References

Mudholkar, G. S., & Srivastava, D. K. (1993). Exponentiated Weibull family for analyzing bathtub failure-rate data. IEEE Transactions on Reliability, 42(2), 299-302.

Examples

set.seed(123)
sim_data <- r_gw(100, delta = 0.5, zeta = 1.2, xi = 1.1)

Random Generation for Inverse Gaussian (IG) Frailty Distribution

Description

Generates random variates from the Inverse Gaussian frailty distribution with mean E[W] = 1.

Usage

r_ig(n, eta)

Arguments

n

Number of observations to generate.

eta

Frailty variance parameter (eta > 0).

Value

Numeric vector of length n.

References

Hougaard, P. (1984). Life table methods for heterogeneous populations: distributions of frailties. Biometrika, 71(1), 75-83.

Examples

set.seed(123)
w_ig <- r_ig(100, eta = 0.5)

Residual Diagnostics and Goodness-of-Fit Tests for MultiFrailty Models

Description

Computes Cox-Snell, Martingale, Deviance, raw, standardized, and studentized residuals along with Kolmogorov-Smirnov (K-S) goodness-of-fit statistics against standard Exponential(1).

Usage

residuals_frailty(fit)

Arguments

fit

A fitted object of class "multifrailty_fit".

Value

A list containing residual vectors, summary accuracy metrics (MSE, RMSE, MAE, R_square, Adj_R_square), and K-S test results (KS_stat, KS_pvalue).

References

Cox, D. R., & Snell, E. J. (1968). A general definition of residuals. Journal of the Royal Statistical Society: Series B (Methodological), 30(2), 248-265.

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples

set.seed(123)
dat <- r_frailty(n = 60, baseline = "weibull", bpar = c(2, 1.5), frailty = "gamma", fpar = c(0.8))
fit <- fit_frailty(time = dat$time, status = dat$status, baseline = "weibull", frailty = "gamma")
res <- residuals_frailty(fit)
print(res$KS_pvalue)

Risk Score Predictions

Description

Risk Score Predictions

Usage

risk_predict(fit, newdata, times = NULL)

Arguments

fit

Fitted multifrailty_fit object.

newdata

Data frame of new observations.

times

Optional time grid.

Value

Numeric vector of risk scores exp(X %*% beta).

Safe exponential function with overflow capping

Description

Safe exponential function with overflow capping

Usage

safe_exp(x)

Arguments

x

Numeric vector.

Value

Numeric vector capped to prevent overflow/underflow.

Safe power evaluation

Description

Safe power evaluation

Usage

safe_pow(base, exponent)

Arguments

base

Positive numeric base.

exponent

Numeric exponent.

Value

Numeric result of base^exponent.

Monte Carlo Simulation Framework for MultiFrailty Models

Description

Evaluates frequentist Maximum Likelihood Estimation performance (bias, relative bias, MSE, empirical coverage) across repeated Monte Carlo replications under user-specified baseline, frailty, and censoring schemes.

Usage

simulate_mle_performance(
  n_sim = 50,
  n = 100,
  baseline = "weibull",
  bpar = c(2, 1.5),
  frailty = "gamma",
  fpar = c(0.8),
  beta = 0.5,
  cen_type = "right",
  cen_rate = 0.2
)

Arguments

n_sim

Number of Monte Carlo simulation replicates. Default is 50.

n

Sample size per replicate. Default is 100.

baseline

Character string for baseline hazard ("weibull" or "gw").

bpar

True baseline parameter vector.

frailty

Character string for frailty family ("none", "gamma", "ig", "gl1", or "gl2").

fpar

True frailty parameter vector.

beta

True regression parameter. Default 0.5.

cen_type

Censoring mechanism. Default "right".

cen_rate

Censoring rate. Default 0.2.

Value

A data frame summarizing parameter true values, mean estimates, bias, relative bias, MSE, and coverage.

References

Pandey, A., Hanagal, D. D., & Tyagi, S. (2022). Shared Frailty Models Based on Cancer Data. International Journal of Statistics and Reliability Engineering, 9(3), 461-474.

Pandey, A., & Tyagi, S. (2021). Comparison of Multiplicative Frailty Models Under Weibull Baseline Distribution. Lobachevskii Journal of Mathematics, 42(13), 3184-3195.

Examples


set.seed(123)
sim_res <- simulate_mle_performance(n_sim = 10, n = 50, baseline = "weibull",
                                    bpar = c(2, 1.5), frailty = "gamma", fpar = c(0.8))
print(sim_res)

Run Comprehensive Simulation Study

Description

Run Comprehensive Simulation Study

Usage

simulation_study(
  sample_sizes = c(25, 50, 100),
  n_sim = 20,
  baseline = "weibull",
  frailty = "gl2"
)

Arguments

sample_sizes

Vector of sample sizes. Default c(25, 50, 100).

n_sim

Number of simulation replicates per setting. Default 20.

baseline

Baseline hazard. Default "weibull".

frailty

Frailty distribution. Default "gl2".

Value

Master data frame of Monte Carlo simulation metrics across sample sizes.

Survival Probability at Specific Time Points

Description

Survival Probability at Specific Time Points

Usage

survival_at(fit, times, newdata = NULL)

Arguments

fit

Fitted multifrailty_fit object.

times

Vector of target time points.

newdata

Optional new data frame.

Value

Matrix of survival probabilities.