library(survival) library(boot) data=NULL lambda=NULL result=NULL
pat=rep(1:102,each=1) trt=rep(c(1,0),51) status=rep(1,102) site=rep(1:51,
each=2) nr.datasets=100 seed=2006 beta=log(1/2) for (i in 1:51) {
lambda[i]=1+((3-1)/50)*(i-1)} lambda1=rep(lambda, each=2)
dummy=rep(c(exp(beta),1),51) elf=lambda1*dummy r=70 #the number of bootstrap
replicates #mymodels=function(dataset.number){ #seed=2006*dataset.number
time=rexp(102, rate=elf) data=cbind(pat,trt,status,site,time)
data1=data.frame(data) #### Cox proportional hazard model stratified on site
model1.cox=coxph(Surv(time,status) ~ trt + strata(site), data=data1,
iter.max=500) model1.surv=survfit(model1.cox) #### Cox proportional hazard
model with dummy covariates for site model2.cox=coxph(Surv(time,status) ~ trt
+ factor(site), data=data1, iter.max=500) model2.surv=survfit(model2.cox)
beta.fun1 = function(dataname) { cox =
coef(coxph(Surv(dataname$time,dataname$status) ~
dataname$trt+strata(dataname$site))) } beta.fun2 = function(dataname) {
cox =coef(coxph(Surv(dataname$time,dataname$status) ~
dataname$trt+factor(dataname$site))) #cox=cox[1] } dataname=data1 ##
Calculate Cox's regression coefficients' standard error and bias using
the ## ordinary Bootstrap # set.seed(1234)#set the random start cox.ord1 =
censboot(data=dataname,statistic=beta.fun1,R=r,F.surv=model1.surv,
cox=model1.cox,sim="ordinary") cox.ord1 cox.ord2 =
censboot(data=dataname,statistic=beta.fun2,R=r,F.surv=model2.surv,
cox=model2.cox,sim="ordinary") cox.ord2 join=list(model1,
model2) #return(join) #} #for (i in 1:2){ #result[i]=mymodels}
I run the above code but for the following part:
cox.ord2 = censboot(data=dataname,statistic=beta.fun2,R=r,F.surv=model2.surv,
cox=model2.cox,sim="ordinary") cox.ord2
I get an error and says that it did not converge. I am not sure whats wrong in
my coding.
I want to calculate the bias. Any help is very much appreciated.
Thanks,
Orlando
---------------------------------
[[alternative HTML version deleted]]
