[R-sig-ME] GLM: Difference between treatment groups for each colony (Tukey posthoc test)

[paul debes]

Hi Sophie,

If the treatment:colony term is NS, what about the colony term?

If the colony term is NS, you could also remove the colony term from the  
model and only leave the treatment term in the model (no need for a  
post-hoc test for the treatment term: it has only two levels: F-test/LRT  
results are sufficient).

If the colony term is also important (but the interaction is NS), the  
treatment effects are the same for each level of colony, and the colony  
effects differ from each other similarly within each treatment level. To  
find what colony levels differ from each other you will need to specify  
the colony term for your post-hoc test.
You can plot the means in two ways, depending on what is your study focus:  
1) two means for the two treatments averaged across colonies AND/OR 3  
means for the three colonies averaged across treatments (in 2 figures). 2)  
six means of the three colonies within each treatment, or vice versa (1  
figure).

To predict model means for many model types, I personally like using the  
library "lsmeans".
This might work (use "summary" to facilitate getting the results as  
data.frame for easier plotting):

library(lsmeans)
model.means.treatment = data.frame(summary(lsmeans(model, ~ treatment)))
model.means.colony = data.frame(summary(lsmeans(model, ~ colony)))
model.means.inter = data.frame(summary(lsmeans(model, ~ treatment:colony)))

Hope this helps,
Paul


On Mon, 07 Dec 2015 14:12:37 +0200, Sophie Waegebaert  
<sophie.waegebaert at gmail.com> wrote:

> Hi Paul,
>
> Thank you very much for the hint!
>
> The treatment:colony term is not significant, so I will not include it in
> the model.
>
> However, I get the same results for the Tukey posthoc test as I had  
> before.
> I am able to say that the DWV treatment significantly reduces the life
> expectancy, but I do not know the difference between the colonies. With  
> the
> results of the Tukey test I want to make a figure (see attachment), but I
> do not know how to do it. The asterisks are based on Tukey posthoc test.
>
> Do you know a way to do this in R?
>
> Thanks a lot!
>
> Cheers,
> Sophie
>
>
>
> 2015-12-07 12:30 GMT+01:00 paul debes <paul.debes at utu.fi>:
>
>> Hi Sophie and List,
>>
>> If you are interested in the treatment contrasts for each level of the
>> Colony factor (with three levels), it may be better to view your  
>> experiment
>> as a factorial and specify a simple linear model rather than the  
>> presently
>> specified linear mixed model. A random Colony term with only three  
>> levels
>> may also not provide the best estimate for the correlation of colony  
>> data
>> between treatments (i.e., as intraclass correlation) that is taken into
>> account when testing the general treatment term.
>>
>> Would this work for you?:
>>
>> fit_life = lm(last.scan ~ 1 + treatment*colony, data = data)
>>
>> In case the treatment:colony term is significant you could conduct the
>> additional pairwise tests of interest.
>>
>> Best,
>> Paul
>>
>>
>>
>> On Mon, 07 Dec 2015 13:02:36 +0200, Sophie Waegebaert <
>> sophie.waegebaert at gmail.com> wrote:
>>
>> Hi,
>>>
>>> I'm still learning how to use R and I have some trouble making using  
>>> Tukey
>>> posthoc tests. I have a dataset with 3 colonies (A, B and C). Each  
>>> colony
>>> is divided into 2 treatments: control and DWV. I want to run a GLM to  
>>> test
>>> wether there is a difference in life expectancy (last.scan) between the
>>> treatment groups for each of the colonies, but I do not know if I am  
>>> using
>>> the right strategy.
>>>
>>> I have taken 'treatment' as a fixed factor and 'colony' as a random
>>> factor:
>>>
>>> fit_life = lmer(last.scan~treatment + (1|colony), data =
>>> data)Anova(fit_life, type = 3) # Type of treatment has a significant
>>> effect on on the life expectancy.
>>>     Response: last.scan
>>>           Chisq Df Pr(>Chisq)
>>>     (Intercept) 106.976  1  < 2.2e-16 ***
>>>     treatment    25.373  1  4.724e-07 ***
>>>
>>> And this is the code I use to do a Tukey posthoc test:
>>>
>>> mcp = glht(fit_life, linfct = mcp(treatment = "Tukey"))
>>> summary(mcp)# DWV treatment significantly changes life expectancy (z =
>>> -9.734, p = < 2e-16)
>>>
>>> Is it possible to find the difference for each colony?
>>>
>>> Thanks a lot for an explanation or hint!
>>>
>>> Cheers,
>>>
>>> Sophie
>>>
>>>         [[alternative HTML version deleted]]
>>>
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>
>>
>> --
>> Paul Debes
>> DFG Research Fellow
>> University of Turku
>> Department of Biology
>> Itäinen Pitkäkatu 4
>> 20520 Turku
>> Finland
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>


-- 
Paul Debes
DFG Research Fellow
University of Turku
Department of Biology
Itäinen Pitkäkatu 4
20520 Turku
Finland