On Sep 15, 2015, at 12:09 PM, Huot, Matthieu wrote: > Hi Tom > > I know the post is over 7-8 years old but I am having the same question. How > to do a post-hoc test like TukeyHSD on coxph type output.
Create a new variable using the `interaction`-function, apply you contrasts to that object, and that should let you side-step all the errors of the person trying to follow Crawley's book. -- David. > > Have you received any info in this matter? > Thanks > Matthieu > > Looking for Post-hoc tests (a la TukeyHSD) or interaction-level independent > contrasts for survival analysis. > Thomas Oliver toliver at stanford.edu > <mailto:r-help%40r-project.org?Subject=Re%3A%20%5BR%5D%20Looking%20for%20Post-hoc%20tests%20%28a%20la%20TukeyHSD%29%20or%20interaction-level%0A%20independent%20contrasts%20for%20survival%20analysis.&In-Reply-To=%3C6.2.5.6.2.20080429131826.04deaeb0%40stanford.edu%3E> > Tue Apr 29 23:12:33 CEST 2008 > > * Previous message: [R] AIC extract and comparison > <https://stat.ethz.ch/pipermail/r-help/2008-April/160916.html> > * Next message: [R] for loop in nls function > <https://stat.ethz.ch/pipermail/r-help/2008-April/160886.html> > * Messages sorted by: [ date > ]<https://stat.ethz.ch/pipermail/r-help/2008-April/date.html#160885> [ thread > ]<https://stat.ethz.ch/pipermail/r-help/2008-April/thread.html#160885> [ > subject > ]<https://stat.ethz.ch/pipermail/r-help/2008-April/subject.html#160885> [ > author ]<https://stat.ethz.ch/pipermail/r-help/2008-April/author.html#160885> > > ________________________________ > > Hello all R-helpers, > > > > I've performed an experiment to test for differential effects of > > elevated temperatures on three different groups of corals. I'm > > currently performing a cox proportional hazards regression with > > censoring on the survivorship (days to mortality) of each individual > > in the experiment with two factors: Temperature Treatment (2 levels: > > ambient and elevated) and experimental group (3 levels: say 1,2,3). > > > > In my experiment, all three groups survived equally well in the > > ambient control treatment, but two of three of the groups succumbed > > to heat stress in the elevated temperature treatment. I can see that > > the third group had a small degree of mortality, but it appears to be > > significantly less than the other two and may be not significantly > > different from the ambient controls. > > > > I would like to ask three questions: 1) Is group 3 different from > > controls? 2) Is group 3 different from group 1 and/or group 2 in the > > elevated treatment? and 3) are groups 1 and 2 different from each > > other in the elevated treatment? > > > > Because I'm testing for differential effects among the elevated > > temperature treatment group, and "I've seen the data" by now, the > > analysis would be easiest for me if I performed a responsible > > multiple comparisons test like TukeyHSD to see how each of the six > > Treatment:Group subgroups compared to each other. TukeyHSD does not > > appear to be defined for outputs from the function coxph -- (see > > survival library). > > > > cph <- coxph(Surv(DayToMort, Censor) ~ Treatment*Group, data=subb) > > > > --> Does anyone know of an implementation of TukeyHSD that would > > work, or of another post-hoc multiple comparison test? > > > > I believe that another responsible tack would be to clearly define > > the contrasts I'd like to make within the interaction term. However > > this has yet to work as fully as I'd like it. > > > > I've successfully set the contrasts matrix for the three-level > > factor "Group" following Crawley's The R Book. > > > > cmat<-cbind(c(-1,1,0),c(0,-1,1)) > > contrasts(subb$Group)<-cmat > > contrasts(subb$Group) > > > > By setting these contrasts and then looking at the interaction terms > > in the coxph model, this allows me to compare groups _within_ each > > separate treatment, and confirms both that #2) that groups 1 and 2 > > are not sig. different in the elevated treatment, and #3) the group3 > > corals survived significantly better than the other groups in the > > elevated treatment. BUT it does not allow me to say if the group3 > > survival is or is not different from its own control.(#1 above). > > > > To make this comparison, I've tried setting the contrast matrix for > > the Treatment:Group interaction term, with no success. Whenever I > > attempt to do so, I run the code below: > > > > cmat<-cbind(c(-1,0,0,1,0,0),c(0,-1,0,0,1,0),c(0,0,-1,0,0,1),c(0,0,0,-1,1,0),c(0,0,0,0,-1,1)) > > #Build a matrix > > rownames(cmat)<-rownames(contrasts(subb$Treatment:subb$Group)) #give > > cmat the correct rownames > > colnames(cmat)<-colnames(contrasts(subb$Treatment:subb$Group)) > > #give cmat the correct colnames > > contrasts(subb$Treatment:subb$Group)<-cmat #try to assign cmat > > > > and I get this error message: > > Error in contrasts(subb$Treatment:subb$Group, ) <- cmat : > > could not find function ":<-" > > > > Alternatively I could run: > > options(contrasts=c("contr.sum","contr.poly")) > > where: > > contrasts(subb$Treatment:subb$Group) > > > > [,1] [,2] [,3] [,4] [,5] > > Con:1 1 0 0 0 0 > > Con:2 0 1 0 0 0 > > Con:3 0 0 1 0 0 > > Exp:1 0 0 0 1 0 > > Exp:2 0 0 0 0 1 > > Exp:3 -1 -1 -1 -1 -1 > > > > But even that doesn't appear to affect the output of : > > > > cph <- coxph(Surv(DayToMort, Censor) ~ Treatment*Group, data=subb) > > > > > > --> Is what I'm trying to do statistically invalid and R is trying to > > quietly save me from statistical destruction, or it is just being a > > pain? Is there a way around it? > > > > --> Any other suggestions? > > > > Many Thanks in Advance, > > > > Tom Oliver > > > > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > 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. David Winsemius Alameda, CA, USA ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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.