Package 'survstan'

Description

The aim of the R package survstan is to provide a toolkit for fitting survival models using Stan. The R package survstan can be used to fit right-censored survival data under independent censoring. The implemented models allow the fitting of survival data in the presence/absence of covariates. All inferential procedures are currently based on the maximum likelihood (ML) approach.

_PACKAGE

References

Stan Development Team (2023). “RStan: the R interface to Stan.” R package version 2.21.8, https://mc-stan.org/.

Lawless JF (2002). Statistical Models and Methods for Lifetime Data, Wiley Series in Probability and Statistics, 2nd Edition edition. John Wiley and Sons. ISBN 9780471372158.

Bennett S (1983). “Analysis of survival data by the proportional odds model.” Statistics in Medicine, 2(2), 273-277. doi:10.1002/sim.4780020223.

Chen YQ, Wang M (2000). “Analysis of Accelerated Hazards Models.” Journal of the American Statistical Association, 95(450), 608-618. doi:10.1080/01621459.2000.10474236.

Demarqui FN, Mayrink VD (2021). “Yang and Prentice model with piecewise exponential baseline distribution for modeling lifetime data with crossing survival curves.” Brazilian Journal of Probability and Statistics, 35(1), 172 – 186. doi:10.1214/20-BJPS471.

Description

Function to fit accelerated failure time (AFT) models.

Usage

aftreg(formula, data, baseline = "weibull", dist = NULL, init = 0, ...)

Arguments

Value

aftreg returns an object of class "aftreg" containing the fitted model.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
summary(fit)

Description

Function to fit accelerated hazard (AH) models.

Usage

ahreg(formula, data, baseline = "weibull", dist = NULL, init = 0, ...)

Arguments

Value

ahreg returns an object of class "ahreg" containing the fitted model.

Examples

library(survstan)
fit <- ahreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
summary(fit)

Description

Akaike information criterion

Usage

## S3 method for class 'survstan'
AIC(object, ..., k = 2)

Arguments

Value

the Akaike information criterion value when a single model is passed to the function; otherwise, a data.frame with the Akaike information criterion values and the number of parameters is returned.

Examples

library(survstan)
fit1 <- aftreg(Surv(futime, fustat) ~ 1, data = ovarian, baseline = "weibull", init = 0)
fit2 <- aftreg(Surv(futime, fustat) ~ rx, data = ovarian, baseline = "weibull", init = 0)
fit3 <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
AIC(fit1, fit2, fit3)

Description

Compute analysis of variance (or deviance) tables for one or more fitted model objects.

Usage

## S3 method for class 'survstan'
anova(...)

Arguments

Value

the ANOVA table.

Examples

library(survstan)
fit1 <- aftreg(Surv(futime, fustat) ~ 1, data = ovarian, baseline = "weibull", init = 0)
fit2 <- aftreg(Surv(futime, fustat) ~ rx, data = ovarian, baseline = "weibull", init = 0)
fit3 <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
anova(fit1, fit2, fit3)

Description

This function is used to allow the user to specify an arbitrary value for the polynomial's degree m. If m = NULL, then m = min(m_max, ceiling(n^0.4)) is used, where m_max = 15.

Usage

bernstein(m = NULL)

Arguments

Value

a list with the baseline name and the polynomial's degree m.

Description

Estimated regression coefficients

Usage

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

Arguments

Value

the estimated regression coefficients

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
coef(fit)

Description

Confidence intervals for the regression coefficients

Usage

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

Arguments

Value

100(1-alpha) confidence intervals for the regression coefficients.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
confint(fit)

Description

Generic S3 method cross_time

Usage

cross_time(object, ...)

Arguments

Value

the crossing survival time

Description

Computes the crossing survival times

Usage

## S3 method for class 'survstan'
cross_time(
  object,
  newdata1,
  newdata2,
  conf.level = 0.95,
  nboot = 1000,
  cores = 1,
  ...
)

