# A tibble: 4 × 2
crab_type n
<fct> <int>
1 female 21
2 intact male 21
3 male major removed 21
4 male minor removed 21
Single-factor ANOVA
Covered last time — Multiple Regression:
fct_relevel()Note
✅ Key idea from last lecture
You already built an ANOVA table for a regression, and you already turned a categorical predictor into dummy variables.
Today every predictor is categorical. The machinery does not change at all — only the name on the front does.
One-way ANOVA, end to end:
Note
📖 Reading
Whitlock & Schluter, Ch. 15 — Comparing Means of More Than Two Groups.
Gotelli & Ellison, Ch. 10 — The Analysis of Variance.
| Continuous X | Categorical X | |
|---|---|---|
| Continuous Y | Regression | ANOVA |
| Categorical Y | Logistic | Contingency |
A male fiddler crab has one hugely enlarged claw. Darnell & Munguia asked whether it doubles as a heat radiator, dumping heat so the crab can stay out foraging in the sun.
Four groups of 21 crabs each, body temperature above ambient measured in the field:
intact male — the controlmale major removed — the big claw taken offmale minor removed — a small claw taken off, the procedural controlfemale — never had a major claw at allIf the big claw sheds heat, removing it should make males hotter.
# A tibble: 4 × 2
crab_type n
<fct> <int>
1 female 21
2 intact male 21
3 male major removed 21
4 male minor removed 21
Note
A balanced design
Equal n in every group — which makes the degrees-of-freedom arithmetic below clean, and the test robust.
Warning
🗳️ Six times today, you make the call
Six slides are marked 🗳️ DECIDE. Each asks a multiple-choice question, and I will not tell you the answer on that slide.
How it works: write your letter down — actually write it — then hold up fingers 1–4. If the room splits, turn to the person next to you, argue for ninety seconds, and we vote again. Nothing here is graded. The point is that you commit before you find out.
Regression
\[Y_i = \beta_0 + \beta_1X_i + \varepsilon_i\]
ANOVA
\[Y_{ij} = \mu + A_i + \varepsilon_{ij}\]
Both:
lm() cannot multiply by the word “female.” Behind the formula it builds a dummy variable for every level except the first — a column of 0s and 1s saying “is this crab in that group?”
(Intercept) crab_typeintact male crab_typemale major removed
1 1 0 0
2 1 0 0
3 1 0 0
crab_typemale minor removed
1 0
2 0
3 0
Four groups, three dummy columns. The level with no column of its own is the reference, and it is carried by the intercept.
Note
✅ This is the whole trick
A one-way ANOVA is a multiple regression whose predictors happen to be 0/1 columns. You met those columns last week under the name dummy variables.
aov() and lm() fit the identical model — aov() just prints it as an ANOVA table by default.
# A tibble: 4 × 3
crab_type group_mean from_coefficients
<fct> <dbl> <dbl>
1 female 1.68 1.68
2 intact male 1.29 1.29
3 male major removed 1.41 1.41
4 male minor removed 1.21 1.21
The intercept is the mean of the reference group. Every other coefficient is that group’s difference from it, so adding them back recovers each group mean exactly.
Important
✅ Why this matters in a moment
If the coefficients are “differences from the reference,” then which level is the reference decides what every coefficient means.
R chose it alphabetically. Part 7 is about taking that choice back.
With four crab groups there are six possible pairwise t-tests. Suppose every group really has the same mean, and you run all six at α = 0.05.
What is the chance you get at least one “significant” result?
(a) 5% — that is what α = 0.05 means
(b) About 12%
(c) About 26%
(d) 30% — 6 × 5%
Warning
✍️ Write your letter down, then vote.
Note
The question to ask
What is the chance of getting no false positives in six independent tries?
Important
✅ The answer is (c) — about 26%.
\[1 - (1 - 0.05)^6 = 0.265\]
Each test has a 95% chance of behaving. Six in a row behave only 73.5% of the time, so you have roughly a one in four chance of at least one false positive — before you have learned anything.
(d) is the common guess. Adding probabilities overshoots, because the ways of being wrong overlap.
ANOVA’s job is to ask one question first: is there any difference at all among these groups?
One test, one α, one p-value. Only if it rejects do you go looking for which groups differ — and then you use a post-hoc method that controls the error rate across that whole family of comparisons.
With 10 groups there are 45 pairs, and the chance of at least one false positive is over 90%.
Even if every group had the same true mean, the sample means would differ, just from sampling error.
ANOVA asks whether the sample means differ more than sampling error alone would produce:
\[F = \frac{MS_{among}}{MS_{within}}\]

