## simulate right censored data from a two state model
set.seed(100)
dat <- SimSurv(100)
# with(dat,plot(Hist(time,status)))
### marginal Kaplan-Meier estimator
kmfit <- prodlim(Hist(time, status) ~ 1, data = dat)
plot(kmfit)
# change time range
plot(kmfit,xlim=c(0,100))
# change scale of y-axis
plot(kmfit,percent=FALSE)
# change axis label and position of ticks
plot(kmfit,
xlim=c(0,100),
axis1.at=seq(0,100,30.25),
axis1.labels=0:(length(seq(0,100,30.25))-1),
xlab="Months",
axis2.las=2,
atrisk.at=seq(0,100,30.25),
atrisk.label="Patients")
# change background color
plot(kmfit,
xlim=c(0,100),
confint.citype="shadow",
col=1,
axis1.at=seq(0,100,30.25),
axis1.labels=0:(length(seq(0,100,30.25))-1),
xlab="Months",
axis2.las=2,
atrisk.at=seq(0,100,30.25),
atrisk.label="Patients",
background=TRUE,
background.fg="white",
background.horizontal=seq(0,1,.25/2),
background.vertical=seq(0,100,30.25),
background.bg=c("gray88"))
# change type of confidence limits
plot(kmfit,
xlim=c(0,100),
confint.citype="dots",
col=4,
background=TRUE,
background.bg=c("white","gray88"),
background.fg="gray77",
background.horizontal=seq(0,1,.25/2),
background.vertical=seq(0,100,30.25))
### Kaplan-Meier in discrete strata
kmfitX <- prodlim(Hist(time, status) ~ X1, data = dat)
plot(kmfitX)
# move legend
plot(kmfitX,legend.x="bottomleft",atRisk.cex=1.3)
## Control the order of strata
## unfortunately prodlim does not obey the order of
## factor levels
dat$group <- factor(cut(dat$X2,c(-Inf,0,0.5,Inf)),labels=c("Low","Intermediate","High"))
kmfitG <- prodlim(Hist(time, status) ~ group, data = dat)
plot(kmfitG)
## we want to re-order the labels such
## that in the legend and in the numbers at-risk
## "Low" comes before "Intermediate" before "High"
dat$group2 <- as.numeric(dat$group)
kmfitG2 <- prodlim(Hist(time, status) ~ group2, data = dat)
plot(kmfitG2,legend.legend=levels(dat$group),atrisk.labels=levels(dat$group))
# add log-rank test to legend
plot(kmfitX,
atRisk.cex=1.3,
logrank=TRUE,
legend.x="topright",
atrisk.title="at-risk")
# change atrisk labels
plot(kmfitX,
legend.x="bottomleft",
atRisk.cex=1.3,
atrisk.title="Patients",
atrisk.labels=c("X1=0","X1=1"))
### Kaplan-Meier in continuous strata
kmfitX2 <- prodlim(Hist(time, status) ~ X2, data = dat)
plot(kmfitX2,xlim=c(0,50))
# specify values of X2 for which to show the curves
plot(kmfitX2,xlim=c(0,50),newdata=data.frame(X2=c(-1.8,0,1.2)))
### Cluster-correlated data
library(survival)
cdat <- cbind(SimSurv(20),patnr=sample(1:5,size=20,replace=TRUE))
kmfitC <- prodlim(Hist(time, status) ~ cluster(patnr), data = cdat)
plot(kmfitC,atrisk.labels=c("Units","Patients"))
## simulate right censored data from a competing risk model
datCR <- SimCompRisk(100)
with(datCR,plot(Hist(time,event)))
### marginal Aalen-Johansen estimator
ajfit <- prodlim(Hist(time, event) ~ 1, data = datCR)
plot(ajfit) # same as plot(ajfit,cause=1)
# cause 2
plot(ajfit,cause=2)
# both in one
plot(ajfit,cause=1)
plot(ajfit,cause=2,add=TRUE,col=2)
### stacked plot
plot(ajfit,cause="stacked")
### conditional Aalen-Johansen estimator
ajfitX <- prodlim(Hist(time, event) ~ X1+X2, data = datCR)
plot(ajfitX,newdata=data.frame(X1=c(1,1,0),X2=c(4,10,10)))
plot(ajfitX,newdata=data.frame(X1=c(1,1,0),X2=c(4,10,10)),cause=2)
## stacked plot
plot(ajfitX,
newdata=data.frame(X1=0,X2=0.1),
cause="stacked",
legend.title="X1=0,X2=0.1",
legend.legend=paste("cause:",getStates(ajfitX$model.response)),
plot.main="Subject specific stacked plot")
Run the code above in your browser using DataLab