Arguments

Value

the crossing survival time

Examples

library(survstan)
data(ipass)
fit <- ypreg(Surv(time, status)~arm, data=ipass, baseline = "weibull")
summary(fit)
newdata1 <- data.frame(arm=0)
newdata2 <- data.frame(arm=1)
tcross <- cross_time(fit, newdata1, newdata2, nboot = 10)
tcross

Description

Function to fit Extended Hazard (EH) models.

Usage

ehreg(formula, data, baseline = "weibull", dist = NULL, init = 0, ...)

Arguments

Value

ehreg returns an object of class "ehreg" containing the fitted model.

Examples

library(survstan)
fit <- ehreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
summary(fit)

Description

Functions required for compatibility of survstan with emmeans. Users are not required to call these functions themselves. Instead, they will be called automatically by the emmeans function of the emmeans package.

Usage

recover_data.survstan(object, ...)

recover_data.ypreg(object, term = c("short", "long"), ...)

recover_data.ehreg(object, term = c("AF", "RH"), ...)

Arguments

Description

Parameters estimates of a survstan model

Usage

estimates(object, ...)

Arguments

Value

the parameters estimates of a given survstan model.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
estimates(fit)

Description

Computes the (generalized) Akaike An Information Criterion for a fitted parametric model.

Usage

## S3 method for class 'survstan'
extractAIC(fit, scale, k = 2, ...)

Arguments

Value

the ANOVA table.

Examples

library(survstan)
fit1 <- aftreg(Surv(futime, fustat) ~ 1, data = ovarian, baseline = "weibull", init = 0)
fit2 <- aftreg(Surv(futime, fustat) ~ rx, data = ovarian, baseline = "weibull", init = 0)
fit3 <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
extractAIC(fit1)
extractAIC(fit2)
extractAIC(fit3)

Description

Data set from a clinical trial conducted by the Gastrointestinal Tumor Study Group (GTSG) in 1982. The data set refers to the survival times of patients with locally nonresectable gastric cancer. Patients were either treated with chemotherapy combined with radiation or chemotherapy alone.

Format

A data frame with 90 rows and 3 variables:

• time: survival times (in days)

• status: failure indicator (1 - failure; 0 - otherwise)

• trt: treatments (1 - chemotherapy + radiation; 0 - chemotherapy alone)

Author(s)

Fabio N. Demarqui [email protected]

References

Gastrointestinal Tumor Study Group. (1982) A Comparison of Combination Chemotherapy and Combined Modality Therapy for Locally Advanced Gastric Carcinoma. Cancer 49:1771-7.

Description

Probability function, distribution function, quantile function and random generation for the distribution with parameters mu, sigma and varphi.

Usage

dggprentice(x, mu, sigma, varphi, log = FALSE)

pggprentice(q, mu = 0, sigma = 1, varphi, lower.tail = TRUE, log.p = FALSE)

qggprentice(p, mu = 0, sigma = 1, varphi, lower.tail = TRUE, log.p = FALSE)

rggprentice(n, mu = 0, sigma = 1, varphi, ...)

Arguments

Details

Probability density function:

$f(x | \mu, \sigma, \varphi) = \begin{cases} \frac{|\varphi|(\varphi^{-2})^{\varphi^{-2}}}{\sigma x\Gamma(\varphi^{-2})}\exp\{\varphi^{-2}[\varphi w - \exp(\varphi w)]\}I_{[0, \infty)}(x), & \varphi \neq 0 \\ \frac{1}{\sqrt{2\pi}x\sigma}\exp\left\{-\frac{1}{2}\left(\frac{log(x)-\mu}{\sigma}\right)^2\right\}I_{[0, \infty)}(x), & \varphi = 0 \end{cases}$

where $w = \frac{\log(x) - \mu}{\sigma}$, for $-\infty < \mu < \infty$, $\sigma>0$ and $-\infty < \varphi < \infty$.

