[R-sig-ME] lme with interaction: how to decide which treatment is more efficaceous? how to deal with p value column of the lme summary?
Thierry Onkelinx
thierry.onkelinx at inbo.be
Tue Mar 22 11:26:59 CET 2016
Dear Martina,
I would reparametrise the model.
lme(var ~ treatment + treatment:period+treatment:side,
method="ML",
random=list(IDlog=~1), na.action=na.omit,
data=changes)
Another option is to use glht() from the multcomp package and specify
specific contrasts.
Best regards,
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium
To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to say
what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey
2016-03-21 16:58 GMT+01:00 Martina Giovannella via R-sig-mixed-models <
r-sig-mixed-models op r-project.org>:
>
> I am using an lme to check the effect of 3 types of treatment (A,B, C) on
> my variable, over time. I do not expect any effect of treatment A, which is
> a fake treatment.
>
> I have measurement over 5 periods,on 2 sides.
>
> My question is: for each treatment and each period, is my variable
> different from the reference point (treatment A, period 1)?
>
> My model is
>
> lme(var ~ treatment*period+side,
> method="ML",
> random=list(IDlog=~1), na.action=na.omit,
> data=changes)
>
> This is the anova of the model
>
> numDF denDF F-value p-value
> (Intercept) 1 473 79.36094 <.0001
> treatment 2 473 3.49473 0.0311
> period 4 473 12.51296 <.0001
> side 1 473 12.16210 0.0005
> treatment:period 8 473 2.02865 0.0416
>
> and the summary for fixed effects
>
> Fixed effects: var ~ treatment * period + hemisphere
> Value Std.Error DF t-value p-value
> (Intercept) 4.622038 3.180983 473 1.453022 0.1469
> treatmentB -1.376755 3.703398 473 -0.371754 0.7102
> treatmentC -1.113021 3.703398 473 -0.300540 0.7639
> period2 3.799946 4.168792 473 0.911522 0.3625
> period3 6.124463 4.168792 473 1.469122 0.1425
> period4 4.309267 4.168792 473 1.033697 0.3018
> period5 5.482068 4.168792 473 1.315026 0.1891
> sideright -4.278672 1.226887 473 -3.487420 0.0005
> treatmentB:period2 4.059350 5.190104 473 0.782133 0.4345
> treatmentC:period2 7.508426 5.190104 473 1.446681 0.1486
> treatmentB:period3 1.965207 5.190104 473 0.378645 0.7051
> treatmentC:period3 5.312525 5.190104 473 1.023587 0.3066
> treatmentB:period4 9.031016 5.190104 473 1.740045 0.0825
> treatmentC:period4 7.819397 5.190104 473 1.506597 0.1326
> treatmentB:period5 13.620365 5.190104 473 2.624295 0.0090
> treatmentC:period5 4.224340 5.215896 473 0.809897 0.4184
>
> The anova states the treatment has an effect, period as well, and there is
> an interaction between treatment and period. Sides are different. This is
> what I expected.
>
> How can I use the p values columns to state for which treatment and which
> period I see a difference from the reference value?
>
> What I would like to retrieve is the numbers I obtain when I test each
> treatment separately. I know how to obtain the mean, what I have to sum up
> in the value column. What about the p values?
>
> Treatment A
>
> Fixed effects: var ~ period + hemisphere
> Value Std.Error DF t-value p-value
> (Intercept) 1.849555 3.768596 94 0.4907810 0.6247
> period2 3.799946 4.174739 94 0.9102236 0.3650
> period3 6.124463 4.174739 94 1.4670291 0.1457
> period4 4.309267 4.174739 94 1.0322244 0.3046
> period5 5.482068 4.174739 94 1.3131524 0.1923
> sideright -1.142925 2.640337 94 -0.4328710 0.6661
>
> Anova
>
> numDF denDF F-value p-value
> (Intercept) 1 94 4.967367 0.0282
> period 4 94 0.654946 0.6248
> side 1 94 0.187377 0.6661
>
> Treatment B
>
> Fixed effects: var ~ period + hemisphere
> Value Std.Error DF t-value p-value
> (Intercept) 4.674863 2.938724 175 1.590780 0.1135
> period2 7.859296 3.241315 175 2.424724 0.0163
> period3 8.089670 3.241315 175 2.495799 0.0135
> period4 13.340283 3.241315 175 4.115700 0.0001
> period5 19.102433 3.241315 175 5.893420 0.0000
> sideright -7.137831 2.049988 175 -3.481889 0.0006
>
> Anova
>
> numDF denDF F-value p-value
> (Intercept) 1 175 34.37766 <.0001
> period 4 175 9.60042 <.0001
> side 1 175 12.12355 6e-04
>
> Treatment C
>
> Fixed effects: var ~ period + hemisphere
> Value Std.Error DF t-value p-value
> (Intercept) 2.936039 2.166417 173 1.355251 0.1771
> period2 11.308372 2.392119 173 4.727346 0.0000
> period3 11.436988 2.392119 173 4.781112 0.0000
> period4 12.128664 2.392119 173 5.070260 0.0000
> period5 9.949632 2.428331 173 4.097314 0.0001
> sideright -3.132715 1.520530 173 -2.060278 0.0409
>
> Anova
>
> numDF denDF F-value p-value
> (Intercept) 1 173 58.13208 <.0001
> period 4 173 8.99083 <.0001
> side 1 173 4.24474 0.0409
>
> From this test it is very clear that for treatment A, var is not bigger in
> any period respect to the reference point, while it is for treatment B and
> C; right side smaller than left.
>
> I want to obtain this information using the first model, with the 3
> treatments.
>
> _______________________________________________
> R-sig-mixed-models op r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list