R help - Feb 2008 - Nested frailty model

Dear R-help,
I am trying to estimate a Cox model with nested effects, or better

h(t,v,w)=v*w*h0(t)*exp(B'x)

where h(t,v,w) is the individual hazard function
       w and v are both frailty terms (gamma or normal distributed)
I have 12 clusters and for each one of them I would like to associate a 
realization of v, while w is a random effect for the whole population. 
At the population level I should integrate out both the frailty 
distributions.
This type of model is usually called nested frailty model.
I would like to know if there is a package which could allow me to 
estimate it.

Thank you very much indeed
Silvia

Hi,

the package 'survival' provides the possibility to estimate Cox PH 
models with a frailty term.
library(survival)
?frailty

If this is not what you are looking for, there is also a package called 
'frailtypack'. But I have no experience with that.

I hope this helps.

Best,
Roland

P.S. Bocconi and survival analysis? Sounds like you work together with 
(or take a course from) Francesco Billari. If this is the case, please 
say hello to him from me.

ah06981 at studbocconi.it wrote:
> Dear R-help,
> I am trying to estimate a Cox model with nested effects, or better
> 
> h(t,v,w)=v*w*h0(t)*exp(B'x)
> 
> where h(t,v,w) is the individual hazard function
>        w and v are both frailty terms (gamma or normal distributed)
> I have 12 clusters and for each one of them I would like to associate a 
> realization of v, while w is a random effect for the whole population. 
> At the population level I should integrate out both the frailty 
> distributions.
> This type of model is usually called nested frailty model.
> I would like to know if there is a package which could allow me to 
> estimate it.
> 
> Thank you very much indeed
> Silvia
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

You want to use coxme, currently found in the kinship package.  It uses lme 
style input, e.g.
 	coxme(Surv(time, status) ~ age + weight, data=mydata,
 		random= ~1 | school/teacher)
 		

  One note: your message seems to imply that you want to use group/subject as 
the random effect.  A Cox model with one random effect per subject is not stable
-- actually an one where df(random effect) approaches the number of events is 
not --- and you will likely get a variance estimate whose confidence interval is
(0, huge number).
  
  	Terry Therneau