aov_car {afex}  R Documentation 
These functions allow convenient specification of any type of ANOVAs (i.e.,
purely withinsubjects ANOVAs, purely betweensubjects ANOVAs, and mixed
betweenwithin or splitplot ANOVAs) for data in the long format
(i.e., one observation per row). If the data has more than one observation
per individual and cell of the design (e.g., multiple responses per
condition), the data will be automatically aggregated. The default settings
reproduce results from commercial statistical packages such as SPSS or SAS.
aov_ez
is called specifying the factors as character vectors,
aov_car
is called using a formula similar to aov
specifying an error strata for the withinsubject factor(s), and aov_4
is called with a lme4like formula (all ANOVA functions return
identical results). The returned object can be passed to e.g., emmeans
for further analysis (e.g., followup tests, contrasts, plotting, etc.).
These functions employ Anova
(from the car package)
to provide test of effects avoiding the somewhat unhandy format of
car::Anova
.
aov_car(
formula,
data,
fun_aggregate = NULL,
type = afex_options("type"),
factorize = afex_options("factorize"),
observed = NULL,
anova_table = list(),
include_aov = afex_options("include_aov"),
return = afex_options("return_aov"),
...
)
aov_4(
formula,
data,
observed = NULL,
fun_aggregate = NULL,
type = afex_options("type"),
factorize = afex_options("factorize"),
return = afex_options("return_aov"),
anova_table = list(),
include_aov = afex_options("include_aov"),
...,
print.formula = FALSE
)
aov_ez(
id,
dv,
data,
between = NULL,
within = NULL,
covariate = NULL,
observed = NULL,
fun_aggregate = NULL,
transformation,
type = afex_options("type"),
factorize = afex_options("factorize"),
return = afex_options("return_aov"),
anova_table = list(),
include_aov = afex_options("include_aov"),
...,
print.formula = FALSE
)
formula 
A formula specifying the ANOVA model similar to

data 
A 
fun_aggregate 
The function for aggregating the data before running the
ANOVA if there is more than one observation per individual and cell of the
design. The default 
type 
The type of sums of squares for the ANOVA. The default is given
by 
factorize 
logical. Should between subject factors be factorized (with
note) before running the analysis. The default is given by

observed 

anova_table 

include_aov 
Boolean. Allows suppressing the calculation of the aov
object. If TRUE the aov model is part of the returned 
return 
What should be returned? The default is given by

... 
Further arguments passed to 
print.formula 

id 

dv 

between 

within 

covariate 

