flexsurv (version 2.2)

coxsnell_flexsurvreg: Cox-Snell residuals from a parametric survival model

Description

Cox-Snell residuals from a parametric survival model

Usage

coxsnell_flexsurvreg(x)

Value

A data frame with a column called est giving the Cox-Snell residual, defined as the fitted cumulative hazard at each data point. fitted cumulative hazard at the given observed data point, and other columns indicating the observation time and covariate values defining the data at this point.

An extra column "(qexp)" gives the equally-spaced quantiles of a standard exponential distribution in the same order as est. To check the fit of the model, "(qexp)" is plotted against est, and the points should form a straight line through the origin with slope 1.

Arguments

x

Object returned by flexsurvreg or flexsurvspline representing a fitted survival model

Examples

  fitg <- flexsurvreg(formula = Surv(futime, fustat) ~ age, data = ovarian, dist = "gengamma")
  cs <- coxsnell_flexsurvreg(fitg)
  
  ## Model doesn't appear to fit well since the cumulative hazards are underestimated.
  ## In this example, this is probably because the dataset is small, 
  ## hence the point estimate is noisy.
  plot(cs$"(qexp)", cs$est, pch=19, xlab="Theoretical quantiles", ylab="Cumulative hazard")
  abline(a=0,b=1,col="red",lwd=2)
  
  ## Alternative way to produce the same plot using "qqplot"
  qy <- qexp(ppoints(nrow(cs),0))
  qqplot(qy, cs$est)
  abline(a=0,b=1, col="red", lwd=2)
  
  ## A log transform may or may not bring out the pattern more clearly
  plot(log(cs$"(qexp)"), log(cs$est), pch=19)
  abline(a=0,b=1, col="red", lwd=2)
  
  ## In the model `fitg`, the fitted cumulative hazard is lower than the true cumulative hazard
  ## Another way to show this is to compare parametric vs nonparametric estimates of 
  ## the cumulative hazard 
  plot(fitg, type="cumhaz", ci=FALSE)
  
  ## Alternative example: where the true model is fitted to simulated data
  ## The model fits well
  y <- rweibull(10000, 2, 2)
  fite <- flexsurvreg(Surv(y) ~ 1, dist="weibull")
  cs <- coxsnell_flexsurvreg(fite)
  plot(cs$"(qexp)", cs$est, pch=19, xlab="Theoretical quantiles", ylab="Cumulative hazard")
  abline(a=0,b=1,col="red",lwd=2)