survplot: Plot Survival Curves and Hazard Functions

Description

Plot estimated survival curves, and for parametric survival models, plot hazard functions. There is an option to print the number of subjects at risk at the start of each time interval. Curves are automatically labeled at the points of maximum separation (using the labcurve function), and there are many other options for labeling that can be specified with the label.curves parameter. For example, different plotting symbols can be placed at constant x-increments and a legend linking the symbols with category labels can automatically positioned on the most empty portion of the plot.

For the case of a two stratum analysis by survfit, survdiffplot plots the difference in two Kaplan-Meier estimates along with approximate confidence bands for the differences, with a reference line at zero. The number of subjects at risk is optionally plotted. This number is taken as the minimum of the number of subjects at risk over the two strata.

Usage

survplot(fit, ...)
## S3 method for class 'rms':
survplot(fit, \dots, xlim,
         ylim=if(loglog) c(-5, 1.5) else if
                 (what == "survival" & missing(fun)) c(0, 1),
         xlab, ylab, time.inc,
         what=c("survival","hazard"),
         type=c("tsiatis","kaplan-meier"),
         conf.type=c("log","log-log","plain","none"),
         conf.int=FALSE, conf=c("bands","bars"),
         add=FALSE, label.curves=TRUE,
         abbrev.label=FALSE, lty, lwd=par("lwd"),
         col=1, col.fill=gray(seq(.95, .75, length=5)),
         adj.subtitle=TRUE, loglog=FALSE, fun,
         n.risk=FALSE, logt=FALSE, dots=FALSE, dotsize=.003,
         grid=FALSE, srt.n.risk=0, sep.n.risk=0.056, adj.n.risk=1, 
         y.n.risk, cex.n.risk=.6, pr=FALSE)
## S3 method for class 'survfit':
survplot(fit, xlim, 
         ylim, xlab, ylab, time.inc,
         conf=c("bands","bars","none"), add=FALSE, 
         label.curves=TRUE, abbrev.label=FALSE,
         lty,lwd=par('lwd'),
         col=1, col.fill=gray(seq(.95, .75, length=5)),
         loglog=FALSE,fun,n.risk=FALSE,logt=FALSE,
         dots=FALSE,dotsize=.003,
         grid=FALSE,
         srt.n.risk=0,sep.n.risk=.056,adj.n.risk=1,
         y.n.risk,cex.n.risk=.6, pr=FALSE, ...)
survdiffplot(fit, order=1:2,
           xlim, ylim, xlab, ylab="Difference in Survival Probability",
           time.inc, conf.int=.95,
           conf=c("shaded", "bands","none"),
           add=FALSE, lty=1, lwd=par('lwd'), col=1,
           n.risk=FALSE, grid=FALSE,
           srt.n.risk=0, adj.n.risk=1,
           y.n.risk, cex.n.risk=.6)

Arguments

fit

result of fit (cph, psm, survfit, survest.psm). For survdiffplot, fit must be the result of survfit.

...

list of factors with names used in model. For fits from survfit, these arguments do not appear - all strata are plotted. Otherwise the first factor listed is the factor used to determine different survival curves. Any other factors are used

xlim

a vector of two numbers specifiying the x-axis range for follow-up time. Default is (0,maxtime) where maxtime was the pretty()d version of the maximum follow-up time in any stratum, stored in fit$maxtime

ylim

y-axis limits. Default is c(0,1) for survival, and c(-5,1.5) if loglog=TRUE. If fun or loglog=TRUE are given and ylim is not, the limits will be computed from the data. For <

xlab

x-axis label. Default is units attribute of failure time variable given to Surv.

ylab

y-axis label. Default is "Survival Probability" or "log(-log Survival Probability)". If fun is given, the default is "". For what="hazard", the default is "Hazard Function".

time.inc

time increment for labeling the x-axis and printing numbers at risk. If not specified, the value of time.inc stored with the model fit will be used.

type

specifies type of estimates, "tsiatis" (the default) or "kaplan-meier". "tsiatis" here corresponds to the Breslow estimator. This is ignored if survival estimates stored with surv=TRUE are being used. Fo

conf.type

specifies the basis for confidence limits. This argument is ignored for fits from survfit.

conf.int

Default is FALSE. Specify e.g. .95 to plot 0.95 confidence bands. For fits from parametric survival models, or Cox models with x=TRUE and y=TRUE specified to the fit, the exact asymptotic formulas will

conf

"bars" for confidence bars at each time.inc time point. If the fit was from cph(..., surv=TRUE), the time.inc used will be that stored with the fit. Use conf="bands" (the default) for bands

what

defaults to "survival" to plot survival estimates. Set to "hazard" or an abbreviation to plot the hazard function (for psm fits only). Confidence intervals are not available for what="hazard".

add

set to TRUE to add curves to an existing plot.

label.curves

default is TRUE to use labcurve to label curves where they are farthest apart. Set label.curves to a list to specify options to labcurve, e.g., label.curves=list(method="arrow", cex=.

abbrev.label

set to TRUE to abbreviate() curve labels that are plotted

lty

vector of line types to use for different factor levels. Default is c(1,3,4,5,6,7,...).

lwd

vector of line widths to use for different factor levels. Default is current par setting for lwd.

col

color for curve, default is 1. Specify a vector to assign different colors to different curves.

col.fill

a vector of colors to used in filling confidence bands

adj.subtitle

set to FALSE to suppress plotting subtitle with levels of adjustment factors not plotted. Defaults to TRUE. This argument is ignored for survfit.

loglog

set to TRUE to plot log(-log Survival) instead of Survival

fun

specifies any function to translate estimates and confidence limits before plotting

logt

set to TRUE to plot log(t) instead of t on the x-axis

n.risk

set to TRUE to add number of subjects at risk for each curve, using the surv.summary created by cph or using the failure times used in fitting the model if y=TRUE was specified to the fit or if the fit w

srt.n.risk

angle of rotation for leftmost number of subjects at risk (since this number may run into the second or into the y-axis). Default is 0.

adj.n.risk

justification for leftmost number at risk. Default is 1 for right justification. Use 0 for left justification, .5 for centered.

sep.n.risk

multiple of upper y limit - lower y limit for separating lines of text containing number of subjects at risk. Default is .056*(ylim[2]-ylim[1]).

y.n.risk

When n.risk=TRUE, the default is to place numbers of patients at risk above the x-axis. You can specify a y-coordinate for the bottom line of the numbers using y.n.risk.

cex.n.risk

character size for number of subjects at risk (when n.risk is TRUE)

dots

set to TRUE to plot a grid of dots. Will be plotted at every time.inc (see cph) and at survival increments of .1 (if d>.4), .05 (if

.2 < d <= .4<="" code="">), or .025
(if d <= .2<="" code="">), where <

dotsize

size of dots in inches

grid

defaults to FALSE. Set to a color shading to plot faint lines. Set to 1 to plot solid lines. Default is .05 if TRUE.

set to TRUE to print survival curve coordinates used in the plots

order

an integer vector of length two specifying the order of groups when computing survival differences. The default of 1:2 indicates that the second group is subtracted from the first. Specify order=2:1 to instead subtract th

Value

list with components adjust (text string specifying adjustment levels) and curve.labels (vector of text strings corresponding to levels of factor used to distinguish curves). For survfit, the returned value is the vector of strata labels, or NULL if there are no strata.

Side Effects

plots. If

par()$mar[4]<4< code="">, issues par(mar=) to increment mar[4] by 2
if n.risk=TRUE and add=FALSE. The user may want to reset par(mar) in
this case to not leave such a wide right margin for plots. You usually
would issue par(mar=c(5,4,4,2)+.1).

Details

survplot will not work for Cox models with time-dependent covariables. Use survest or survfit for that purpose.

Use ps.slide, win.slide, gs.slide to set up nice defaults for plotting. These also set a system option mgp.axis.labels to allow x and y-axes to have differing mgp graphical parameters (see par). This is important when labels for y-axis tick marks are to be written horizontally (par(las=1)), as a larger gap between the labels and the tick marks are needed. You can set the axis-specific 2nd component of mgp using mgp.axis.labels(c(xvalue,yvalue)).

Examples

Run this code

# Simulate data from a population model in which the log hazard
# function is linear in age and there is no age x sex interaction
n <- 1000
set.seed(731)
age <- 50 + 12*rnorm(n)
label(age) <- "Age"
sex <- factor(sample(c('male','female'), n, TRUE))
cens <- 15*runif(n)
h <- .02*exp(.04*(age-50)+.8*(sex=='Female'))
dt <- -log(runif(n))/h
label(dt) <- 'Follow-up Time'
e <- ifelse(dt <= cens,1,0)
dt <- pmin(dt, cens)
units(dt) <- "Year"
dd <- datadist(age, sex)
options(datadist='dd')
S <- Surv(dt,e)

# When age is in the model by itself and we predict at the mean age,
# approximate confidence intervals are ok

f <- cph(S ~ age, surv=TRUE)
survplot(f, age=mean(age), conf.int=.95)
g <- cph(S ~ age, x=TRUE, y=TRUE)
survplot(g, age=mean(age), conf.int=.95, add=TRUE, col='red', conf='bars')

# Repeat for an age far from the mean; not ok
survplot(f, age=75, conf.int=.95)
survplot(g, age=75, conf.int=.95, add=TRUE, col='red', conf='bars')


#Plot stratified survival curves by sex, adj for quadratic age effect
# with age x sex interaction (2 d.f. interaction)

f <- cph(S ~ pol(age,2)*strat(sex), x=TRUE, y=TRUE)
#or f <- psm(S ~ pol(age,2)*sex)
Predict(f, sex=., age=c(30,50,70))
survplot(f, sex=., n.risk=TRUE)           #Adjust age to median
survplot(f, sex=., logt=TRUE, loglog=TRUE)   #Check for Weibull-ness (linearity)
survplot(f, sex=c("male","female"), age=50)
                                        #Would have worked without datadist
                                        #or with an incomplete datadist
survplot(f, sex=., label.curves=list(keys=c(2,0), point.inc=2))
                                        #Identify curves with symbols


survplot(f, sex=., label.curves=list(keys=c('m','f')))
                                        #Identify curves with single letters


#Plots by quintiles of age, adjusting sex to male
options(digits=3)
survplot(f, age=quantile(age,(1:4)/5), sex="male")


#Plot survival Kaplan-Meier survival estimates for males
f <- survfit(S ~ 1, subset=sex=="male")
survplot(f)


#Plot survival for both sexes
f <- survfit(S ~ sex)
survplot(f)
#Check for log-normal and log-logistic fits
survplot(f, fun=qnorm, ylab="Inverse Normal Transform")
survplot(f, fun=function(y)log(y/(1-y)), ylab="Logit S(t)")

#Plot the difference between sexes
survdiffplot(f)

options(datadist=NULL)

Run the code above in your browser using DataCamp Workspace