[R] coxph models for insects
therneau at mayo.edu
Tue Oct 9 15:07:40 CEST 2007
You have an interesting problem, and a serious (reliable) consultation would
take more time than I have to give at the moment. Which is to say that you
should take these comments with a grain of salt.
First, I don't think that you have censored data. You have 2 subdistribution
functions F1(t) and F2(t), F1(t) + F2(t) = F(t) = the "time to endpoint"
distribution. With censored data you would have some insects whose endpoint had
not been observed, e.g., it's time to write the paper and some of the durn
things have neither emerged nor died yet.
How best to model this data is a larger and harder question. What you have
done with coxph, creating an artificial censoring at the time of the competing
event, is a first step in competing risks modeling. You might look at the
appropriate chapter in Kalbfleisch and Prentice or a similar survival textbook
to get a feeling about how to think about the fitted coefficients. I would
guess that there is a lot of work on competing risks outside of the survival
(censored) data literature, in financial markets for instance, that is
appropriate; but I am not familiar with that literature. If they use linear
models, as I suspect they might, then that leads to appropriate answers to your
question about blocking.
For a coxph model, the code for grouped variance structures a la GEE models is
not very mature. One could allow for a single level, eg correlation within
chamber by adding "+ cluster(chamber)" to the model line, but more complex
correlation structures are not supported. As an alternative one could use a
mixed effects Cox model via the coxme function (currently distributed as a part
of the kinship package, but will be added to the survival package "Real Soon
Now"), with blocks as random effects.
Your use of the weight argument in coxph is correct.
I'm attempting to run a Cox proportional hazards function on survival data
from insects and I have a few questions.
My current command that I'm using to call the model is as follows (using
coxph() from the survival library):
coxph(Surv(day, censor) ~ treatment + room + chamber %in% treatment, data =
Day indicates which day a particular observation occurred, and censor indicates
the type of observation; that is, mortality or emergence.
My first question is: how do I code my censor variable? I am interested in
treatment effects on mortality, and the mortality observations are currently
coded as a 1 in the censor variable. Emergence is another mechanism by which
the insects in my study may leave the system, similar to a patient opting
out in the middle of a clinical study, and are coded as a 0.
I also have a nested and block variable that I'm unsure of how to deal with.
Chamber is nested within treatment, where each chamber has a maximum of 4
observations or insects, all within 1 treatment. The block variable is room,
where each room contains ~40 chambers each of which has been randomly
assigned a treatment. I believe that I have chamber %in% treatment right in
representing my data, but how do I account for room as a block?
More information about the R-help