# PRACTICAL EXERCISE 8 # ==================== # Read data and disregard the individuals who smoke pipe or cigar: norw.death=read.table("http://folk.uio.no/borgan/abg-2008/data/causes_death.txt", header=T) norw.death=norw.death[norw.death\$smkgr!=6,] # We will use the survival package, so this has to be loaded library(survival) # Question a # ---------- # We fit a Cox model with sex and smoking habits as categorical covariates # and systolic blodd pressure and body mass index as numeric covariates: fit.b=coxph(Surv(agestart,agestop,dead)~factor(sex)+factor(smkgr)+sbp+bmi, data=norw.death) summary(fit.b) # The model assumes a log-linear effect of the numeric covariates and # proportionality of all covariates # Question b # ---------- # We fit a model with smoothing spline for sbp: fit.b=coxph(Surv(agestart,agestop,dead)~factor(sex)+factor(smkgr)+pspline(sbp)+bmi, data=norw.death) print(fit.b) par(mfrow=c(2,2)) termplot(fit.b,se=T) # The effect of systolic blood pressure is fairly log-linear. # Question c # ---------- # We fit a model with smoothing spline for bmi: fit.c=coxph(Surv(agestart,agestop,dead)~factor(sex)+factor(smkgr)+sbp+pspline(bmi), data=norw.death) print(fit.c) par(mfrow=c(2,2)) termplot(fit.c,se=T) # The effect of body mass index in not log-linear. It seems difficult to find # a useful trasformation of body mass index, so we will work with grouped bmi. # We make a categorical variable for groupe bmi (cf. problem text) and fit a model with this # covariate using the igroup [20,25) as reference: norw.death\$gbmi=cut(norw.death\$bmi,c(0,20,25,30,60),right=F) norw.death\$gbmi=relevel(norw.death\$gbmi,ref="[20,25)") fit.c2=coxph(Surv(agestart,agestop,dead)~factor(sex)+factor(smkgr)+sbp+factor(gbmi), data=norw.death) summary(fit.c2) # Question d # ---------- # We perform a test and make plots to check proportionality for the model with grouped body mass index cox.zph(fit.c2,transform='log') par(mfrow=c(3,3)) plot(cox.zph(fit.c2,transform='log')) # The assumption of proportionality seems reasonable for all covariates # (There may be some problem with proportionality for smoking group 5, but this # is not serious enough to "force us" to modify the model.)