transformation 
In 
aov_ez
will concatenate
all betweensubject factors using *
(i.e., producing all main effects
and interactions) and all covariates by +
(i.e., adding only the main
effects to the existing betweensubject factors). The withinsubject factors
do fully interact with all betweensubject factors and covariates. This is
essentially identical to the behavior of SPSS's glm
function.
The formula
s for aov_car
or aov_4
must contain a single
Error
term specifying the ID
column and potential
withinsubject factors (you can use mixed
for running
mixedeffects models with multiple error terms). Factors outside the
Error
term are treated as betweensubject factors (the withinsubject
factors specified in the Error
term are ignored outside the
Error
term; in other words, it is not necessary to specify them
outside the Error
term, see Examples).
Suppressing the intercept
(i.e, via 0 +
or  1
) is ignored. Specific specifications of
effects (e.g., excluding terms with 
or using ^
) could be okay
but is not tested. Using the I
or poly
function
within the formula is not tested and not supported!
To run an ANCOVA you need to set factorize = FALSE
and make sure that
all variables have the correct type (i.e., factors are factors and numeric
variables are numeric and centered).
Note that the default behavior is to include calculation of the effect size
generalized etasquared for which all nonmanipluated (i.e.,
observed) variables need to be specified via the observed
argument to
obtain correct results. When changing the effect size to "pes"
(partial etasquared) or "none"
via anova_table
this becomes
unnecessary.
Factor contrasts will be set to "contr.sum"
for all betweensubject
factors if default contrasts are not equal to "contr.sum"
or
attrib(factor, "contrasts") != "contr.sum"
. (withinsubject factors
are hardcoded "contr.sum"
.)
Type 3 sums of squares are default in afex. While some authors argue that socalled type 3 sums of squares are dangerous and/or problematic (most notably Venables, 2000), they are the default in many commercial statistical application such as SPSS or SAS. Furthermore, statisticians with an applied perspective recommend type 3 tests (e.g., Maxwell and Delaney, 2004). Consequently, they are the default for the ANOVA functions described here. For some more discussion on this issue see here.
Note that lower order effects (e.g., main effects) in type 3 ANOVAs are only
meaningful with
effects
coding. Therefore, contrasts are set to contr.sum
which
ensures meaningful results. For a discussion of the other (nonrecommended)
coding schemes see
here.
The S3 object returned
per default can be directly passed to emmeans::emmeans
for further
analysis. This allows to test any type of contrasts that might be of interest
independent of whether or not this contrast involves betweensubject
variables, withinsubject variables, or a combination thereof. The general
procedure to run those contrasts is the following (see Examples for a full
example):
Estimate an afex_aov
object with the function returned here. For example: x < aov_car(dv ~ a*b + (id/c), d)
Obtain a emmGridclass
object by running emmeans
on the afex_aov
object from step 1 using the factors involved in the contrast. For example: r < emmeans(x, ~a:c)
Create a list containing the desired contrasts on the reference grid object from step 2. For example: con1 < list(a_x = c(1, 1, 0, 0, 0, 0), b_x = c(0, 0, 0.5, 0.5, 0, 1))
Test the contrast on the reference grid using contrast
. For example: contrast(r, con1)
To control for multiple testing pvalue adjustments can be specified. For example the BonferroniHolm correction: contrast(r, con1, adjust = "holm")
Note that emmeans allows for a variety of advanced settings and
simplifications, for example: all pairwise comparison of a single factor
using one command (e.g., emmeans(x, "a", contr = "pairwise")
) or
advanced control for multiple testing by passing objects to multcomp.
A comprehensive overview of the functionality is provided in the
accompanying vignettes (see
here).
Since version 1.0, afex per default uses the multivariate
model
(i.e., the lm
slot of the afex_aov
object) for followup tests
with emmeans. Compared to the univariate
model (i.e., the
aov
slot), this can handle unbalanced data and addresses sphericity
better. To use the older (and not recommended) model = "univariate"
make sure to set include_aov = TRUE
when estimating the ANOVA.
Starting with afex version 0.22, emmeans is not
loaded/attached automatically when loading afex. Therefore,
emmeans now needs to be loaded by the user via
library("emmeans")
or require("emmeans")
.
afex_aov
Objects A full overview over the
methods provided for afex_aov
objects is provided in the corresponding
help page: afex_aovmethods
. The probably most important ones
for endusers are summary
, anova
, and nice
.
The summary
method returns, for ANOVAs containing withinsubject
(repeatedmeasures) factors with more than two levels, the complete
univariate analysis: Results without dfcorrection, the GreenhouseGeisser
corrected results, the HyunhFeldt corrected results, and the results of the
Mauchly test for sphericity.
The anova
method returns a data.frame
of class "anova"
containing the ANOVA table in numeric form (i.e., the one in slot
anova_table
of a afex_aov
). This method has arguments such as
correction
and es
and can be used to obtain an ANOVA table with
different correction than the one initially specified.
The nice
method also returns a data.frame
, but rounds
most values and transforms them into characters for nice printing. Also has
arguments like correction
and es
which can be used to obtain an
ANOVA table with different correction than the one initially specified.
aov_car
, aov_4
, and aov_ez
are wrappers for
Anova
and aov
, the return value is
dependent on the return
argument. Per default, an S3 object of class
"afex_aov"
is returned containing the following slots:
"anova_table"
An ANOVA table of class c("anova",
"data.frame")
.
"aov"
aov
object returned from aov
(should not be used to evaluate significance of effects, but can be passed
to emmeans
for posthoc tests).
"Anova"
object returned from Anova
, an
object of class "Anova.mlm"
(if withinsubjects factors are present)
or of class c("anova", "data.frame")
.
"lm"
the object fitted with lm
and passed to
Anova
(i.e., an object of class "lm"
or "mlm"
). Also
returned if return = "lm"
.
"data"
a list containing: (1) long
(the possibly
aggregated data in long format used for aov
), wide
(the data
used to fit the lm
object), and idata
(if withinsubject
factors are present, the idata
argument passed to
car::Anova
). Also returned if return = "data"
.
In addition, the object has the following attributes: "dv"
,
"id"
, "within"
, "between"
, and "type"
.
The print method for afex_aov
objects
(invisibly) returns (and prints) the same as if return
is
"nice"
: a nice ANOVA table (produced by nice
) with the
following columns: Effect
, df
, MSE
(meansquared
errors), F
(potentially with significant symbols), ges
(generalized etasquared), p
.
aov_4()
: Allows definition of ANOVAmodel using
lme4::lmer
like Syntax (but still fits a standard ANOVA).
aov_ez()
: Allows definition of ANOVAmodel using character strings.
Calculation of ANOVA models via aov
(which is done per default)
can be comparatively slow and produce comparatively large objects for
ANOVAs with many withinsubjects factors or levels. To avoid this
calculation set include_aov = FALSE
. You can also disable this
globally with: afex_options(include_aov = FALSE)
The id variable and variables entered as withinsubjects (i.e.,
repeatedmeasures) factors are silently converted to factors. Levels of
withinsubject factors are converted to valid variable names using
make.names(...,unique=TRUE)
. Unused factor levels are
silently dropped on all variables.
Contrasts attached to a factor as an attribute are probably not preserved and not supported.
The workhorse is aov_car
. aov_4
and aov_ez
only
construe and pass an appropriate formula to aov_car
. Use
print.formula = TRUE
to view this formula.
In contrast to aov
aov_car
assumes that all factors to
the right of /
in the Error
term are belonging together.
Consequently, Error(id/(a*b))
and Error(id/a*b)
are identical
(which is not true for aov
).
Henrik Singmann
The design of these functions was influenced by ezANOVA
from package ez.
Cramer, A. O. J., van Ravenzwaaij, D., Matzke, D., Steingroever, H., Wetzels, R., Grasman, R. P. P. P., ... Wagenmakers, E.J. (2015). Hidden multiplicity in exploratory multiway ANOVA: Prevalence and remedies. Psychonomic Bulletin & Review, 18. doi:10.3758/s1342301509135
Maxwell, S. E., & Delaney, H. D. (2004). Designing Experiments and Analyzing Data: A ModelComparisons Perspective. Mahwah, N.J.: Lawrence Erlbaum Associates.
Venables, W.N. (2000). Exegeses on linear models. Paper presented to the SPlus User's Conference, Washington DC, 89 October 1998, Washington, DC. Available from: http://www.stats.ox.ac.uk/pub/MASS3/Exegeses.pdf
Various methods for objects of class afex_aov
are available:
afex_aovmethods
nice
creates the nice ANOVA tables which is by default printed.
See also there for a slightly longer discussion of the available effect
sizes.
mixed
provides a (formula) interface for obtaining pvalues for
mixedmodels via lme4. The functions presented here do not estimate
mixed models.
##########################
## 1: Specifying ANOVAs ##
##########################
# Example using a purely withinsubjects design
# (Maxwell & Delaney, 2004, Chapter 12, Table 12.5, p. 578):
data(md_12.1)
aov_ez("id", "rt", md_12.1, within = c("angle", "noise"),
anova_table=list(correction = "none", es = "none"))
# Default output
aov_ez("id", "rt", md_12.1, within = c("angle", "noise"))
# examples using obk.long (see ?obk.long), a long version of the OBrienKaiser
# dataset (car package). Data is a splitplot or mixed design: contains both
# within and betweensubjects factors.
data(obk.long, package = "afex")
# estimate mixed ANOVA on the full design:
aov_car(value ~ treatment * gender + Error(id/(phase*hour)),
data = obk.long, observed = "gender")
aov_4(value ~ treatment * gender + (phase*hourid),
data = obk.long, observed = "gender")
aov_ez("id", "value", obk.long, between = c("treatment", "gender"),
within = c("phase", "hour"), observed = "gender")
# the three calls return the same ANOVA table:
# Anova Table (Type 3 tests)
#
# Response: value
# Effect df MSE F ges p.value
# 1 treatment 2, 10 22.81 3.94 + .198 .055
# 2 gender 1, 10 22.81 3.66 + .115 .085
# 3 treatment:gender 2, 10 22.81 2.86 .179 .104
# 4 phase 1.60, 15.99 5.02 16.13 *** .151 <.001
# 5 treatment:phase 3.20, 15.99 5.02 4.85 * .097 .013
# 6 gender:phase 1.60, 15.99 5.02 0.28 .003 .709
# 7 treatment:gender:phase 3.20, 15.99 5.02 0.64 .014 .612
# 8 hour 1.84, 18.41 3.39 16.69 *** .125 <.001
# 9 treatment:hour 3.68, 18.41 3.39 0.09 .002 .979
# 10 gender:hour 1.84, 18.41 3.39 0.45 .004 .628
# 11 treatment:gender:hour 3.68, 18.41 3.39 0.62 .011 .641
# 12 phase:hour 3.60, 35.96 2.67 1.18 .015 .335
# 13 treatment:phase:hour 7.19, 35.96 2.67 0.35 .009 .930
# 14 gender:phase:hour 3.60, 35.96 2.67 0.93 .012 .449
# 15 treatment:gender:phase:hour 7.19, 35.96 2.67 0.74 .019 .646
# 
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1
#
# Sphericity correction method: GG
# "numeric" variables are per default converted to factors (as long as factorize
# = TRUE):
obk.long$hour2 < as.numeric(as.character(obk.long$hour))
# gives same results as calls before
aov_car(value ~ treatment * gender + Error(id/phase*hour2),
data = obk.long, observed = c("gender"))
# ANCOVA: adding a covariate (necessary to set factorize = FALSE)
aov_car(value ~ treatment * gender + age + Error(id/(phase*hour)),
data = obk.long, observed = c("gender", "age"), factorize = FALSE)
aov_4(value ~ treatment * gender + age + (phase*hourid),
data = obk.long, observed = c("gender", "age"), factorize = FALSE)
aov_ez("id", "value", obk.long, between = c("treatment", "gender"),
within = c("phase", "hour"), covariate = "age",
observed = c("gender", "age"), factorize = FALSE)
# aggregating over one withinsubjects factor (phase), with warning:
aov_car(value ~ treatment * gender + Error(id/hour), data = obk.long,
observed = "gender")
aov_ez("id", "value", obk.long, c("treatment", "gender"), "hour",
observed = "gender")
# aggregating over both withinsubjects factors (again with warning),
# only betweensubjects factors:
aov_car(value ~ treatment * gender + Error(id), data = obk.long,
observed = c("gender"))
aov_4(value ~ treatment * gender + (1id), data = obk.long,
observed = c("gender"))
aov_ez("id", "value", obk.long, between = c("treatment", "gender"),
observed = "gender")
# only withinsubject factors (ignoring betweensubjects factors)
aov_car(value ~ Error(id/(phase*hour)), data = obk.long)
aov_4(value ~ (phase*hourid), data = obk.long)
aov_ez("id", "value", obk.long, within = c("phase", "hour"))
### changing defaults of ANOVA table:
# no dfcorrection & partial etasquared:
aov_car(value ~ treatment * gender + Error(id/(phase*hour)),
data = obk.long, anova_table = list(correction = "none", es = "pes"))
# no dfcorrection and no MSE
aov_car(value ~ treatment * gender + Error(id/(phase*hour)),
data = obk.long,observed = "gender",
anova_table = list(correction = "none", MSE = FALSE))
# add pvalue adjustment for all effects (see Cramer et al., 2015, PB&R)
aov_ez("id", "value", obk.long, between = "treatment",
within = c("phase", "hour"),
anova_table = list(p_adjust_method = "holm"))
###########################
## 2: Followup Analysis ##
###########################
# use data as above
data(obk.long, package = "afex")
# 1. obtain afex_aov object:
a1 < aov_ez("id", "value", obk.long, between = c("treatment", "gender"),
within = c("phase", "hour"), observed = "gender")
if (requireNamespace("ggplot2") & requireNamespace("emmeans")) {
# 1b. plot data using afex_plot function, for more see:
## vignette("afex_plot_introduction", package = "afex")
## default plot uses multivariate modelbased CIs
afex_plot(a1, "hour", "gender", c("treatment", "phase"))
a1b < aov_ez("id", "value", obk.long, between = c("treatment", "gender"),
within = c("phase", "hour"), observed = "gender",
include_aov = TRUE)
## you can use a univariate model and CIs if you refit the model with the aov
## slot
afex_plot(a1b, "hour", "gender", c("treatment", "phase"),
emmeans_arg = list(model = "univariate"))
## in a mixed betweenwithin designs, no errorbars might be preferrable:
afex_plot(a1, "hour", "gender", c("treatment", "phase"), error = "none")
}
if (requireNamespace("emmeans")) {
library("emmeans") # package emmeans needs to be attached for followup tests.
# 2. obtain reference grid object (default uses multivariate model):
r1 < emmeans(a1, ~treatment +phase)
r1
# 3. create list of contrasts on the reference grid:
c1 < list(
A_B_pre = c(rep(0, 6), 0, 1, 1), # A versus B for pretest
A_B_comb = c(0.5, 0.5, 0, 0.5, 0.5, 0, 0, 0, 0), # A vs. B for post and followup combined
effect_post = c(0, 0, 0, 1, 0.5, 0.5, 0, 0, 0), # control versus A&B post
effect_fup = c(1, 0.5, 0.5, 0, 0, 0, 0, 0, 0), # control versus A&B followup
effect_comb = c(0.5, 0.25, 0.25, 0.5, 0.25, 0.25, 0, 0, 0) # control versus A&B combined
)
# 4. test contrasts on reference grid:
contrast(r1, c1)
# same as before, but using BonferroniHolm correction for multiple testing:
contrast(r1, c1, adjust = "holm")
# 2. (alternative): all pairwise comparisons of treatment:
emmeans(a1, "treatment", contr = "pairwise")
}
#######################
## 3: Other examples ##
#######################
data(obk.long, package = "afex")
# replicating ?Anova using aov_car:
obk_anova < aov_car(value ~ treatment * gender + Error(id/(phase*hour)),
data = obk.long, type = 2)
# in contrast to aov you do not need the withinsubject factors outside Error()
str(obk_anova, 1, give.attr = FALSE)
# List of 5
# $ anova_table:Classes ‘anova’ and 'data.frame': 15 obs. of 6 variables:
# $ aov : NULL
# $ Anova :List of 14
# $ lm :List of 13
# $ data :List of 3
obk_anova$Anova
# Type II Repeated Measures MANOVA Tests: Pillai test statistic
# Df test stat approx F num Df den Df Pr(>F)
# (Intercept) 1 0.96954 318.34 1 10 6.532e09 ***
# treatment 2 0.48092 4.63 2 10 0.0376868 *
# gender 1 0.20356 2.56 1 10 0.1409735
# treatment:gender 2 0.36350 2.86 2 10 0.1044692
# phase 1 0.85052 25.61 2 9 0.0001930 ***
# treatment:phase 2 0.68518 2.61 4 20 0.0667354 .
# gender:phase 1 0.04314 0.20 2 9 0.8199968
# treatment:gender:phase 2 0.31060 0.92 4 20 0.4721498
# hour 1 0.93468 25.04 4 7 0.0003043 ***
# treatment:hour 2 0.30144 0.35 8 16 0.9295212
# gender:hour 1 0.29274 0.72 4 7 0.6023742
# treatment:gender:hour 2 0.57022 0.80 8 16 0.6131884
# phase:hour 1 0.54958 0.46 8 3 0.8324517
# treatment:phase:hour 2 0.66367 0.25 16 8 0.9914415
# gender:phase:hour 1 0.69505 0.85 8 3 0.6202076
# treatment:gender:phase:hour 2 0.79277 0.33 16 8 0.9723693
# 
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1