
Lecture: Nested ANOVA
Hierarchical designs
Where We Left Off
Covered last time — two-factor ANOVA:
- Example, linear model, analysis of variance
- Null hypotheses
- Interactions and main effects
- Unequal sample size
- Assumptions
✅ Key idea from last lecture
In a factorial design, every level of B appears under every level of A. Today: what happens when B is only ever nested inside one specific level of A?

Part 1 · Multifactor Designs
Two-Factor (2-Way) ANOVA, Recap
- We often consider more than 1 factor (independent categorical variables) to reduce unexplained variance and look at interactions
- 2-factor designs (2-way ANOVA) are very common in ecology
- Can have more factors (e.g., 3-way ANOVA) — interpretation gets tricky
- Most multifactor designs are: factorial, nested, or mixed

Factorial Versus Nested Designs
Consider two factors, A and B:
- Factorial/crossed: every level of B in every level of A
- Nested/hierarchical: levels of B occur only in 1 level of A
- both fixed — rare
- one fixed, one random
- both random — rare in ecology, common in genetics

Part 2 · Factorial Designs — Examples
Factorial Designs Overview
Factorial designs: both factors are typically fixed (but not always).

Factorial Design Example: Seedling Growth
Effects of light level on growth of seedlings of different size:
- 3 light levels (Factor A)
- 3 size classes (Factor B)
- 5 replicate seedlings in each cell

Factorial Design Example: Salamander Growth
Effects of food level and tadpole presence on larval salamander growth:
- 2 food levels (Factor A)
- Presence/absence of tadpoles (Factor B)
- 8 replicates in each cell

Factorial Design Example: Limpet Fecundity
Effect of season and density on limpet fecundity:
- 2 seasons (Factor A)
- 4 density treatments (Factor B)
- 3 replicates in each cell

Part 3 · Nested Designs — Examples
Nested Designs Overview
Nested designs cover: nested design examples, the linear model, analysis of variance, null hypotheses, unbalanced designs, assumptions, mixed-model ANOVA.
Nested designs: Factor A is usually fixed; Factor B is usually random.

Nested Design Example: Limpet Growth
Study on effects of enclosure size on limpet growth:
- 2 enclosure sizes (Factor A)
- 5 replicate enclosures (Factor B)
- 5 replicate limpets per enclosure

Nested Design Example: Reef Fish
Study on reef fish recruitment:
- 5 sites (Factor A)
- 6 transects at each site (Factor B)
- Replicate observations along each transect

Nested Design Example: Sea Urchin Grazing
Effects of sea urchin grazing on biomass of filamentous algae:
- 4 urchin densities: 0%, 33%, 66%, 100% of natural density
- 4 patches of rocky bottom (3–4 m²) nested in each grazing level
- 5 replicate quadrats per patch
This is our worked example for the rest of the lecture.

Part 4 · The Nested Design Model
Nested Design: Linear Model Structure
Consider a nested design with:
- p levels of factor A (i = 1…p) — e.g., 4 grazing levels
- q levels of factor B (j = 1…q), nested within each level of A — e.g., 4 different patches per grazing level
- n replicates (k = 1…n) in each combination of A and B — e.g., 5 replicate quadrats in each patch in each grazing level

Calculating Means in a Nested Design
We can calculate several means:
- Overall mean (across all levels of A and B) = \(\bar{y}\)
- A mean for each level of A (across all levels of B in that A) = \(\bar{y}_i\)
- A mean for each level of B within each A = \(\bar{y}_{j(i)}\)

Nested Design Means, Visualized

