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).  
I would like to ask: 


*

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

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

*How can I perform post hoc analysis taking on the Prodotto term? 
 A: 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.