Every point’s distance from the grand mean splits into two pieces, exactly as in regression:
\[\underbrace{(Y_{ij} - \bar{Y})}_{\text{total}} = \underbrace{(\bar{Y}_i - \bar{Y})}_{\text{among}} + \underbrace{(Y_{ij} - \bar{Y}_i)}_{\text{within}}\]
Square and sum over every observation:
\[SS_{total} = SS_{among} + SS_{within}\]
The dashed line is the grand mean; the coloured bars are the group means.

\[SS_{total} = \sum_{i}\sum_{j}(Y_{ij} - \bar{Y})^2\]
\[SS_{among} = \sum_{i} n_i(\bar{Y}_i - \bar{Y})^2\]
\[SS_{within} = \sum_{i}\sum_{j}(Y_{ij} - \bar{Y}_i)^2\]
with a groups and \(n_i\) observations in group i:
Tip
🖐 Note the \(n_i\)
Each group mean’s deviation is weighted by that group’s own sample size. With a balanced design every \(n_i\) is the same and it simplifies to \(n\sum(\bar{Y}_i - \bar{Y})^2\) — but only then.
Our design is balanced at n = 21, so both forms agree here.
grand_mean <- mean(crab_df$body_temp_c)
crab_parts <- crab_df %>%
group_by(crab_type) %>%
mutate(group_mean = mean(body_temp_c)) %>%
ungroup()
ss_total <- sum((crab_parts$body_temp_c - grand_mean)^2)
ss_among <- sum((crab_parts$group_mean - grand_mean)^2)
ss_within <- sum((crab_parts$body_temp_c - crab_parts$group_mean)^2)
c(total = ss_total, among = ss_among, within = ss_within) total among within
6.108929 2.641310 3.467619
[1] 6.108929
2.641 + 3.468 = 6.109 — exactly, as it must.
Important
✅ The same identity as regression
Different names, same partition, same arithmetic.
Every number is one you just built:
F = 20.31, p = 7.00e-10.
The ANOVA returns F = 20.31, p = 7.00e-10.
What have you just learned?
(a) All four crab types have different mean temperatures
(b) At least one crab type differs from at least one other
(c) Removing the major claw raises body temperature
(d) The major claw explains 43.2% of the variation in temperature, so the radiator hypothesis is supported
Warning
✍️ One letter. Write it down.
Note
The question to ask
What was the null hypothesis, exactly?
A test can only tell you about the hypothesis it tested.
Important
✅ The answer is (b).
H₀ was “all four means are equal.” Rejecting it says at least one differs — not which, not how many, and not in which direction.
(a) claims far more than the test supports. (c) is a directional claim about one specific pair; that needs a post-hoc comparison, and Part 8 shows it is not even true. (d) misreads R²: 0.432 is the variance explained by crab type, which bundles sex and two surgeries together, not by the claw alone.
Warning
⚠️ The most common ANOVA mistake
Reporting a significant F and then describing whichever group means look farthest apart, as though the F test had said so.
It did not. The F test is an entry ticket, not an answer.
\[F = \frac{MS_{among}}{MS_{within}}\]
Here F = 20.31 — the group means are 20 times further apart than within-group noise alone can explain.
Tip
🖐 Why F is never negative
It is a ratio of two variances, and a variance is a sum of squares. The direction of a difference disappears the moment you square it — which is why a significant F cannot tell you which way anything went.
Take just two of the four groups — intact males and females — and test them both ways:
t = 5.958, so t² = 35.5 — and F = 35.5.
Identical, with identical p-values.
Important
✅ ANOVA generalises the t-test
With two groups they are the same test. With three or more, the t-test has nowhere to go and ANOVA does.
This is the same t² = F identity you proved on a regression slope in Week 6.
\[R^2 = \frac{SS_{among}}{SS_{total}} = \frac{2.64}{6.11} = 0.432\]
Crab type accounts for about 43.2% of the variation in body temperature. The rest is individual crabs being individual crabs.
Note
Significant and partial, at the same time
p = 7.00e-10 says the group differences are not noise. R² = 0.43 says they explain under half of what is going on.
Both are true, and a results section should carry both.
One-way ANOVA assumes:
ANOVA is fairly robust to mild departures from 2 and 3 when the design is balanced, which is one more reason to plan for equal n.
Warning
⚠️ Normality of what?
Not body_temp_c. Not each group separately. The residuals of the fitted model — exactly as in regression.
A response that looks bimodal because the groups differ is perfectly fine; that is the effect you are testing.
The same four panels you read for a regression:

Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 3 0.1587 0.9238
80
H₀ is “the variances are equal”, so a large p-value is the good outcome — the same backwards logic as Shapiro-Wilk.
p = 0.924 — nothing to worry about.
Tip
🖐 Levene vs Bartlett
Bartlett’s test is more powerful but assumes normality itself. Levene’s works on absolute deviations from the group median, so it is far more robust. Use Levene.
Shapiro-Wilk normality test
data: residuals(crab_model)
W = 0.97566, p-value = 0.112
Note what goes inside: residuals(crab_model).
p = 0.112 — the residuals are consistent with normal.
Warning
⚠️ Do not let the test overrule the plot
At small n it cannot detect real problems; at large n it flags departures too small to matter. The QQ plot is the primary evidence.
Read the file without thinking about it and the four crab types become a factor in alphabetical order:
So the reference is female — and every coefficient in the model becomes “compared with a female crab.”
Warning
⚠️ That is the wrong baseline for this experiment
The question is what happens when you remove a claw from a male. The comparison that answers it is against an intact male, not against a female.
The model fits identically either way. The output only answers your question if you tell it what your question was.
You are about to set the reference level deliberately.
Which level should it be, and why?
(a) female — keep R’s default; the model fit is identical anyway
(b) intact male — the untouched version of the animals that were manipulated
(c) male minor removed — it is the procedural control, so it accounts for the surgery
(d) It makes no difference, because the F-test and R² do not change
Warning
✍️ Careful — two of these are true statements. Only one is the right decision.
Note
The question to ask
Which single comparison do I most want to read straight off the coefficient table?
fct_relevel() · Saying What You MeantImportant
✅ The answer is (b).
is true — F, R², and the ANOVA table are identical whichever level is the reference. That is not a reason to accept a default chosen by spelling.
is tempting and it is the right comparison to make eventually — but as the baseline it hides the very contrast the experiment was built around.
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.286 0.045 28.300 0.000
crab_typefemale 0.390 0.064 6.077 0.000
crab_typemale major removed 0.129 0.064 2.001 0.049
crab_typemale minor removed -0.076 0.064 -1.186 0.239
Now every row reads as “degrees hotter than an intact male”:
female +0.39 — clearly hottermale major removed +0.13male minor removed -0.08Tip
🖐 The habit
Set the reference in the same line that creates the factor. fct_relevel(x, "control") moves one level to the front and leaves the rest alone.
The F-test said something differs. Post-hoc tests say what, while keeping the family-wide error rate at α.
Tukey’s HSD is the standard for all pairwise comparisons with equal n.
Warning
⚠️ Deciding after you look is not a plan
Choosing which comparisons to make because they looked big re-introduces exactly the multiple-comparisons problem the correction exists to solve.
Tukey multiple comparisons of means
95% family-wise confidence level
Fit: aov(formula = crab_model)
$crab_type
diff lwr upr
female-intact male 0.39047619 0.22189182 0.55906056
male major removed-intact male 0.12857143 -0.04001294 0.29715579
male minor removed-intact male -0.07619048 -0.24477484 0.09239389
male major removed-female -0.26190476 -0.43048913 -0.09332040
male minor removed-female -0.46666667 -0.63525103 -0.29808230
male minor removed-male major removed -0.20476190 -0.37334627 -0.03617754
p adj
female-intact male 0.0000002
male major removed-intact male 0.1962212
male minor removed-intact male 0.6375171
male major removed-female 0.0006106
male minor removed-female 0.0000000
male minor removed-male major removed 0.0108312
Each row is one pair: the difference in means, a confidence interval, and a p-value already adjusted for the whole family of six comparisons.
A pair differs when its interval excludes zero.
The radiator hypothesis predicts that a male whose major claw was removed should run hotter than an intact male.
Tukey gives male major removed − intact male a difference of +0.129 with an adjusted p = 0.196.
What should you conclude?
(a) The hypothesis is supported — the difference is in the predicted direction
(b) The hypothesis fails — removing the major claw did not significantly change temperature
(c) The hypothesis is supported, because females are significantly hotter than every male group
(d) Nothing can be concluded, because the overall ANOVA was significant
Warning
✍️ One letter, then argue it with your neighbour.
Note
The question to ask
Which comparison does the hypothesis actually make a prediction about?
Important
✅ The answer is (b).
Removing the major claw moved temperature by 0.13 °C, adjusted p = 0.2 — not distinguishable from zero. The direction is right and the evidence is not there.
(c) is the trap: females are hotter than every male group, but females differ from males in many ways besides the claw. That comparison cannot test a claw-removal hypothesis.
What did differ (adjusted p):
| Comparison | p |
|---|---|
| female − intact male | 0 |
| major removed − intact | 0.196 |
| minor removed − intact | 0.638 |
| major − minor removed | 0.011 |
The big effect in this dataset is sex, not surgery.
emmeans gives the same comparisons with more control — any contrast you like, any adjustment method:
contrast estimate SE df t.ratio p.value
intact male - female -0.3905 0.0643 80 -6.077 <0.0001
intact male - male major removed -0.1286 0.0643 80 -2.001 0.1962
intact male - male minor removed 0.0762 0.0643 80 1.186 0.6375
female - male major removed 0.2619 0.0643 80 4.076 0.0006
female - male minor removed 0.4667 0.0643 80 7.263 <0.0001
male major removed - male minor removed 0.2048 0.0643 80 3.187 0.0108
P value adjustment: tukey method for comparing a family of 4 estimates
Tip
🖐 Why emmeans and not just TukeyHSD()
TukeyHSD() only does all-pairs on a balanced aov. emmeans handles unbalanced designs, covariates, several factors and planned contrasts — and it is what you will use from Week 9 onward.
A results sentence needs the test, the statistic with both df, the p-value, and the direction — then the post-hoc detail.
Body temperature differed among crab types (one-way ANOVA, F3,80 = 20.31, p < 0.001, R² = 0.43). Females were warmer than intact males by 0.39 °C (Tukey HSD, adjusted p < 0.001), but removing the major claw did not detectably change male temperature (0.13 °C, adjusted p = 0.2).
Note
✅ What that sentence does
The crab data behaved. Most data do not.
Here are 200 nematodes, 50 in each of four treatments, exposed to a parasite as larvae, as adults, as both, or to water as a control. The response is lifespan in days.

