# Anemia data
require(survival)
anemia <- read.delim("http://myweb.uiowa.edu/pbreheny/data/anemia.txt")
head(anemia)

# Surv objects
S.gvhd <- with(anemia, Surv(Time_gvhd, Status_gvhd))
S.mort <- with(anemia, Surv(Time, Status))
head(S.mort)

# survfit
Trt <- anemia$Trt
fit <- survfit(S.mort ~ Trt)
fit

# survfit plot 1
plot(fit, xlab="Days on study", ylab="Survival")

# survfit plot 2
col <- c("#FF4E37FF", "#008DFFFF")
plot(fit, xlab="Days on study", ylab="Survival", mark.time=FALSE,
     col=col, lwd=2, las=1, bty="n")
text(1000, 0.9, "MTX+CSP")
text(800, 0.4, "MTX")

# survfit plot 3
plot(fit, xlab="Days on study", ylab="Survival", mark.time=FALSE, conf.int=TRUE,
     col=col, las=1)
legend("bottomleft", legend=unique(Trt), col=col, lty=1)

# survfit plot 4
source("http://myweb.uiowa.edu/pbreheny/5710/f16/labs/ciband.R")
source("http://myweb.uiowa.edu/pbreheny/5710/f16/labs/nrisk.R")
col <- c("#FF4E37FF", "#008DFFFF")
plot(fit, xlab="Days on study", ylab="Survival", mark.time=FALSE,
     col=col, lwd=2, las=1, bty="n")
ciband(fit, col=c("#FF4E3740", "#008DFF40"))
nrisk(fit)
text(1400, 0.87, "MTX+CSP", xpd=TRUE)
text(1300, 0.47, "MTX")

# survdiffplot
require(rms)
fit.mort <- npsurv(S.mort ~ Trt)
survdiffplot(fit.mort, xlab="Days on Study", order=2:1)

# survdiff
survdiff(S.mort ~ Trt)
survdiff(S.gvhd ~ Trt)