The Nested Design Linear Model
\[y_{ijk} = \mu + \alpha_i + \beta_{j(i)} + \varepsilon_{ijk}\]
- \(y_{ijk}\) is the response variable — the value of the k-th replicate in the j-th level of B in the i-th level of A (e.g., algal biomass in the 3rd quadrat, in the 2nd patch, in the low-grazing treatment)
- \(\mu\) is the overall mean (overall average algal biomass)
🖐 Notice
This looks almost identical to the two-way ANOVA model — the difference is entirely in the subscript: \(\beta_{j(i)}\), not \(\beta_j\). B only exists inside a specific A.
Fixed and Random Effects in the Nested Model
\[y_{ijk} = \mu + \alpha_i + \beta_{j(i)} + \varepsilon_{ijk}\]
- \(\alpha_i\) is the fixed effect of factor A — the difference between the average biomass in all low-grazing-level quadrats and the overall mean
- \(\beta_{j(i)}\) is the random effect of factor B nested within factor A — usually a random variable measuring variance among all possible levels of B within each level of A (variance among all possible patches that could have been used in the low-grazing treatment)
✅ Key idea
B being random is the whole point of nesting: you don’t care about these specific patches — you care about patch-to-patch variability in general.
The Error Term in the Nested Model
\[y_{ijk} = \mu + \alpha_i + \beta_{j(i)} + \varepsilon_{ijk}\]
- \(\varepsilon_{ijk}\) is the error term — the unexplained variance associated with the k-th replicate in the j-th level of B in the i-th level of A
- \(\alpha_i\) is the effect of the i-th level of A: \(\mu_i - \mu\)
- \(\varepsilon_{ijk}\): the difference between the observed algal biomass in the 3rd quadrat in the 2nd patch in the low-grazing treatment, and the predicted biomass (the average biomass in that 2nd patch in the low-grazing treatment)
🖐 Notice
The residual is measured against the cell mean (this specific patch), not the treatment mean — that’s what makes it “nested,” arithmetically.
Part 5 · Partitioning Variance
Where the Variation Goes
One quadrat’s distance from the grand mean splits into three pieces, not two:
\[\underbrace{(Y_{ijk} - \bar{Y})}_{\text{total}} = \underbrace{(\bar{Y}_i - \bar{Y})}_{A} + \underbrace{(\bar{Y}_{ij} - \bar{Y}_i)}_{B(A)} + \underbrace{(Y_{ijk} - \bar{Y}_{ij})}_{\text{residual}}\]
Square and sum over all 80 quadrats:
\[SS_{total} = SS_A + SS_{B(A)} + SS_{resid}\]
On the urchin data: 54 781 = 14 429 + 21 242 + 19 110 — the same three numbers aov() prints in Part 8.
The new row is the whole point
A one-way ANOVA has no \(SS_{B(A)}\) row — that variation lands in the residual. Here it is bigger than the treatment effect, and it is what Factor A gets tested against.
Analysis of Variance: Residual and Total
- \(SS_{resid}\) is the difference between each observation and the mean for its level of factor B, summed over all observations
- \(SS_{total} = SS_A + SS_B + SS_{resid}\)
- SS can be turned into MS by dividing by the appropriate df

Analysis of Variance: SSA
As before, we partition the variance in the response variable using SS. \(SS_A\) is the SS of differences between the mean in each level of A and the overall mean.

Analysis of Variance: SSB
\(SS_B\) is the SS of the difference between the mean in each level of B and the mean of the corresponding level of A, summed across all levels of A.

The Nested ANOVA Table

Part 6 · Hypotheses and Assumptions
Null Hypotheses: Factor A
Two hypotheses are tested on values of MS. 1. No effects of factor A — assuming A is fixed:
- \(H_0(A)\): \(\mu_1 = \mu_2 = ... = \mu_i = \mu\)
- Same as in 1-factor ANOVA, using means from B factors nested within each level of A
- (no difference in algal biomass across all grazing levels: \(\mu_{none} = \mu_{low} = \mu_{med} = \mu_{high}\))

Null Hypotheses: Factor B
2. No effects of factor B nested in A — assuming B is random:
- \(H_0(B)\): \(\sigma^2_\beta = 0\) (no variance added due to differences between all possible levels of B)
- (no variance added due to differences between patches)

Conclusions from the Analysis
Conclusion:
“Significant variation between replicate patches within each treatment, but no significant difference in the amount of filamentous algae between treatments.”

Unbalanced Nested Designs
Unequal sample sizes can arise from:
- An uneven number of B levels within each A
- An uneven number of replicates within each level of B
The second kind is the one that bites. Patch means stop being equally precise, so:
aov(Error(...))loses its exact F-tests- averaging to patch means stops being equivalent to the nested test
- you can no longer weight patches by their precision without \(\sigma^2_{within}\)
The principled fix is a mixed model — Week 11.

Nested Design Assumptions
As usual, we assume: equal variance, normality, no outliers, independence of observations.
Equal variance + normality need to be assessed at both levels:
- Since means for each level of B within each A are used for the H-test about A, need to assess whether those means meet normality and equal variance
- Examine residuals for the H-test about B
- Transformations can be used
⚠️ Watch out!
This is the part people most often get wrong: checking assumptions on the raw data isn’t enough for a nested design — you need to check them at the level each hypothesis test actually operates on.
Part 7 · The Sea Urchin Example — Data
The Nested Design, the Hard Way
This analysis examines the effects of varying sea urchin densities on the percentage cover of filamentous algae. The experiment was designed with:
- Four urchin density treatments — 0% (urchins removed), 33%, 66%, and 100% of natural density (the unmanipulated control)
- Four patches nested inside each treatment — 16 patches in all, each one receiving a single treatment
- Five replicate quadrats scored inside every patch, so 80 quadrats
📖 Reference
Quinn & Keough (2002), Experimental Design and Data Analysis for Biologists — this dataset and analysis follow their treatment of the sea urchin/algae example.
Data Overview
The data frame contains: patch (random patches 1–16 where treatments were applied), treat (urchin density treatment), quad (replicate quadrats within each treatment-patch combination), algae (percentage cover of filamentous algae — the response variable).
## Data
## # A tibble: 6 × 4
## treat patch quad algae
## <fct> <fct> <dbl> <dbl>
## 1 100% 1 1 0
## 2 100% 1 2 0
## 3 100% 1 3 0
## 4 100% 1 4 6
## 5 100% 1 5 2
## 6 100% 2 1 0
##
##
## Treatment levels:
## [1] "0%" "33%" "66%" "100%"
## Summary statistics by treatment:
## # A tibble: 4 × 7
## treat n mean sd se min max
## <fct> <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0% 20 39.2 28.7 6.41 0 83
## 2 33% 20 19 25.7 5.74 0 71
## 3 66% 20 21.6 25.1 5.62 0 79
## 4 100% 20 1.3 3.18 0.711 0 13
A First Look at the Data

Part 8 · Fitting It — the Wrong Way, Then the Right One
Manual Nested ANOVA — Setting Up
In this experimental design, patch is nested within treat because each patch received only one treatment level — a hierarchical design where the effect of patches must be considered within each treatment (Quinn & Keough 2002).
We’ll use a traditional nested ANOVA. You could make this easy by averaging algae per patch and running a one-way ANOVA on the 16 patch means. On a perfectly balanced design that gives the identical F, df and p — but it throws away the within-patch variance, and the moment the design stops being balanced that shortcut starts costing you. We come back to both halves of this at the end of the lecture.
This first attempt is really a regular factorial ANOVA — not appropriate, since it’s pseudoreplicated. Also, the F-value for Treat is using \(MS_{resid}\) when it should use \(MS_{treat:patch}\).
u_fact_model <- aov(algae ~ treat + treat:patch, data = u_df)
cat("Factorial model summary:\n\n")
## Factorial model summary:
summary(u_fact_model)
## Df Sum Sq Mean Sq F value Pr(>F)
## treat 3 14429 4810 16.108 0.0000000658 ***
## treat:patch 12 21242 1770 5.928 0.0000008323 ***
## Residuals 64 19110 299
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Why the Factorial ANOVA Is Not Correct
- We can correct this by specifying the error term
- The MS needed for the F-value is really MS(patch nested in treatment), because those are the true replicates
- We can specify that — the problem is this approach can’t really handle unbalanced designs
# Error(patch) -- not Error(treat:patch). Patch IDs are already unique
# (1-16) and each patch sits in exactly one treatment, so naming the patch
# alone is enough. Writing treat:patch asks for 4 x 16 = 64 combinations
# when only 16 exist, which makes R warn "Error() model is singular".
# Both give identical numbers; this one does it without the warning.
u_nested_model <- aov(algae ~ treat + Error(patch), data = u_df)
cat("Correctly specified nested ANOVA:\n\n")
## Correctly specified nested ANOVA:
summary(u_nested_model)
##
## Error: patch
## Df Sum Sq Mean Sq F value Pr(>F)
## treat 3 14429 4810 2.717 0.0913 .
## Residuals 12 21242 1770
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 64 19110 298.6Same Numerator, Different Denominator
Both models used the same 80 quadrats and produced the same sums of squares. Only the F-test’s denominator changed:
| MS treat | ÷ MS error | = F | p | |
|---|---|---|---|---|
| Naive factorial | 4810 | ÷ 299 (quadrats) | 16.11 | 0.00000007 |
| Nested | 4810 | ÷ 1770 (patches) | 2.72 | 0.091 |
The p-value moved by six orders of magnitude.
The naive model tested the treatment against variation between quadrats within a patch — but quadrats are not what we assigned treatments to. Patches are. 80 quadrats, but only 16 patches: 12 denominator df, not 64.
✅ How to find the right denominator
Ask: what did I randomly assign the treatment to?
That unit — patch here — supplies the error term. Everything measured inside it is a subsample, and subsamples tell you about measurement noise, not about treatment effects.
🖐 Notice MS(patch) is large
1770 vs 299 — patches genuinely differ from one another. That variability is real biology, and the nested model makes you pay for it. The naive model pretended it was not there.
Where You Have Seen This Before
This is the pine needle problem again, one level up:
| Our pine study | This urchin study |
|---|---|
| 8 needles per side | 5 quadrats per patch |
| within one tree | within one patch |
| 4 trees are the replicates | 16 patches are the replicates |
| 64 rows, n = 4 | 80 rows, n = 16 |
In Week 1 we solved it by averaging each tree to one number. That works, and it is honest.
✅ Key idea
Averaging throws the subsamples away. Nesting uses them — to estimate how noisy your measurements are — while still refusing to count them as replicates. On a balanced design the treatment test comes out the same either way; the variance components never do, and once the design is unbalanced neither does the test.
⚠️ The limit of aov(Error(...))
This syntax needs a (near-)balanced design. Ours is perfectly balanced — 4 patches per treatment, 5 quadrats per patch — so it works.
Real data usually is not. That is the problem mixed models solve, and it is where we go next.
How Do I Know If My Design Is Nested?
Three questions, in order. Ask them of your own data before you fit anything.
1. What did I apply the treatment to? → Whole patches. Not quadrats.
2. Does every level of B appear under every level of A? → Patch 7 only ever received the 66% density. It cannot appear under 100%. → No ⇒ nested. (Yes ⇒ crossed/factorial.)
3. How many of those units do I have per treatment? → 4 patches, not 20 quadrats. That is your n.
The answer to Q1 is your error term.
aov(algae ~ treat + Error(patch), data = u_df)
# ^^^^^^^^^^^^
# the thing treatments were assigned to🖐 A quick tell
If two rows of your data share a value for the thing you randomized, they are not independent replicates of it.
Two quadrats in patch 3 share patch 3. Two needles on one tree share that tree.
⚠️ Nested factors are random
We are not asking “is patch 7 different from patch 12?” — those specific patches are of no interest. They are a sample of possible patches, so patch enters as a random effect.
The factor above it (treat) is fixed: those four densities are exactly the ones we care about.
Part 9 · Post-Hoc, and What the Test Really Is
Post-Hoc on a Nested Model
u_emm <- emmeans(u_nested_model, ~ treat)
u_emm
## treat emmean SE df lower.CL upper.CL
## 0% 39.2 9.41 12 18.70 59.7
## 33% 19.0 9.41 12 -1.50 39.5
## 66% 21.6 9.41 12 1.05 42.0
## 100% 1.3 9.41 12 -19.20 21.8
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
multcomp::cld(u_emm, Letters = letters)
## treat emmean SE df lower.CL upper.CL .group
## 100% 1.3 9.41 12 -19.20 21.8 a
## 33% 19.0 9.41 12 -1.50 39.5 a
## 66% 21.6 9.41 12 1.05 42.0 a
## 0% 39.2 9.41 12 18.70 59.7 a
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 4 estimates
## significance level used: alpha = 0.05
## NOTE: If two or more means share the same grouping symbol,
## then we cannot show them to be different.
## But we also did not show them to be the same.Means from 1.3 % to 39.2 % — a thirty-fold spread — and every one of them shares the letter a.
The SE is 9.41 on every mean, and it comes from \(MS_{B(A)}\), the patch stratum. With only 4 patches per treatment there is simply not enough replication to separate them.
⚠️ emmeans will warn you
“EMMs are biased unless design is perfectly balanced.” Ours is balanced — 4 patches per treatment, 5 quadrats per patch — so read past it. The warning fires on any aovlist object.
The Nested F Is a One-Way ANOVA on the Patch Means — If Balanced
# Collapse all 80 quadrats to the 16 patch means ...
patch_df <- u_df %>%
group_by(treat, patch) %>%
summarise(mean_algae = mean(algae), .groups = "drop")
# ... and run an ordinary one-way ANOVA on those 16 rows
summary(aov(mean_algae ~ treat, data = patch_df))
## Df Sum Sq Mean Sq F value Pr(>F)
## treat 3 2886 961.9 2.717 0.0913 .
## Residuals 12 4248 354.0
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1F = 2.717 on 3 and 12 df, p = 0.091 — identical to the nested model, to every digit. The SS are exactly one fifth, and dividing numerator and denominator by the same 5 leaves F untouched.
✅ What the subsamples bought you
Not df for the treatment test — that was 12 either way. They made each patch’s value more precise, and let you estimate the within-patch variance. For power against the treatment you need more patches: twenty quadrats in four patches is still n = 4.
⚠️ All of this depends on BALANCE — 5 quadrats in every patch. Change that and the equality breaks. Next slide.
…And Why Averaging Breaks When You Are Not
# a storm scours all but one quadrat in four patches
thin_df <- u_df %>%
filter(!(patch %in% c(2, 6, 10, 14) & quad > 1))
# nested ANOVA on the thinned data
summary(aov(algae ~ treat + Error(patch), thin_df))[[1]]
## Df Sum Sq Mean Sq F value Pr(>F)
## treat 3 8837 2946 2.25 0.135
## Residuals 12 15713 1309
# the patch-mean shortcut, same data
pm <- thin_df %>% group_by(treat, patch) %>%
summarise(m = mean(algae), .groups = "drop")
anova(aov(m ~ treat, pm))
## Analysis of Variance Table
##
## Response: m
## Df Sum Sq Mean Sq F value Pr(>F)
## treat 3 3174.1 1058.0 2.9628 0.07493 .
## Residuals 12 4285.2 357.1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1nested aov |
patch means | |
|---|---|---|
| balanced | F=2.72, p=0.091 | F=2.72, p=0.091 |
| unbalanced | F=2.25, p=0.135 | F=2.96, p=0.075 |
Same data — and the shortcut is now the most optimistic.
⚠️ The patch means are no longer equally precise
\[\operatorname{Var}(\bar y_j) = \sigma^2_{patch} + \frac{\sigma^2_{within}}{n_j}\]
A patch with 1 quadrat is far noisier than one with 5 — but the ANOVA on the means treats them as interchangeable, violating its own equal-variance assumption by construction.
To weight them you need \(\sigma^2_{within}\) — which averaging already threw away. That is the cost of losing sight of the variance inside the nested groups.
Effect Size — How Big, Not Just Whether
p = 0.091 says “we could not resolve it.” It does not say “there is nothing there.” The effect size says how big the thing we could not resolve is.
effectsize takes the aovlist straight, and the Group column tells you which error stratum it used — patch, the correct denominator.
eta_squared(u_nested_model)
## # Effect Size for ANOVA (Type I)
##
## Group | Parameter | Eta2 (partial) | 95% CI
## -------------------------------------------------
## patch | treat | 0.40 | [0.00, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].
omega_squared(u_nested_model)
## # Effect Size for ANOVA (Type I)
##
## Group | Parameter | Omega2 (partial) | 95% CI
## ---------------------------------------------------
## patch | treat | 0.24 | [0.00, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].| Model | partial eta2 | p |
|---|---|---|
| Nested (correct) | 0.40 | 0.091 |
| Naive factorial (wrong) | 0.43 | 0.00000007 |
🖐 Pseudoreplication corrupts the p-value, not the effect size
The wrong model and the right model agree the effect is large — 0.40 against 0.43. What they disagree about, by seven orders of magnitude, is how confident you are allowed to be.
So the honest reading of this experiment is not “urchins do not matter.” It is “a large effect that 4 patches per treatment could not resolve.” That is a power problem, and the effect size is what lets you say so.
Reporting a Nested ANOVA
How to report this result
“Algal cover did not differ significantly among urchin density treatments (nested ANOVA, F₃,₁₂ = 2.72, p = 0.091), with patches treated as a random factor nested within treatment. Variation among patches within treatments was substantial, accounting for approximately 50% of the total variance (σ²_patch = 294, σ²_within = 299; n = 4 patches per treatment, 5 quadrats per patch).”
The three things people get wrong:
- The df. It is F(3, 12) — from the patch stratum, not the 64 on the “Within” line. F(3, 64) is the wrong test.
- Silence about the nesting. Say “nested” and name the random factor; “F = 2.72” alone tells a reader nothing.
- Throwing away the patch variance. It is a result, not a nuisance — here it is half the story.
# Variance components, straight from
# the two mean squares
MS_patch <- 1770.2 # Error: patch
MS_within <- 298.6 # Error: Within
n_quad <- 5 # quadrats per patch
var_patch <- (MS_patch - MS_within) / n_quad
var_within <- MS_within
cat("among patches:", round(var_patch, 1), "
")
## among patches: 294.3
cat("within patch: ", round(var_within, 1), "
")
## within patch: 298.6
cat("% among patches:",
round(100 * var_patch /
(var_patch + var_within), 1), "%
")
## % among patches: 49.6 %🖐 Half the variance is between patches
Ignore that and you are pretending 80 quadrats are 80 independent observations. They are closer to 16.
Wrap-up — What You Did Today
- Told a nested design from a crossed one, and learned the three questions that decide it
- Found the right error term by asking what the treatment was applied to
- Watched the same 80 quadrats give F = 16.11 (p < 0.0000001) the wrong way and F = 2.72 (p = 0.091) the right way
- Proved the nested F is just a one-way ANOVA on the 16 patch means — and watched that equality fall apart the moment the design went unbalanced
- Estimated variance components — half the variation lives between patches
- Put an effect size next to the p-value, and found a large effect we were not powered to resolve
- Reported it with the correct degrees of freedom
⚠️ The mistake this lecture exists to prevent
Using the residual mean square as your denominator when the treatment was applied to something coarser than a single measurement.
It inflates F, shrinks p, and it is the most common error in ecological field experiments (Hurlbert 1984).
Wrap-up — This Has Been the Plan All Along
| Week | Same idea, new tool |
|---|---|
| 1 | 64 needles, but 4 trees |
| 3 | two-stage summarise() |
| 4 | the paired t-test |
| 5 | power: 44 trees vs 7 |
| 9 | nested ANOVA |
Today it finally got its own name and its own model.
Sums of squares were never the problem — the denominator was.
🖐 Before next class
Take a dataset you actually care about — yours, or one from a paper — and answer the three questions from “How Do I Know If My Design Is Nested?”.
Write down what your true n is. If it is smaller than the number of rows, you have found a nested design.
Up next — Linear Mixed Models
These same data come back with patch as a random effect instead of an error term. Error() is the 1960s way to get the right denominator; a random effect is the modern way — and it handles unbalanced designs, which Error() does not.
Readings
- Gotelli & Ellison, Ch. 10 — The Analysis of Variance (the nested / hierarchical designs section).
- Quinn & Keough (2002), Experimental Design and Data Analysis for Biologists, Ch. 9 — Multifactor Analysis of Variance (§9.1, nested designs). This is the source of today’s worked example and of the ANOVA table on the hypothesis slides.
- Hurlbert, S. H. (1984). Pseudoreplication and the design of ecological field experiments. Ecological Monographs, 54(2), 187–211.
📖 On the dataset
The sea urchin grazing data are Andrew & Underwood’s, as analysed in Quinn & Keough Ch. 9. The textbook PDFs are on Canvas.
🖐 Why no Whitlock & Schluter here
W&S Ch. 18 covers multiple explanatory variables but does not treat nested designs. For this topic, Gotelli and Quinn & Keough are the references.