Note that the water control sorts last alphabetically — another reason to relevel on purpose.
Shapiro-Wilk normality test
data: residuals(worm_model)
W = 0.97626, p-value = 0.001808
Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 3 22.481 1.493e-12 ***
196
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Look back at the boxplots: group SDs run from 3.3 to 8.7 days.
Warning
⚠️ Nothing here is simulated
This is the published dataset, unmodified. Real data fails these tests routinely.
An ANOVA on these data returns F = 15.25, p = 5.87e-09 — but the residuals are non-normal and the variances are unequal.
What is the right next step?
(a) Report the ANOVA — it is significant, and ANOVA is robust
(b) Drop the group with the biggest variance and re-run
(c) Try a transformation, and if that fails use a rank-based test such as Kruskal-Wallis
(d) Report nothing — the data cannot be analysed
Warning
✍️ Write your letter down, then vote.
Note
The question to ask
Is the problem with the data, or with the test I chose for it?
Important
✅ The answer is (c).
Transform first — a log or square root often fixes both problems at once, and keeps you in a framework that extends to more factors. If it does not, move to a rank-based test.
(b) is data deletion to make a test pass, which is not a repair.
Kruskal-Wallis replaces every value by its rank and asks whether the rank totals differ more than chance — so it needs no normality and is not thrown by unequal spread.
χ² = 29.46, df = 3, p = 1.79e-06.
The conclusion survives. What changes is that you can now defend it.
Kruskal-Wallis has the same limit as F: it says something differs, not what. Dunn’s test is its post-hoc.
Comparison Z P.unadj P.adj
1 water - adult -3.120 1.78e-03 2.67e-03
2 water - larva 0.917 3.59e-01 4.31e-01
3 adult - larva 4.040 5.32e-05 1.60e-04
4 water - larva+adult -3.510 4.56e-04 9.11e-04
5 adult - larva+adult -0.381 7.03e-01 7.03e-01
6 larva - larva+adult -4.420 9.76e-06 5.85e-05
Note
Which adjustment?
"bonferroni" is the strictest and costs the most power. "bh" (Benjamini-Hochberg) controls the false discovery rate instead, and is the usual choice when you have more than a handful of comparisons.
pairwise.wilcox.test() does a similar job; Dunn’s is preferred because it uses the same ranking as the Kruskal-Wallis it follows.
| Situation | Parametric route | Non-parametric route |
|---|---|---|
| Compare 2 groups | t-test | Mann-Whitney U / Wilcoxon rank-sum |
| Compare 3+ groups | one-way ANOVA | Kruskal-Wallis |
| Post-hoc, all pairs | Tukey HSD, emmeans |
Dunn’s test, pairwise Wilcoxon |
| What it compares | means | mean ranks (≈ medians) |
| Needs | normal residuals, equal variance | neither |
| Costs you | nothing, when assumptions hold | power, and the original units |
Try a transformation before you give up the parametric route — it keeps your effect sizes interpretable and extends to the multi-factor designs coming in Week 9.
Nematode lifespan differed among parasite-exposure treatments (Kruskal-Wallis χ²3 = 29.46, p < 0.001). A Kruskal-Wallis test was used because the residuals of the equivalent ANOVA were non-normal (Shapiro-Wilk p = 0.00181) and variances were unequal (Levene p = 1.49e-12). Dunn’s tests with Benjamini-Hochberg adjustment showed that larva-exposed worms outlived worms exposed as adults or as both.
Important
✅ Say why you went non-parametric
A reader cannot tell whether you chose a rank test because your data required it or because it gave you the answer you wanted — unless you show the assumption checks that made the decision.
Report the numbers that drove the choice, not just the choice.
Today you:
lm(), the same partition, the same Fanova()fct_relevel() so the coefficients answered the experiment’s questionImportant
⚠️ A significant F is an entry ticket
It says at least one group differs. It never says which, how many, or in which direction. Everything past that comes from post-hoc comparisons you planned or corrected for.
Note
✅ Assumptions live in the residuals
Except independence, which lives in the design — and no plot will ever show it to you.
Important
✅ Set the reference level
fct_relevel(x, "control"), in the line that makes the factor. The fit never changes; the readability always does.
Note
Up next — Two-Way ANOVA
Two factors at once, and the thing that only exists when you have two: the interaction.
When code breaks — and it will, that is normal:
residuals(model) to shapiro.test() rather than the raw response? Is car loaded for leveneTest(), FSA for dunnTest()??function_name opens the help pageNote
✅ Key idea
Every working scientist googles error messages daily. Getting stuck is not failing — it is the job.