tests/bladder.R

options(na.action=na.exclude) # preserve missings
options(contrasts=c('contr.treatment', 'contr.poly'),
                     show.signif.stars=FALSE) #ensure constrast type
library(survival)

#
# Fit the models found in Wei et. al.
#
wfit <- coxph(Surv(stop, event) ~ (rx + size + number)* strata(enum),
		 cluster=id, bladder, ties='breslow')
wfit

# Check the rx coefs versus Wei, et al, JASA 1989
rx <- c(1,4,5,6)  # the treatment coefs above
cmat <- diag(4); cmat[1,] <- 1;          #contrast matrix
wfit$coef[rx] %*% cmat                   # the coefs in their paper (table 5)
t(cmat) %*% wfit$var[rx,rx] %*% cmat  # var matrix (eqn 3.2)

# Anderson-Gill fit
fita <- coxph(Surv(start, stop, event) ~ rx + size + number, cluster=id,
		  bladder2,  ties='breslow')
summary(fita)

# Prentice fits.  Their model 1 a and b are the same
fit1p  <- coxph(Surv(stop, event) ~ rx + size + number, bladder2,
		subset=(enum==1), ties='breslow')
fit2pa <- coxph(Surv(stop, event) ~ rx + size + number, bladder2,
		subset=(enum==2), ties='breslow')
fit2pb <- coxph(Surv(stop-start,  event) ~ rx + size + number, bladder2,
		   subset=(enum==2), ties='breslow')
fit3pa <- coxph(Surv(stop, event) ~ rx + size + number, bladder2,
		subset=(enum==3), ties='breslow')
 #and etc.
fit1p
fit2pa
fit2pb
fit3pa
rm(rx, cmat, wfit, fita, fit1p, fit2pa, fit2pb, fit3pa)

Try the survival package in your browser

Any scripts or data that you put into this service are public.

survival documentation built on Aug. 24, 2021, 5:06 p.m.