Distribution function:

$F(x|\mu, \sigma, \varphi) = \begin{cases} F_{G}(y|1/\varphi^2, 1), & \varphi > 0 \\ 1-F_{G}(y|1/\varphi^2, 1), & \varphi < 0 \\ F_{LN}(x|\mu, \sigma), & \varphi = 0 \end{cases}$

where $y = \displaystyle\left(\frac{x}{\sigma}\right)^\varphi$, $F_{G}(\cdot|\nu, 1)$ is the distribution function of a gamma distribution with shape parameter $1/\varphi^2$ and scale parameter equals to 1, and $F_{LN}(x|\mu, \sigma)$ corresponds to the distribution function of a lognormal distribution with location parameter $\mu$ and scale parameter $\sigma$.

Value

dggprentice gives the (log) probability function, pggprentice gives the (log) distribution function, qggprentice gives the quantile function, and rggprentice generates random deviates.

Description

Generic S3 method ggresiduals

Usage

ggresiduals(object, ...)

Arguments

Details

Generic method to plot residuals of survival models.

Value

the desired residual plot.

Description

ggresiduals method for survstan models

Usage

## S3 method for class 'survstan'
ggresiduals(object, type = c("coxsnell", "martingale", "deviance"), ...)

Arguments

Details

This function produces residuals plots of Cox-Snell residuals, martingale residuals and deviance residuals.

Value

the desired residual plot.

Examples

library(survstan)
ovarian$rx <- as.factor(ovarian$rx)
fit <- aftreg(Surv(futime, fustat) ~ age + rx, data = ovarian, baseline = "weibull", init = 0)
ggresiduals(fit, type = "coxsnell")
ggresiduals(fit, type = "martingale")
ggresiduals(fit, type = "deviance")

Description

Probability function, distribution function, quantile function and random generation for the distribution with parameters alpha, gamma and kappa.

Usage

dggstacy(x, alpha, gamma, kappa, log = FALSE)

pggstacy(q, alpha, gamma, kappa, log.p = FALSE, lower.tail = TRUE)

qggstacy(
  p,
  alpha = 1,
  gamma = 1,
  kappa = 1,
  log.p = FALSE,
  lower.tail = TRUE,
  ...
)

rggstacy(n, alpha = 1, gamma = 1, kappa = 1, ...)

Arguments

Details

Probability density function:

$f(x|\alpha, \gamma, \kappa) = \frac{\kappa}{\gamma^{\alpha}\Gamma(\alpha/\kappa)}x^{\alpha-1}\exp\left\{-\left(\frac{x}{\gamma}\right)^{\kappa}\right\}I_{[0, \infty)}(x),$

for $\alpha>0$, $\gamma>0$ and $\kappa>0$.

Distribution function:

$F(t|\alpha, \gamma, \kappa) = F_{G}(x|\nu, 1),$

where $x = \displaystyle\left(\frac{t}{\gamma}\right)^\kappa$, and $F_{G}(\cdot|\nu, 1)$ corresponds to the distribution function of a gamma distribution with shape parameter $\nu = \alpha/\gamma$ and scale parameter equals to 1.

Value

dggstacy gives the (log) probability function, pggstacy gives the (log) distribution function, qggstacy gives the quantile function, and rggstacy generates random deviates.

Description

Probability function, distribution function, quantile function and random generation for the distribution with parameters alpha and gamma.

Usage

dgompertz(x, alpha = 1, gamma = 1, log = FALSE, ...)

pgompertz(q, alpha = 1, gamma = 1, lower.tail = TRUE, log.p = FALSE, ...)

qgompertz(p, alpha = 1, gamma = 1, lower.tail = FALSE, log.p = FALSE, ...)

rgompertz(n, alpha = 1, gamma = 1, ...)

Arguments

Details

Probability density function:

$f(x|\alpha, \gamma) = \alpha\gamma \exp\{\gamma x - \alpha(e^{\gamma x} - 1)\}I_{[0, \infty)}(x),$

