jeudi 28 mai 2015

R, Pairwise comparison of random variable in mixed model

We measured temperatures of a pond repeatedly every day at each hour for a month at two different depths (i.e., top and bottom). We want to see if the temperatures at the top of the pond are significantly different from the bottom and if so at what hours. Initially I did a two-way anova in R:

aov.result <- aov(temp ~ depth * hour, data = pondtemp)

followed by post hoc tukey hsd test:

tukey.result <- TukeyHSD(aov.result)

The pairwise comparison of the depth*hour interaction term is what I need to see which hours have significantly different temperatures between top and bottom. This worked out well but someone pointed out that since it is a repeated measure it does not satisfy the assumption of independence. Therefore I tried using a linear mixed model. I took the depth as the fixed variable and figured the hour (since it has multiple observations) in a month should be the random variable:

pondmdl <- lmer(temp ~ depth + (1+variable|hour), data = pondtemp)

And used glht package for post-hoc:

summary(glht(pondmdl, mcp(depth = "Tukey")))

However, this does not allow me to do the pair wise comparison I want to do (i.e., comparing Top Hour 1 to Bottom Hour 0 -23)

I found one way by introducing an interaction factor:

pondtemp$depth.hour <- interaction (pondtemp$depth, pondtemp$hour)

and then using this in my model and glht function:

pondmdl <- lmer(temp~depth.hour + (1+depth|hour), data = pondtemp)
summary(glht(pondmdl, mcp(depth.hour = "Tukey")))

However, I'm not sure if I can allow fixed and random variables interact like that and still use the same random variable in the random variable error term.

Please advise what is my best option. Thank you !




Aucun commentaire:

Enregistrer un commentaire