I have the following between-within anova: aovn1 <- aov(amplitude ~ stereo*site*stimulus + Error(subject/(site*stimulus)), stereon1) This works fine. BUT I need to do Tukey HSD multiple comparisons, and the aov() approach won't work. So I am trying the method posted on r-help: lmen1 <- lme(amplitude ~ stereo*site*stimulus, random ~1|subject/(site*stimulus), stereon1) This doesn't work:> lmen1 <- lme(amplitude ~ stereo*site*stimulus, random = ~1|subject/(site*stimulus), stereon1)Error in getGroups.data.frame(dataMix, groups) : Invalid formula for groups Please tell me what to do to get the lme() fit to work. After I get lmen1, I will do: anova(lmen1) summary(lmen1) summary(glht(lmen1, linfct=mcp(V="Tukey"))) Please tell me if that sounds right. Thanks very much for any help! Bill ============= This is the example I am following (with modifications for my expt design) You want to use lme() in package nlme, then glht() in the multcomp package. This will give you multiplicity adjusted p-values and confidence intervals. ## Example require(MASS) ## for oats data set require(nlme) ## for lme() require(multcomp) ## for multiple comparison stuff Aov.mod <- aov(Y ~ N + V + Error(B/V), data = oats) Lme.mod <- lme(Y ~ N + V, random = ~1 | B/V, data = oats) summary(Aov.mod) anova(Lme.mod) summary(Lme.mod) summary(glht(Lme.mod, linfct=mcp(V="Tukey")))