for $\alpha>0$ and $\gamma>0$.

Distribution function:

$F(x|\alpha, \gamma) = 1 - \exp\{- \alpha(e^{\gamma x} - 1)\},$

for $x>0$, $\alpha>0$ and $\gamma>0$.

Value

dgompertz gives the (log) probability function, pgompertz gives the (log) distribution function, qgompertz gives the quantile function, and rgompertz generates random deviates.

Description

Reconstructed IPASS clinical trial data reported in Argyropoulos and Unruh (2015). Although reconstructed, this data set preserves all features exhibited in references with full access to the observations from this clinical trial. The data base is related to the period of March 2006 to April 2008. The main purpose of the study is to compare the drug gefitinib against carboplatin/paclitaxel doublet chemotherapy as first line treatment, in terms of progression free survival (in months), to be applied to selected non-small-cell lung cancer (NSCLC) patients.

Format

A data frame with 1217 rows and 3 variables:

• time: progression free survival (in months)

• status: failure indicator (1 - failure; 0 - otherwise)

• arm: (1 - gefitinib; 0 - carboplatin/paclitaxel doublet chemotherapy)

Author(s)

Fabio N. Demarqui [email protected]

References

Argyropoulos, C. and Unruh, M. L. (2015). Analysis of time to event outcomes in randomized controlled trials by generalized additive models. PLOS One 10, 1-33.

Description

Extracts the log-likelihood function for a fitted parametric model.

Usage

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

Arguments

Value

the log-likelihood value when a single model is passed to the function; otherwise, a data.frame with the log-likelihood values and the number of parameters is returned.

Examples

library(survstan)
fit1 <- aftreg(Surv(futime, fustat) ~ 1, data = ovarian, baseline = "weibull", init = 0)
fit2 <- aftreg(Surv(futime, fustat) ~ rx, data = ovarian, baseline = "weibull", init = 0)
fit3 <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
logLik(fit1, fit2, fit3)

Description

Reconstruct the model matrix for a survstan model.

Usage

## S3 method for class 'survstan'
model.matrix(object, ...)

Arguments

Value

The model matrix (or matrices) for the fit.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
model.matrix(fit)

Description

Hazard and cumulative hazard functions of the PE distribution

Usage

hpexp(x, rho, rates)

Hpexp(x, rho, rates)

Arguments

Value

hpexp gives the hazard function and Hpexp gives the cumulative hazard function of the PE distribution.

Description

Probability function, distribution function, quantile function and random generation for the Piecewise Exponential (PE) distribution.

Usage

dpexp(x, rho, rates, log = FALSE)

ppexp(q, rho, rates, lower.tail = TRUE, log.p = FALSE)

qpexp(p, rho, rates, lower.tail = TRUE, log.p = FALSE)

rpexp(n, rho, rates)

Arguments

Value

dpexp gives the (log) probability function, ppexp gives the (log) distribution function, qpexp gives the quantile function, and rpexp generates random deviates.

Examples

n <- 10
rho <- c(0, 1, 3, 7, Inf)
rates <- c(0.5, 4, 0.8, 0.1)
x <- sort(rpexp(n, rho=rho, rates=rates))
Fx <- ppexp(x, rho, rates)
y <- qpexp(Fx, rho, rates)
# checking:
x==y

Description

Function to fit proportional hazards (PH) models.

Usage

phreg(formula, data, baseline = "weibull", dist = NULL, init = 0, ...)

Arguments

Value

phreg returns an object of class "phreg" containing the fitted model.

Examples

library(survstan)
fit <- phreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
summary(fit)

Description

This function is used to allow the user to specify the piecewise exponential baseline with arbitrary time grid/number of intervals.

Usage

piecewise(rho = NULL, m = NULL)

Arguments

Value

a list with the baseline name, and the provided time grid and number of intervals.

Description

