ANOVA Exercises in R: 15 One-Way and Two-Way Practice Problems

Explanation: aov() is a thin wrapper around lm() whose summary.aov method prints the classical ANOVA table instead of a regression coefficient table. The F-statistic is the ratio of between-group mean square to within-group mean square; under the null of equal group means it follows an F distribution with (k-1, N-k) degrees of freedom. A p-value of 0.016 rejects equal means at the 5% level, but does not tell you which pair differs (that needs a post-hoc test like Tukey HSD).

Explanation: The F value of 119 is extremely large because between-species variance dominates within-species variance for sepal length: setosa is much shorter than virginica on average. With three species the numerator df is 2 and the denominator df is 150 minus 3 = 147. A p-value reported as <2e-16 is at the machine-precision floor; you can pull the exact number with summary(ex_1_2)[[1]][["Pr(>F)"]][1]. Always pair an "obviously significant" F test with effect size (eta-squared) and a plot, because significance alone says nothing about practical magnitude.

Explanation: The single-factor model lumps all dose levels together, which inflates the residual sum of squares because dose is a real (and large) source of variance left unmodelled. That is why supplement looks only marginally significant here (p = 0.060) but jumps to p = 0.0004 once you add dose as a second factor in Exercise 3.1. Treat a borderline p-value from a one-way ANOVA as a hint that an important covariate may be missing from the right-hand side.

Explanation: summary.aov returns a list whose first element is a data frame with columns Df, Sum Sq, Mean Sq, F value, Pr(>F). The first row is the factor effect and the second row is the residual (which has no F or p). Indexing with [[1]] and then a column name gives a numeric vector you can subset. For programmatic pipelines, broom::tidy(ex_1_1) returns the same numbers as a tibble with stable column names and is friendlier than positional indexing.

Explanation: Tukey HSD inverts the studentized range distribution to control the family-wise error rate at 5%, so the three pairwise p-values are already adjusted (you do not apply Bonferroni on top). The only significant pair is trt2 vs trt1; both treatments versus control are non-significant. This is a common pattern with three-arm trials: the two treatments together "move the needle" enough to flag the omnibus F test, but only the most extreme pair survives a controlled pairwise comparison.

Explanation: Bonferroni multiplies each raw p-value by the number of comparisons (here 3) and caps at 1. It is the most conservative of the common corrections, which is why the trt2 vs trt1 adjusted p of 0.013 is slightly larger than Tukey's 0.012. With more comparisons (six groups: 15 pairs) Bonferroni quickly becomes underpowered; Tukey HSD or Holm are usually preferable for many comparisons. Use the pooled-SD variant when you trust the equal-variance assumption from Levene's test; otherwise pass pool.sd = FALSE for a Welch-style version.

Explanation: Once dose is added to the right-hand side, the residual mean square drops from 56 (one-way) to 14.7 because most of that variance was really dose, not noise. Supplement now has an F of 14 and a p of 0.0004 instead of being marginal: a clean illustration of why omitting a strong covariate inflates within-group variance and hides real effects. Wrap dose in factor() so it is treated as three discrete levels (0.5, 1, 2) with two degrees of freedom; leaving it numeric would force a linear trend with only one degree of freedom.

Explanation: The interaction row is significant (p = 0.022), so the supplement effect is NOT the same at every dose. From the cell means you will see OJ outperforms VC at 0.5 and 1 mg, but the two are nearly tied at 2 mg. When the interaction is significant, you should interpret main effects with care or pivot to simple-effect tests within each dose level rather than reporting one aggregate "supp effect". * expands to all main effects plus the interaction; (supp + factor(dose))^2 would give the same model for two factors.

Explanation: The OJ-minus-VC gap is +5.25 at 0.5 mg, +5.93 at 1 mg, but only -0.08 at 2 mg. That collapsing gap is exactly the interaction the ANOVA flagged. Cell means tables are the most important diagnostic before reading any interaction p-value, because a "significant interaction" can be ordinal (lines stay in the same rank order) or disordinal (lines cross), and the interpretation differs. For tidyverse users, the equivalent is ToothGrowth |> group_by(supp, dose) |> summarise(mean = mean(len)).

Explanation: Levene's test runs a one-way ANOVA on the absolute deviations from each group's center (median by default, which is the robust Brown-Forsythe variant). A non-significant result (p = 0.34) means there is not enough evidence to reject equal variances; classical aov() is fine. If Levene's were significant you would either transform weight (often log()), fit a Welch-style ANOVA via oneway.test(..., var.equal = FALSE), or move to a non-parametric test like Kruskal-Wallis.

Explanation: Test the model's residuals, not the raw response, because ANOVA assumes Normal(0, sigma^2) errors conditional on the predictors. With n = 30 the Shapiro test has modest power, so a non-significant p (0.44) is weak evidence of normality, not proof; always pair it with a QQ-plot plot(ex_1_1, which = 2). ANOVA is also robust to mild non-normality when sample sizes are equal across groups, so a borderline Shapiro result on balanced data rarely overturns the substantive conclusion.

Explanation: Welch's ANOVA is the F-test analogue of Welch's t-test: it weights each group by its own variance and shrinks the denominator degrees of freedom so the test stays calibrated when variances differ. Notice the denom df here is a non-integer 17.13, smaller than the classical 27, which costs some power but buys robustness. For pairwise follow-ups under unequal variances, use pairwise.t.test(..., pool.sd = FALSE); the Games-Howell post-hoc (from package rstatix or userfriendlyscience) is the Welch-style analogue of Tukey HSD.

Explanation: Eta-squared is the proportion of total variance in the outcome attributable to the factor, the ANOVA analogue of R-squared. Cohen's rough benchmarks call 0.01 small, 0.06 medium, and 0.14 large; the PlantGrowth effect of 0.26 is large. Eta-squared is upward-biased for small samples; partial eta-squared and omega-squared correct for that and are reported by effectsize::eta_squared() and effectsize::omega_squared(). Always report an effect size next to the F value because significance scales with sample size while effect size does not.

Explanation: Assigning a contrast matrix to a factor changes what the lm coefficients estimate without refitting any new data. The two columns must be orthogonal and sum to zero within each column; c(2, -1, -1) and c(0, 1, -1) qualify. group1 is one-third of the ctrl-versus-mean-of-treatments contrast (the +2/-1/-1 coding scales it), and its non-significant p (0.71) says the two treatments together do NOT differ from control on average. group2 tests trt1 versus trt2 and is significant (p = 0.035), matching the Tukey HSD finding from Exercise 2.2. Custom contrasts let you test one targeted scientific question per degree of freedom instead of all pairs.

Explanation: All three group-difference tests (classical, Welch, Kruskal-Wallis) agree on the conclusion at the 5% level, which is the reassuring case: the decision is robust to choice of test. When tests disagree, prefer the one whose assumptions are satisfied: classical F if Levene and Shapiro both pass, Welch's F if Levene flags unequal variances, Kruskal-Wallis if Shapiro flags clearly non-normal residuals (especially with small or unbalanced groups). Reporting all three is overkill for a paper, but during analysis it is a useful sensitivity check.