If you want to stick with the aov()
function you can use the emmeans
package which can handle aovlist
(and many other) objects.
library("emmeans")
# set orthogonal contrasts
options(contrasts = c("contr.sum", "contr.poly"))
aov_velocity <- aov(Velocity ~ Material + Error(Subject / Material), data = scrd)
After creating an emmGrid
object as follows
emm <- emmeans(aov_velocity, ~ Material)
it is very easy to get all (post hoc) pairwise comparisons using the pairs()
function or any desired contrast using the contrast()
function of the emmeans
package. Multiple-testing adjustments can be achieved via the adjust
argument of these functions:
pairs(emm) # adjust argument not specified -> default p-value adjustment in this case is "tukey"
For more information on this I found the detailed emmeans vignettes and the documentation to be very helpful.
Also, you can find a complete (reproducible) example including a description on how to get the correct contrast weights in my answer here.
Note, however, that using a univariate model for the post hoc tests can result in anti-conservative p-values if sphericity is violated.