Function to fit proportional odds (PO) models.

Usage

poreg(formula, data, baseline = "weibull", dist = NULL, init = 0, ...)

Arguments

Value

poreg returns an object of class "poreg" containing the fitted model.

Examples

library(survstan)
fit <- poreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
summary(fit)

Description

Produces a printed summary of a fitted survstan model.

Usage

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

Arguments

Value

No return value, called for side effects.

Description

Rank a collection of survstan models

Usage

rank_models(formula, data, survreg, baseline, dist = NULL, ...)

Arguments

Value

a tibble containing the fitted models ranked according to their AICs.

Examples

library(survstan)
library(dplyr)

veteran <- veteran %>%
  mutate(across(c(trt, prior, celltype), as.factor))
fits <- rank_models(
  formula = Surv(time, status) ~ celltype+karno,
  data = veteran,
  survreg = c("aftreg", "ahreg", "phreg", "poreg", "ypreg", "ehreg"),
  baseline = c("exponential", "weibull", "lognormal", "loglogistic")
)

Description

residuals method for survstan models

Usage

## S3 method for class 'survstan'
residuals(object, type = c("coxsnell", "martingale", "deviance"), ...)

Arguments

Details

This function extracts the residuals, martingale residuals and deviance residuals of a survstan object.

Value

a vector containing the desired residuals.

Examples

library(survstan)
ovarian$rx <- as.factor(ovarian$rx)
fit <- aftreg(Surv(futime, fustat) ~ age + rx, data = ovarian, baseline = "weibull", init = 0)
residuals(fit, type = "coxsnell")
residuals(fit, type = "martingale")
residuals(fit, type = "deviance")

Description

Generic S3 method se

Usage

se(object, ...)

Arguments

Value

the standard errors associated with a set of parameter estimators for a given model.

Description

Estimated standard errors

Usage

## S3 method for class 'survstan'
se(object, ...)

Arguments

Value

a vector with the standard errors.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
se(fit)

Description

Summary for a survstan object

Usage

## S3 method for class 'survstan'
summary(object, conf.level = 0.95, ...)

Arguments

Value

an object of the class summary.survstan containing a summary of the fitted model.

Description

Computes the predicted survivor function for a phpe model.

Usage

## S3 method for class 'survstan'
survfit(formula, newdata = NULL, ...)

Arguments

Value

a data.frame containing the estimated survival probabilities.

Examples

library(survstan)
library(ggplot2)
data(ipass)
ipass$arm <- as.factor(ipass$arm)
fit <- ypreg(Surv(time, status)~arm, data=ipass, baseline = "weibull")
summary(fit)
newdata <- data.frame(arm=as.factor(0:1))
surv <- survfit(fit, newdata)
ggplot(surv, aes(x=time, y=surv, color = arm)) +
  geom_line()

Description

Tidy a survstan object

Usage

## S3 method for class 'survstan'
tidy(x, conf.int = FALSE, conf.level = 0.95, ...)

Arguments

Details

Convert a fitted model into a tibble.

Value

a tibble with a summary of the fit.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
tidy(fit)

Description

Time grid

Usage

time_grid(time, event, m = NULL)

Arguments

Value

Time grid.

Description

This function extracts and returns the variance-covariance matrix associated with the regression coefficients when the maximum likelihood estimation approach is used in the model fitting.

Usage

## S3 method for class 'survstan'
vcov(object, all = FALSE, ...)

Arguments

Value

the variance-covariance matrix associated with the parameters estimators.

Examples

library(survstan)
fit <- aftreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull", init = 0)
vcov(fit)

Description

Function to fit Yang and Prentice (YP) models.

Usage

ypreg(formula, data, baseline = "weibull", dist = NULL, init = 0, ...)

Arguments

Value

ypreg returns an object of class "ypreg" containing the fitted model.

Examples

library(survstan)
fit <- ypreg(Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian, baseline = "weibull")
summary(fit)