-
-
Save PhilipPallmann/86b5ee998d203097376b to your computer and use it in GitHub Desktop.
Draw ggplot2-style Kaplan-Meier curves (slightly modified version of mattbaggott/ggsurvival.R).
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
ggsurv <- function(f.survfit, f.CI="default", f.shape=3){ | |
require(ggplot2, quietly=TRUE) | |
# initialise frame variable | |
f.frame <- NULL | |
# check if more then one strata | |
if(length(names(f.survfit$strata)) == 0){ | |
# create data.frame with data from survfit | |
f.frame <- data.frame(time=f.survfit$time, n.risk=f.survfit$n.risk, n.event=f.survfit$n.event, | |
n.censor = f.survfit$n.censor, surv=f.survfit$surv, upper=f.survfit$upper, | |
lower=f.survfit$lower) | |
# create first two rows (start at 1) | |
f.start <- data.frame(time=c(0, f.frame$time[1]), n.risk=c(f.survfit$n, f.survfit$n), n.event=c(0,0), | |
n.censor=c(0,0), surv=c(1,1), upper=c(1,1), lower=c(1,1)) | |
# add first row to dataset | |
f.frame <- rbind(f.start, f.frame) | |
# remove temporary data | |
rm(f.start) | |
} | |
else { | |
# create vector for strata identification | |
f.strata <- NULL | |
for(f.i in 1:length(f.survfit$strata)){ | |
# add vector for one strata according to number of rows of strata | |
f.strata <- c(f.strata, rep(unlist(strsplit(names(f.survfit$strata)[f.i], "="))[2], f.survfit$strata[f.i])) | |
} | |
# create data.frame with data from survfit (create column for strata) | |
f.frame <- data.frame(time=f.survfit$time, n.risk=f.survfit$n.risk, n.event=f.survfit$n.event, | |
n.censor = f.survfit$n.censor, surv=f.survfit$surv, upper=f.survfit$upper, | |
lower=f.survfit$lower, strata=factor(f.strata)) | |
# remove temporary data | |
rm(f.strata) | |
# create first two rows (start at 1) for each strata | |
for(f.i in 1:length(f.survfit$strata)){ | |
# take only subset for this strata from data | |
f.subset <- subset(f.frame, strata==unlist(strsplit(names(f.survfit$strata)[f.i], "="))[2]) | |
# create first two rows (time: 0, time of first event) | |
f.start <- data.frame(time=c(0, f.subset$time[1]), n.risk=rep(f.survfit[f.i]$n, 2), | |
n.event=c(0,0), n.censor=c(0,0), surv=c(1,1), upper=c(1,1), lower=c(1,1), | |
strata=rep(unlist(strsplit(names(f.survfit$strata)[f.i], "="))[2],2)) | |
# add first two rows to dataset | |
f.frame <- rbind(f.start, f.frame) | |
# remove temporary data | |
rm(f.start, f.subset) | |
} | |
# reorder data | |
f.frame <- f.frame[order(f.frame$strata, f.frame$time), ] | |
# rename row.names | |
rownames(f.frame) <- NULL | |
} | |
# use different plotting commands dependig whether or not strata's are given | |
if("strata" %in% names(f.frame) == FALSE){ | |
# confidence intervals are drawn if not specified otherwise | |
if(f.CI=="default" | f.CI==TRUE ){ | |
# create plot with 4 layers (first 3 layers only events, last layer only censored) | |
# hint: censoring data for multiple censoring events at timepoint are overplotted | |
# (unlike in plot.survfit in survival package) | |
ggplot(data=f.frame) + geom_step(aes(x=time, y=surv), direction="hv") + | |
geom_step(aes(x=time, y=upper), directions="hv", linetype=2) + | |
geom_step(aes(x=time,y=lower), direction="hv", linetype=2) + | |
geom_point(data=subset(f.frame, n.censor>=1), aes(x=time, y=surv), shape=f.shape) | |
} | |
else { | |
# create plot without confidence intervalls | |
ggplot(data=f.frame) + geom_step(aes(x=time, y=surv), direction="hv") + | |
geom_point(data=subset(f.frame, n.censor>=1), aes(x=time, y=surv), shape=f.shape) | |
} | |
} | |
else { | |
if(f.CI=="default" | f.CI==FALSE){ | |
# without CI | |
ggplot(data=f.frame, aes(group=strata, colour=strata)) + | |
geom_step(aes(x=time, y=surv), direction="hv") + | |
geom_point(data=subset(f.frame, n.censor>=1), aes(x=time, y=surv), shape=f.shape) | |
} | |
else { | |
# with CI (hint: use alpha for CI) | |
ggplot(data=f.frame, aes(colour=strata, group=strata)) + | |
geom_step(aes(x=time, y=surv), direction="hv") + | |
geom_step(aes(x=time, y=upper), directions="hv", linetype=2, alpha=0.5) + | |
geom_step(aes(x=time,y=lower), direction="hv", linetype=2, alpha=0.5) + | |
geom_point(data=subset(f.frame, n.censor>=1), aes(x=time, y=surv), shape=f.shape) | |
} | |
} | |
} | |
# A short example | |
library(survival) | |
sf <- survfit(Surv(time, status) ~ ph.karno, data=lung) | |
ggsurv(sf) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment