# lme versus aov repeated measures variance specification

I am performing a repeated measures ANOVA in R, where the experimental design is composed of a between-subjects factor (Prodotto, three levels), the within-subjects factor is time (Tempo, four levels&mdahs;an ordered factor) and each Product * Time is repeated on eleven patients (idPaziente).

The aov() syntax would be:

myAov<-aov(valore ~ Tempo*Prodotto + Error(idPaziente/Tempo),data=mydata)


whilst the lme() (from the nlme package) syntax would be

myModel <- lme(valore ~  Tempo*Prodotto,data=mydata, random= ~ 1 | idPaziente/Tempo)


In Finch (2014) Multilevel Modeling with R, I find a specification of a longitudinal model that would be written as:

myModel <- lme(valore ~  Tempo*Prodotto,data=mydata, random= ~ 1 | idPaziente)


I guess the equivalent aov version would have Error(idPaziente) instead of Error(idPaziente/Tempo).

1. What is the correct specification for my repeated measures model, both for aov() and lme()?

2. Why are (1 | idPaziente/Tempo) and Error(idPaziente/Tempo) more appropriate than (1 | idPaziente) and Error(idPaziente)?

3. How can I perform post hoc analysis taking on the Prodotto term?

• This is borderline, but IMO it may have just enough statistical content to be considered on topic here. Dec 26 '15 at 2:30
• Hi Giorgio - Just out of curiosity, did my answer worked for the problem you were experiencing? Or can I update it somehow? Dec 26 '15 at 21:37

This depends on your nesting structure (example dataset):

set.seed(123)
df <- data.frame(
idPaziente = rep(c(1:12), each = 11),
Tempo = rep(c(1:4),33),
Prodotto = rep(letters[1:3], 11, each = 4),
valore = rnorm(132, 6, 1)
)


Using random = ~ 1 | idPaziente/Tempo:

fit <- lme(valore ~  Tempo * Prodotto, data = df, random = ~ 1 | idPaziente/Tempo)


To get the ANOVA table you could do:

library(car)
Anova(fit)


If your idPaziente and Tempo is explicitly nested,

df\$newID <- with(df, paste(idPaziente, Tempo, sep = "-"))


you only need to specify random = ~ 1 | newID

fit2 <- lme(valore ~  Tempo * Prodotto, data = df, random = ~ 1 | newID)

Anova(fit2)


then:

summary(fit) and summary(fit2) will give you the same results. Also have a look at the Number of Groups: statement. This is a good check to see whether your nesting worked out correctly and reflects your data structure.

If Prodotto is significant, you could follow up with:

library(lsmeans)
lsmeans(fit, pairwise ~ Prodotto)
lsmeans(fit2, pairwise ~ Prodotto)


And the same should hold for the aov() function, too.