R help - Mar 2008 - coxme - fitting random treatment effect nested within centre

Dear all,

I am using "coxme" function in Kinship library to fit random treatment
effect nested within centre. I got 3 treatments (0,1,2) and 3 centres. I used
following commands, but got an error.
> ugroup=paste(rep(1:3,each=3),rep(0:2,3),sep='/')
>
mat1=bdsmatrix(rep(c(1,1,1,1,1,1,1,1,1),3),blocksize=rep(3,3),dimnames=list(ugroup,ugroup))
>
mat2=bdsmatrix(rep(c(0,0,0,0,0,0,0,0,1),3),blocksize=rep(3,3),dimnames=list(ugroup,ugroup))
> group=paste(dat1$centre,dat1$treat,sep='/')
> coxme(Surv(time,status) ~ as.factor(treat), data=dat1,random=
~1|group,varlist=list(mat1,mat2),rescale=F,pdcheck=FALSE)
Error in coxme.fit(X, Y, strats, offset, init, control, weights = weights, : 
Random effects variance is not spd
 
Could anyone help me correcting this error? 
Many thanks in advance.
Ruwanthi


     
____________________________________________________________________________________
Looking for last minute shopping deals?

---- included message ----
 Dear all,

I am using "coxme" function in Kinship library to fit random treatment
effect
nested within centre. I got 3 treatments (0,1,2) and 3 centres. I used following
commands, but got an error. 
> ugroup=paste(rep(1:3,each=3),rep(0:2,3),sep='/')
> 
mat1=bdsmatrix(rep(c(1,1,1,1,1,1,1,1,1),3),blocksize=rep(3,3),dimnames=list(ugro
up,ugroup))
> 
mat2=bdsmatrix(rep(c(0,0,0,0,0,0,0,0,1),3),blocksize=rep(3,3),dimnames=list(ugro
up,ugroup))
> group=paste(dat1$centre,dat1$treat,sep='/')
> coxme(Surv(time,status) ~ as.factor(treat), data=dat1,random= 
~1|group,varlist=list(mat1,mat2),rescale=F,pdcheck=FALSE)

Error in coxme.fit(X, Y, strats, offset, init, control, weights = weights, : 
Random effects variance is not spd
 
Could anyone help me correcting this error? 
Many thanks in advance.
Ruwanthi

------- end inclusion -----

 The "build your own bdsmatrix" style in coxme is for variance
structures that
are not built in, like pedigree data.  Your problem is one that coxme can do 
directly, so it is easier to call the routine simply:
 
  fit <- coxme(Surv(time, status) ~ factor(treat), data=data1,
  	         random = ~1 | centre/treat)
  	       
	Terry Therneau