# R commands for Section 7.2
# Only the PBC data are
considered, not the UNOS data.
# Read the PBC data into a data
frame (needs only to be done once if the workspace is saved)
pbc=read.table("http://www.uio.no/studier/emner/matnat/math/STK4900/v11/pbc.txt", sep="\t",header=T,na.strings=".")
# Cox model
for treatment and bilirubin (cf
Tables 7.4 and 7.5, pages 219 and 220):
library(survival)??? # load survival library
pbc.fit=coxph(Surv(years,status)~rx+bilirubin, data=pbc,method="breslow")?
????? # default in R is method="efron" which gives slightly different results when
there are tied survival times
summary(pbc.fit)
# Cox model
with histology as a categorical covariate?
(cf Table 7.6, page 222):
pbc.hist=coxph(Surv(years,status)~factor(histol), data=pbc,method="breslow")?
summary(pbc.hist)
# Cox model
with age in one-year units (cf Table 7.9, page 224):
pbc.age=coxph(Surv(years,status)~age,
data=pbc,method="breslow")?
summary(pbc.age)
# Cox model
with age in five-year units (cf Table 7.10, page 225):
pbc.age5=coxph(Surv(years,status)~I(age/5),
data=pbc,method="breslow")?
summary(pbc.age5)
# Unadjusted
Cox model for bilrubin (cf
Table 7.11, page 226):
pbc.bil=coxph(Surv(years,status)~bilirubin, data=pbc,method="breslow")?
summary(pbc.bil)
# Adjusted Cox
model for bilrubin (cf
Table 7.12, page 226):
pbc.biladj=coxph(Surv(years,status)~bilirubin+edema+hepatom+spiders,
data=pbc,method="breslow")?
summary(pbc.biladj)
# Cox model
with interaction (cf Table 7.13, page 227):
pbc.int=coxph(Surv(years,status)~rx+hepatom+rx:hepatom, data=pbc,method="breslow")?
summary(pbc.int)
# Predicted
survival curve for PBC covariate pattern (cf Figure
7.6, page 231):
pbc.bilhept=coxph(Surv(years,status)~bilirubin+hepatom, data=pbc,method="breslow")?
plot(survfit(pbc.bilhept,
newdata=data.frame(bilirubin=4.5,hepatom=1)),conf.int=F,xlim=c(0,15))
??????????????????????????????????????????????????????????????????????????????????????