Main effects and interactions
A main effect is an effect of a factor which does not depend on the level of the other factor.
An interaction is any dependence on the effect of one factor on the level of the other factor.
STAT 516 Lec 06
Two-way factorial design (balanced)
Karl Gregory
2026-10-26
Tensile strength data from Kuehl (2000)
Tensile strength
Side-by-side boxplots of responses at all factor level combinations:
Tensile strength
A dotplot is better since there are only a small number of EUs per treatment group.
The two-way factorial experimental design
- Two factors of interest.
- Each factor comprehends a set of treatments, called its levels.
- EUs randomly assigned to treatments.
- Each treatment is a unique combinations of factor levels.
- If Factor A has \(a\) levels and Factor B has \(b\) levels.
- There are \(ab\) treatment groups.
We want to make inferences about:
- The effects of each factor.
- The interactions between the two factors.
- Various differences in treatment group means.
Discuss: Give factors and their levels in the tensile strength experiment.
Two-way treatment effects model
Suppose the responses arise as \[ Y_{ijk} = \mu + \tau_i + \gamma_j + (\tau\gamma)_{ij} + \varepsilon_{ijk} \] for \(i = 1,\dots,a\), \(j = 1,\dots,b\), and \(k = 1,\dots,n_{ij}\), where
- \(Y_{ijk}\) is the response for EU \(k\) under level \(i\) of A and level \(j\) of B.
- \(\mu\) represents a baseline or overall mean
- The \(\tau_i\) are the main effects for Factor A.
- The \(\gamma_j\) are the main effects for Factor B.
- The \((\tau\gamma)_{ij}\) are the interaction effects between A and B.
- The \(\varepsilon_{ijk}\) are \(\mathcal{N}(0,\sigma^2)\) error terms.
Assume for now
- a balanced design, i.e. \(n_{ij} = n\) for all \(i,j\)
- with replication, i.e. \(n \geq 2\).
Cell means model representation
Assume \[ Y_{ijk} = \mu_{ij} + \varepsilon_{ijk} \] for \(i = 1,\dots,a\), \(j = 1,\dots,b\), and \(k = 1,\dots,n_{ij}\), where
- \(Y_{ijk}\) is the response for EU \(k\) under level \(i\) of A and level \(j\) of B.
- \(\mu_{ij}\) is the response mean under level \(i\) of A and level \(j\) of B.
- The \(\varepsilon_{ijk}\) are \(\mathcal{N}(0,\sigma^2)\) error terms.
We have merely written \(\mu_{ij} = \mu + \tau_i + \gamma_j + (\tau\gamma)_{ij}\) for all \(i,j\).
Goals in two-way factorial experiments
In the two-way treatment effects model \[ Y_{ijk} = \mu + \tau_i + \gamma_j + (\tau\gamma)_{ij} + \varepsilon_{ijk}, \] where \(\varepsilon_{ijk}\overset{\operatorname{ind}}{\sim}\mathcal{N}(0,\sigma^2)\), we wish to
- Visualize the data.
- Estimate the parameters \(\mu\), \(\tau_i\), \(\gamma_j\), \((\tau\gamma)_{ij}\).
- Estimate the error term variance \(\sigma^2\).
- Check whether the model assumptions are satisfied.
- Decompose the variation in the \(Y_{ijk}\) into its sources.
- Test whether there is any effect of the factors on the response.
- Test for main effects and interaction effects.
- Do multiple comparisons.
Important quantities
Treatment group means: \(\displaystyle \bar Y_{ij.} = \frac{1}{n} \sum_{k=1}^n Y_{ijk}\) for \(i=1,\dots,a\), \(j=1,\dots,b\).
Marginal means: \[\begin{align*} \bar Y_{i..} &= \frac{1}{na}\sum_{j=1}^b \sum_{k=1}^n Y_{ijk}, \quad i = 1,\dots,a\\ \bar Y_{.j.} &= \frac{1}{nb}\sum_{i=1}^a \sum_{k=1}^n Y_{ijk}, \quad j = 1,\dots,b. \end{align*}\]
Treatment group standard deviations: \(\displaystyle S_{ij}^2 = \frac{1}{n-1}\sum_{k=1}^n(Y_{ijk} - \bar Y_{ij.})^2\)
Overall mean \(\displaystyle \bar Y_{...} = \frac{1}{nab}\sum_{i=1}^a\sum_{j=1}^b \sum_{k=1}^{n}Y_{ijk}\).
Tensile strength
Consider how to compute the quantities on the previous slide:
Estimating the mean and treatment effects parameters
We want estimators \(\hat \mu\), \(\hat \tau_i\), \(\hat \gamma_j\), and \(\hat{(\tau\gamma)_{ij}}\) which satisfy \[\begin{align*} \hat \mu + \hat \tau_1 + \hat \gamma_1 + \hat{(\tau\gamma)}_{11} &= \bar Y_{11.}\\ &~\vdots\\ \hat \mu + \hat \tau_1 + \hat \gamma_b + \hat{(\tau\gamma)}_{1b} &= \bar Y_{1b.}\\ &~\vdots\\ \hat \mu + \hat \tau_a + \hat \gamma_1 + \hat{(\tau\gamma)}_{a1} &= \bar Y_{a1.}\\ &~\vdots\\ \hat \mu + \hat \tau_a + \hat \gamma_b + \hat{(\tau\gamma)}_{ab} &= \bar Y_{ab.} \end{align*}\]
Above has \(ab\) equations and \(a + b + ab + 1\) unknowns! Can fix by setting \[ \hat \tau_1 = 0, \quad \hat \gamma_1 = 0, \quad \text{ and }\quad \hat{(\tau\gamma)}_{1j} = \hat{(\tau\gamma)}_{i1} = 0 \text{ for all } i,j. \] Then \(\hat \mu\) estimates the mean of the baseline group \(i=1\), \(j=1\).
Tensile strength
We have \(a = 2\) and \(b=4\). So the complete system of equations is
\[\begin{align*} \hat \mu + \hat \tau_1 + \hat \gamma_1 + \hat{(\tau\gamma)}_{11} &= \bar Y_{11.}\\ \hat \mu + \hat \tau_1 + \hat \gamma_2 + \hat{(\tau\gamma)}_{12} &= \bar Y_{12.}\\ \hat \mu + \hat \tau_1 + \hat \gamma_3 + \hat{(\tau\gamma)}_{13} &= \bar Y_{13.}\\ \hat \mu + \hat \tau_1 + \hat \gamma_4 + \hat{(\tau\gamma)}_{14} &= \bar Y_{14.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_1 + \hat{(\tau\gamma)}_{21} &= \bar Y_{21.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_2 + \hat{(\tau\gamma)}_{22} &= \bar Y_{22.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_3 + \hat{(\tau\gamma)}_{23} &= \bar Y_{23.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_4 + \hat{(\tau\gamma)}_{24} &= \bar Y_{24.} \end{align*}\]
Tensile strength
After setting to zero every term with an index equal to 1, we have
\[\begin{align*} \hat \mu &= \bar Y_{11.}\\ \hat \mu + \hat \gamma_2 &= \bar Y_{12.}\\ \hat \mu + \hat \gamma_3 &= \bar Y_{13.}\\ \hat \mu + \hat \gamma_4 &= \bar Y_{14.}\\ \hat \mu + \hat \tau_2 &= \bar Y_{21.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_2 + \hat{(\tau\gamma)}_{22} &= \bar Y_{22.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_3 + \hat{(\tau\gamma)}_{23} &= \bar Y_{23.}\\ \hat \mu + \hat \tau_2 + \hat \gamma_4 + \hat{(\tau\gamma)}_{24} &= \bar Y_{24.} \end{align*}\]
Tensile strength
After re-arranging to solve for the parameter estimates we obtain
\[\begin{align*} \hat \mu &= \bar Y_{11.}\\ \hat \gamma_2 &= \bar Y_{12.} - \bar Y_{11.}\\ \hat \gamma_3 &= \bar Y_{13.} - \bar Y_{11.}\\ \hat \gamma_4 &= \bar Y_{14.} - \bar Y_{11.}\\ \hat \tau_2 &= \bar Y_{21.} - \bar Y_{11.}\\ \hat{(\tau\gamma)}_{22} &= \bar Y_{22.} - \bar Y_{11.} - (\bar Y_{12.} - \bar Y_{11.}) - (\bar Y_{21.} - \bar Y_{11.})\\ \hat{(\tau\gamma)}_{23} &= \bar Y_{23.} - \bar Y_{11.} - (\bar Y_{13.} - \bar Y_{11.}) - (\bar Y_{21.} - \bar Y_{11.})\\ \hat{(\tau\gamma)}_{24} &= \bar Y_{24.} - \bar Y_{11.} - (\bar Y_{14.} - \bar Y_{11.}) - (\bar Y_{21.} - \bar Y_{11.}) \end{align*}\]
Tensile strength
Call:
lm(formula = y ~ agg + comp + agg:comp)
Coefficients:
(Intercept) aggS compr compst compvl aggS:compr
97.3 -36.7 31.7 -32.0 -40.0 18.7
aggS:compst aggS:compvl
39.0 21.0 See how the parameter estimates build the treatment means:
Hypotheses of interest in the two-way factorial experiment
\(H_0\): Factor A has no main effect. \[ \text{$H_0$: $\tau_1 = \dots = \tau_a$}. \]
\(H_0\): Factor B has no main effect. \[ \text{$H_0$: $\gamma_1 = \dots = \gamma_b$}. \]
\(H_0\): There is no interaction between Factor A and Factor B. \[ \text{$H_0$: $(\tau\gamma)_{11} = \dots = (\tau\gamma)_{ab}$.} \] (Factor A effects are same at all levels of Factor B and vice versa.)
Estimating the error term variance
In the two-way factorial model the
- fitted values are the treatment group means, i.e. \(\hat Y_{ijk} = \bar Y_{ij.}\).
- residuals are the deviations from the group means \[ \hat \varepsilon_{ijk} = Y_{ijk} - \hat Y_{ijk} = Y_{ijk} - \bar Y_{ij.} \].
An unbiased estimator of the error term variance \(\sigma^2\) is given by \[ \hat \sigma^2 = \frac{1}{ab}\sum_{i=1}^a\sum_{j=1}^b (n - 1)S_{ij}^2=\frac{1}{ab(n-1)} \sum_{i=1}^a\sum_{j=1}^b\sum_{k=1}^n\varepsilon_{ijk}^2. \] This estimate pools the variances in all the treatment groups.
Checking model assumptions
Validity of the methods in these slides depends on these assumptions:
The responses are normally distributed around the treatment means (Check QQ plot of residuals).
The response has the same variance in all treatment groups (Check residuals vs fitted values plot).
The response values are independent of each other (No way to check; must trust experimental design).
Tensile strength
Normal quantile-quantile plot of the residuals:
Tensile strength
Residuals versus fitted values plot:
Sums of squares in the two-way factorial experiment
| SS | Symbol | Formula |
|---|---|---|
| Total | \(\operatorname{SS}_{\operatorname{Tot}}\) | \(\sum_{i=1}^a\sum_{j=1}^b\sum_{k=1}^n(Y_{ijk} - \bar Y_{...})^2\) |
| Treatments | \(\operatorname{SS}_{\operatorname{Trt}}\) | \(\sum_{i=1}^a \sum_{j=1}^b n(\bar Y_{ij.} - \bar Y_{...})^2\) |
| Error | \(\operatorname{SS}_{\operatorname{Error}}\) | \(\sum_{i=1}^a\sum_{j=1}^b\sum_{k=1}^n(Y_{ijk} - \bar Y_{ij.})^2\) |
- We have \(\operatorname{SS}_{\operatorname{Tot}}= \operatorname{SS}_{\operatorname{Trt}}+ \operatorname{SS}_{\operatorname{Error}}\).
- We again define \(\displaystyle R^2 = \frac{\operatorname{SS}_{\operatorname{Trt}}}{\operatorname{SS}_{\operatorname{Tot}}}\).
- Can further decompose \(\operatorname{SS}_{\operatorname{Trt}}\) into SS from main effects and interactions.
ANOVA table for overall F test
| Source | Df | SS | MS | F value | p value |
|---|---|---|---|---|---|
| Treatments | \(ab-1\) | \(\operatorname{SS}_{\operatorname{Trt}}\) | \(\operatorname{MS}_{\operatorname{Trt}}\) | \(\operatorname{F}_{\operatorname{test}}\) | \(P(F > \operatorname{F}_{\operatorname{test}})\) |
| Error | \(ab(n-1)\) | \(\operatorname{SS}_{\operatorname{Error}}\) | \(\operatorname{MS}_{\operatorname{Error}}\) | ||
| Total | \(abn-1\) | \(\operatorname{SS}_{\operatorname{Tot}}\) |
The F statistic \(\displaystyle \operatorname{F}_{\operatorname{test}}= \frac{\operatorname{MS}_{\operatorname{Trt}}}{\operatorname{MS}_{\operatorname{Error}}}\) tests \[\begin{align*} &\text{$H_0$: All the $\mu_{ij}$ are the same.}\\ &\text{$H_1$: The $\mu_{ij}$ are not all the same.} \end{align*}\] where \(\mu_{ij} = \mu + \tau_i + \gamma_i + (\tau\gamma)_{ij}\).
Reject \(H_0\) if \(\operatorname{F}_{\operatorname{test}}> F_{ab-1,ab(n-1),\alpha}\).
Tensile strength
Exercise: Make ANOVA table for the tensile strength data. Interpret!
Tensile strength
| Source | Df | SS | MS | F value | p value |
|---|---|---|---|---|---|
| Treatments | 7 | 19122.50 | 2731.79 | 287.5564 | 0.0000 |
| Error | 16 | 152.00 | 9.50 | ||
| Total | 23 | 19274.50 |
Further decomposition of treatments sum of squares
Define main effect and interaction sums of squares \(\operatorname{SS}_{\operatorname{A}}\), \(\operatorname{SS}_{\operatorname{B}}\), and \(\operatorname{SS}_{\operatorname{AB}}\):
- \(\operatorname{SS}_{\operatorname{A}}= bn\sum_{i=1}^a(\bar Y_{i..} - \bar Y_{...})^2\)
- \(\operatorname{SS}_{\operatorname{B}}= an\sum_{j=1}^b(\bar Y_{.j.} - \bar Y_{...})^2\)
- \(\operatorname{SS}_{\operatorname{AB}}= n\sum_{i=1}^a\sum_{j=1}^b(\bar Y_{ij.} - (\bar Y_{i..} + \bar Y_{.j.} - \bar Y_{...}))^2\)
In the balanced design we have \[ \operatorname{SS}_{\operatorname{Trt}}= \operatorname{SS}_{\operatorname{A}}+ \operatorname{SS}_{\operatorname{B}}+ \operatorname{SS}_{\operatorname{AB}}. \]
Full sums of squares in two-way factorial experiment
| SS | Symbol | Formula |
|---|---|---|
| Total | \(\operatorname{SS}_{\operatorname{Tot}}\) | \(\sum_{i=1}^a\sum_{j=1}^b\sum_{k=1}^n(Y_{ijk} - \bar Y_{...})^2\) |
| A | \(\operatorname{SS}_{\operatorname{A}}\) | \(bn\sum_{i=1}^a (\bar Y_{i..} - \bar Y_{...})^2\) |
| B | \(\operatorname{SS}_{\operatorname{B}}\) | \(an\sum_{j=1}^b (\bar Y_{.j.} - \bar Y_{...})^2\) |
| AB | \(\operatorname{SS}_{\operatorname{AB}}\) | \(n\sum_{i=1}^a \sum_{j=1}^b (Y_{ij.} - (\bar Y_{i..} + \bar Y_{.j.} - \bar Y_{...}))^2\) |
| Error | \(\operatorname{SS}_{\operatorname{Error}}\) | \(\sum_{i=1}^a\sum_{j=1}^b\sum_{k=1}^n(Y_{ijk} - \bar Y_{ij.})^2\) |
- We have \(\operatorname{SS}_{\operatorname{Tot}}= \operatorname{SS}_{\operatorname{A}}+ \operatorname{SS}_{\operatorname{B}}+ \operatorname{SS}_{\operatorname{AB}}+ \operatorname{SS}_{\operatorname{Error}}\).
Full ANOVA table for balanced two-way factorial design
| Source | Df | SS | MS | F value |
|---|---|---|---|---|
| A | \(a-1\) | \(\operatorname{SS}_{\operatorname{A}}\) | \(\operatorname{MS}_{\operatorname{A}}\) | \(F_{\operatorname{A}} = \operatorname{MS}_{\operatorname{A}}/ \operatorname{MS}_{\operatorname{Error}}\) |
| B | \(b-1\) | \(\operatorname{SS}_{\operatorname{B}}\) | \(\operatorname{MS}_{\operatorname{B}}\) | \(F_{\operatorname{B}} = \operatorname{MS}_{\operatorname{B}}/\operatorname{MS}_{\operatorname{Error}}\) |
| AB | \((a-1)(b-1)\) | \(\operatorname{SS}_{\operatorname{AB}}\) | \(\operatorname{MS}_{\operatorname{AB}}\) | \(F_{\operatorname{AB}} = \operatorname{MS}_{\operatorname{AB}}/ \operatorname{MS}_{\operatorname{Error}}\) |
| Error | \(ab(n-1)\) | \(\operatorname{SS}_{\operatorname{Error}}\) | \(\operatorname{MS}_{\operatorname{Error}}\) | |
| Total | \(abn-1\) | \(\operatorname{SS}_{\operatorname{Tot}}\) |
- Reject \(H_0\): no Factor A main effect if \(F_{\operatorname{A}} > F_{a-1,ab(n-1),\alpha}\).
- Reject \(H_0\): no Factor B main effect if \(F_{\operatorname{B}} > F_{b-1,ab(n-1),\alpha}\).
- Reject \(H_0\): no A and B interaction if \(F_{\operatorname{AB}} > F_{(a-1)(b-1),ab(n-1),\alpha}\).
Tensile strength
Obtain ANOVA table with anova() function on lm() output.
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
agg 1 1734 1734 182.5 3.6e-10 ***
comp 3 16243 5414 570.0 < 2e-16 ***
agg:comp 3 1145 382 40.2 1.1e-07 ***
Residuals 16 152 10
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Important: anova() function only appropriate for a balanced design.
Interaction is significant. Now what?
When you find a significant interaction:
- Make interaction plots (next slides).
- Be very cautious about interpreting main effects, even when these are statistically significant.
Tensile strength
Use the interaction.plot() function to visualize an interaction:
Tensile strength
Interactions appear as crossing lines or differing slopes.
Population marginal means and overall mean
Define the population marginal means and the overall mean as
- \(\displaystyle \bar \mu_{i.} = \frac{1}{b}\sum_{j=1}^b \mu_{ij}, \quad i = 1,\dots,a\).
- \(\displaystyle \bar \mu_{.j} = \frac{1}{a}\sum_{i=1}^a \mu_{ij}, \quad j = 1,\dots,b\).
- \(\displaystyle \bar \mu_{..} = \frac{1}{ab}\sum_{i=1}^a\sum_{j=1}^b\mu_{ij}\).
Tensile strength
Interpret \(\bar \mu_{i.}\), \(\bar \mu_{.j}\), and \(\bar \mu_{..}\) for the tensile data:
Estimates of cell and marginal means in the balanced case
The estimators of the cell and marginal means are given by
- \(\hat \mu_{ij} = \bar Y_{ij.}\), \(i = 1,\dots,a\), \(j = 1,\dots,b\)
- \(\hat{\bar \mu}_{i.} = \bar Y_{i..}\), \(i = 1,\dots,a\).
- \(\hat{\bar \mu}_{.j} = \bar Y_{.j.}\), \(j = 1,\dots,b\).
We estimate \(\hat{\bar \mu}_{i.}\) with \(\bar Y_{i..}\) and \(\hat{\bar \mu}_{.j}\) with \(\bar Y_{.j.}\) only when \(n_{ij} =n~\forall ij\).
Some CI formulas (without familywise adjustment)
These CI formulas are for the balanced design \(n_{ij} = n~\forall ij\).
| Target | \((1-\alpha)100\%\) confidence interval |
|---|---|
| \(\mu_{ij}\) | \(\bar Y_{ij.} \pm t_{ab(n-1),\alpha/2}\hat \sigma \sqrt{\frac{1}{n}}\) |
| \(\mu_{ij} - \mu_{i'j'}\) | \(\bar Y_{ij.} - \bar Y_{i'j'.} \pm t_{ab(n-1),\alpha/2}\hat \sigma \sqrt{\frac{2}{n}}\) |
| \(\bar \mu_{i.}\) | \(\bar Y_{i..} \pm t_{ab(n-1),\alpha/2}\hat \sigma\sqrt{\frac{1}{bn}}\) |
| \(\bar \mu_{.j}\) | \(\bar Y_{.j.} \pm t_{ab(n-1),\alpha/2}\hat \sigma\sqrt{\frac{1}{an}}\) |
| \(\bar \mu_{i.} - \bar \mu_{i'.}\) | \(\bar Y_{i..} - \bar Y_{i'..} \pm t_{ab(n-1),\alpha/2}\hat \sigma \sqrt{\frac{2}{bn}}\) |
| \(\bar \mu_{.j} - \bar \mu_{.j'}\) | \(\bar Y_{.j.} - \bar Y_{.j'.}\pm t_{ab(n-1),\alpha/2}\hat \sigma \sqrt{\frac{2}{an}}\) |
In the above \(\hat \sigma = \sqrt{\operatorname{MS}_{\operatorname{Error}}}\).
Comparing means at all factor level combinations
- Tukey’s for comparing all pairs among \(\mu_{ij}\), \(i=1,\dots,a\), \(j = 1,\dots,b\): \[ \bar Y_{ij.} - \bar Y_{i'j'.} \pm q_{ab,ab(n-1),\alpha} ~ \hat \sigma \frac{1}{\sqrt{n}}, \quad (i,j)\neq (i',j'). \]
- Dunnett’s for comparing all means \(\mu_{ij}\) to a baseline \(\mu_{ab}\): \[ \bar Y_{ij.} - \bar Y_{ab.} \pm d_{ab,ab(n-1),\alpha} ~ \hat \sigma \sqrt{\frac{2}{n}}, \quad (i,j)\neq (a,b). \] Use \(\hat \sigma = \sqrt{\operatorname{MS}_{\operatorname{Error}}}\).
These may make more comparisons than are of interest
Tensile strength
Tukey multiple comparisons of means
95% family-wise confidence level
Fit: aov(formula = lm(y ~ agg:comp))
$`agg:comp`
diff lwr upr p adj
S:l-B:l -36.67 -45.38 -27.954 0.000
B:r-B:l 31.67 22.95 40.380 0.000
S:r-B:l 13.67 4.95 22.380 0.001
B:st-B:l -32.00 -40.71 -23.287 0.000
S:st-B:l -29.67 -38.38 -20.954 0.000
B:vl-B:l -40.00 -48.71 -31.287 0.000
S:vl-B:l -55.67 -64.38 -46.954 0.000
B:r-S:l 68.33 59.62 77.046 0.000
S:r-S:l 50.33 41.62 59.046 0.000
B:st-S:l 4.67 -4.05 13.380 0.596
S:st-S:l 7.00 -1.71 15.713 0.168
B:vl-S:l -3.33 -12.05 5.380 0.877
S:vl-S:l -19.00 -27.71 -10.287 0.000
S:r-B:r -18.00 -26.71 -9.287 0.000
B:st-B:r -63.67 -72.38 -54.954 0.000
S:st-B:r -61.33 -70.05 -52.620 0.000
B:vl-B:r -71.67 -80.38 -62.954 0.000
S:vl-B:r -87.33 -96.05 -78.620 0.000
B:st-S:r -45.67 -54.38 -36.954 0.000
S:st-S:r -43.33 -52.05 -34.620 0.000
B:vl-S:r -53.67 -62.38 -44.954 0.000
S:vl-S:r -69.33 -78.05 -60.620 0.000
S:st-B:st 2.33 -6.38 11.046 0.979
B:vl-B:st -8.00 -16.71 0.713 0.084
S:vl-B:st -23.67 -32.38 -14.954 0.000
B:vl-S:st -10.33 -19.05 -1.620 0.015
S:vl-S:st -26.00 -34.71 -17.287 0.000
S:vl-B:vl -15.67 -24.38 -6.954 0.000Tensile strength
Easiest way to do Dunnett’s is to convert the design to a one-way:
[1] "B_l" "B_r" "B_st" "B_vl" "S_l" "S_r" "S_st" "S_vl"
Dunnett's test for comparing several treatments with a control :
95% family-wise confidence level
$B_st
diff lwr.ci upr.ci pval
B_l-B_st 32.00 24.64 39.356 2.7e-10 ***
B_r-B_st 63.67 56.31 71.023 < 2e-16 ***
B_vl-B_st -8.00 -15.36 -0.644 0.0303 *
S_l-B_st -4.67 -12.02 2.690 0.3267
S_r-B_st 45.67 38.31 53.023 < 2e-16 ***
S_st-B_st 2.33 -5.02 9.690 0.8882
S_vl-B_st -23.67 -31.02 -16.310 1.0e-07 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Comparing factor level means at fixed level of other factor
Fix Factor A at level \(i\) and compare means across factor B:
Tukey’s for comparing all pairs among \(\bar \mu_{i1},\dots,\bar \mu_{ib}\): \[ \bar Y_{ij.} - \bar Y_{ij'.} \pm q_{b,ab(n-1),\alpha} ~ \hat \sigma \frac{1}{\sqrt{n}}, \quad 1 \leq j < j' \leq b. \]
Dunnett’s for comparing \(\bar \mu_{ij},\dots,\bar \mu_{i,b-1}\) to baseline \(\mu_{ib}\): \[ \bar Y_{ij.} - \bar Y_{ib.} \pm d_{b,ab(n-1),\alpha} ~ \hat \sigma \sqrt{\frac{2}{n}}, \quad j = 1,\dots,b-1. \]
To do this at all levels \(i=1,\dots,a\), divide \(\alpha\) by \(a\) (à la Bonferroni)!
Likewise if comparing means of factor B at fixed levels of factor A.
Tensile strength
For each aggregate type, compare all pairs of compaction method means.
Ensure that the familywise coverage probability is \(0.95\).
Use \(q_{4,16,0.05/2} =\) qtukey(1-0.05/2,4,16) \(= 4.54763\).
y11. <- mean(y[agg == "B" & comp == "l"])
y12. <- mean(y[agg == "B" & comp == "r"])
y13. <- mean(y[agg == "B" & comp == "st"])
y14. <- mean(y[agg == "B" & comp == "vl"])
y21. <- mean(y[agg == "S" & comp == "l"])
y22. <- mean(y[agg == "S" & comp == "r"])
y23. <- mean(y[agg == "S" & comp == "st"])
y24. <- mean(y[agg == "S" & comp == "vl"])
alpha <- 0.05
me <- qtukey(1-alpha/a,b,a*b*(n-1)) * sqrt(MSE) / sqrt(n)Tensile strength
ttab <- rbind(c(y11. - y12. - me, y11. - y12. + me),
c(y11. - y13. - me, y11. - y13. + me),
c(y11. - y14. - me, y11. - y14. + me),
c(y12. - y13. - me, y12. - y13. + me),
c(y12. - y14. - me, y12. - y14. + me),
c(y13. - y14. - me, y13. - y14. + me),
c(y21. - y22. - me, y21. - y22. + me),
c(y21. - y23. - me, y21. - y23. + me),
c(y21. - y24. - me, y21. - y24. + me),
c(y22. - y23. - me, y22. - y23. + me),
c(y22. - y24. - me, y22. - y24. + me),
c(y23. - y24. - me, y23. - y24. + me))
rownames(ttab) <- c("B:l-r","B:l-st","B:l-vl","B:r-st","B:r-vl","B:st-vl",
"S:l-r","S:l-st","S:l-vl","S:r-st","S:r-vl","S:st-vl")
colnames(ttab) <- c("lower","upper")Tensile strength
Tukey/Bonferroni-adjusted confidence intervals comparing all pairs of compaction methods for each aggregate type:
Tensile strength
Identify which differences are not significant in the interaction plot:
Tensile strength
For only the basalt aggregate type, compare all compaction method means to that of the static method.
Use \(d_{4,16,0.05} = 2.59\) from the table on the next slide.
Table A.5 from Mohr, Wilson, and Freund (2021)
Tensile strength
Dunnett’s comparison of compaction method means to the static method when the aggregate type is basalt:
Interaction not significant. Then what?
If the interaction is not significant:
- We can interpret main effects.
- We can make meaningful comparisons among marginal means.
Comparing marginal means in the absence of interaction
For making comparisons among the marginal means of Factor A:
Tukey’s for comparing all pairs among \(\bar \mu_{1.},\dots,\bar \mu_{a.}\): \[ \bar Y_{i..} - \bar Y_{i'..} \pm q_{a,ab(n-1),\alpha} ~ \hat \sigma \frac{1}{\sqrt{bn}}, \quad 1 \leq i < i' \leq a. \]
Dunnett’s for comparing \(\bar \mu_{1.},\dots,\bar \mu_{a-1.}\) to a control mean \(\bar \mu_{a.}\): \[ \bar Y_{i..} - \bar Y_{a..} \pm d_{a,ab(n-1),\alpha} ~ \hat \sigma\sqrt{\frac{2}{bn}}, \quad i = 1,\dots, a-1. \]
Still use \(\hat \sigma = \sqrt{\operatorname{MS}_{\operatorname{Error}}}\).
Do likewise for making comparisons among \(\bar \mu_{.1},\dots,\bar \mu_{.b}\) of Factor B.
Serum glucose example from Kuehl (2000)
Two methods for measuring serum glucose level at three glucose levels.
- Is there an interaction between the method and the glucose level?
- If not, can we describe the main effect of the method?
Serum glucose
Serum glucose
Variance appears smaller at lower glucose level. Try using \(\log(Y_{ijk})\).
Serum glucose
This looks better.
Serum glucose
Normality check looks okay.
Serum glucose
Analysis of Variance Table
Response: log(y)
Df Sum Sq Mean Sq F value Pr(>F)
method 1 0.01 0.01 78.11 1.3e-06 ***
glc 2 7.19 3.60 19670.48 < 2e-16 ***
method:glc 2 0.00 0.00 3.06 0.085 .
Residuals 12 0.00 0.00
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1There is only weak evidence of interaction. Check interaction plot.
Serum glucose
Serum glucose
It appears safe to ignore the interaction and report on main effects.
Serum glucose
Call the method Factor A; build a CI for \(\bar \mu_{1.} - \bar \mu_{2.}\) (just one comparison).
Since \(a = 2\), \(b = 3\), and \(n = 3\), use \(\displaystyle \bar Y_{1..} - \bar Y_{2..} \pm q_{2,2\cdot 3(3-1),0.05} ~ \hat \sigma \frac{1}{\sqrt{3\cdot 3}}\).
[1] 0.0424 0.0702Since \(a=2\), \(q_{a,ab(n-1),\alpha} = \sqrt{2} \cdot t_{ab(n-1),\alpha/2}\), so it is just a \(t\)-interval.
Possible workflow for factorial experiments
flowchart TD
A[Overall F test of significance] -->|Fail to reject null| B(Stop)
A[Overall F test of significance] -->|Reject null| C(Test for interaction effect)
C -->|Fail to reject null| E[Test each main effect]
C -->|Reject null| D[Compare factor level means at each level of the other factor]
E -->|Fail to reject null| G[Stop]
E -->|Reject null| F[Compare factor level means averaging over levels of the other factor]
Three-way factorial experiment
Suppose we have factors A, B, and C with \(a\), \(b\), and \(c\) levels, so that \[ Y_{ijkl} = \mu_{ijk} + \varepsilon_{ijkl} \] for \(i = 1,\dots,a\), \(j = 1,\dots,b\), \(k = 1,\dots,c\), \(l=1,\dots,n_{ij}\).
Treatment effects representation of the mean \(\mu_{ijk}\) is \[ \mu_{ijk} = \tau_i + \gamma_j + \eta_k + (\tau\gamma)_{ij} + (\tau\eta)_{ik} + (\gamma\eta)_{jk} + (\tau\gamma\eta)_{ijk}. \]
- \(Y_{ijkl}\) is response for EU \(k\) under levels \(i,j,k\) of factors A, B, and C.
- \(\mu\) represents a baseline or overall mean
- The \(\tau_i\), \(\gamma_j\), and \(\eta_k\) are the main effects for Factors A, B, and C.
- The \((\tau\gamma)_{ij}\), \((\tau\eta)_{ik}\), and \((\gamma\eta)_{jk}\) are two-way interaction effects between the pairs of factors AB, AC, and BC.
- The \((\tau\gamma\eta)_{ijk}\) are three-way interaction effects between A, B, and C.
- The \(\varepsilon_{ijkl}\) are \(\mathcal{N}(0,\sigma^2)\) error terms.
Shrimp data example
Four-week weight gain of shrimp in different aquarium settings:
Temperature (Celcius), Density (# shrimp per 40 liters), and Salinity.
Data from Kuehl (2000), image from Dr. Longnecker’s notes.
Shrimp
wg <- c(86,52,73,544,371,482,390,290,397,
53,73,86,393,398,208,249,265,243,
439,436,349,249,245,330,247,277,205,
324,305,364,352,267,316,188,223,281)
temp <- as.factor(c(rep(25,18),rep(35,18)))
dens <- as.factor(c(1,1) %x% c(rep(80,9),rep(160,9))) # kronecker product
salt <- as.factor(c(1,1,1,1) %x% c(rep(10,3),rep(25,3),rep(40,3)))
aggregate(wg ~ temp + dens + salt, FUN = mean) temp dens salt wg
1 25 80 10 70.3
2 35 80 10 408.0
3 25 160 10 70.7
4 35 160 10 331.0
5 25 80 25 465.7
6 35 80 25 274.7
7 25 160 25 333.0
8 35 160 25 311.7
9 25 80 40 359.0
10 35 80 40 243.0
11 25 160 40 252.3
12 35 160 40 230.7Shrimp
Shrimp
temp dens wg
1 25 80 298
2 35 80 309
3 25 160 219
4 35 160 291 temp salt wg
1 25 10 70.5
2 35 10 369.5
3 25 25 399.3
4 35 25 293.2
5 25 40 305.7
6 35 40 236.8 dens salt wg
1 80 10 239
2 160 10 201
3 80 25 370
4 160 25 322
5 80 40 301
6 160 40 242Shrimp
Shrimp
Shrimp
ANOVA table for balanced three-way factorial design
| Source | Df | SS | MS | F value |
|---|---|---|---|---|
| A | \(a-1\) | \(\operatorname{SS}_{\operatorname{A}}\) | \(\operatorname{MS}_{\operatorname{A}}\) | \(F_{\operatorname{A}} = \operatorname{MS}_{\operatorname{A}}/ \operatorname{MS}_{\operatorname{Error}}\) |
| B | \(b-1\) | \(\operatorname{SS}_{\operatorname{B}}\) | \(\operatorname{MS}_{\operatorname{B}}\) | \(F_{\operatorname{B}} = \operatorname{MS}_{\operatorname{B}}/\operatorname{MS}_{\operatorname{Error}}\) |
| C | \(c-1\) | \(\operatorname{SS}_{\operatorname{C}}\) | \(\operatorname{MS}_{\operatorname{C}}\) | \(F_{\operatorname{C}} = \operatorname{MS}_{\operatorname{C}}/\operatorname{MS}_{\operatorname{Error}}\) |
| AB | \((a-1)(b-1)\) | \(\operatorname{SS}_{\operatorname{AB}}\) | \(\operatorname{MS}_{\operatorname{AB}}\) | \(F_{\operatorname{AB}} = \operatorname{MS}_{\operatorname{AB}}/ \operatorname{MS}_{\operatorname{Error}}\) |
| AC | \((a-1)(c-1)\) | \(\operatorname{SS}_{\operatorname{AC}}\) | \(\operatorname{MS}_{\operatorname{AC}}\) | \(F_{\operatorname{AC}} = \operatorname{MS}_{\operatorname{AC}}/ \operatorname{MS}_{\operatorname{Error}}\) |
| BC | \((b-1)(c-1)\) | \(\SSBC\) | \(\MSBC\) | \(F_{\operatorname{BC}} = \MSBC / \operatorname{MS}_{\operatorname{Error}}\) |
| ABC | \((a-1)(b-1)(c-1)\) | \(\SSABC\) | \(\MSABC\) | \(F_{\operatorname{ABC}} = \MSABC / \operatorname{MS}_{\operatorname{Error}}\) |
| Error | \(abc(n-1)\) | \(\operatorname{SS}_{\operatorname{Error}}\) | \(\operatorname{MS}_{\operatorname{Error}}\) | |
| Total | \(abcn-1\) | \(\operatorname{SS}_{\operatorname{Tot}}\) |
Rejection rules in the three-way factorial design
With reference to the ANOVA table on the preceding slide:
- Reject \(H_0\): no A main effect if \(F_{\operatorname{A}} > F_{a-1,abc(n-1),\alpha}\).
- Reject \(H_0\): no B main effect if \(F_{\operatorname{B}} > F_{b-1,abc(n-1),\alpha}\).
- Reject \(H_0\): no C main effect if \(F_{\operatorname{C}} > F_{c-1,abc(n-1),\alpha}\).
- Reject \(H_0\): no A and B interaction if \(F_{\operatorname{AB}} > F_{(a-1)(b-1),abc(n-1),\alpha}\).
- Reject \(H_0\): no A and C interaction if \(F_{\operatorname{AC}} > F_{(a-1)(c-1),abc(n-1),\alpha}\).
- Reject \(H_0\): no B and C interaction if \(F_{\operatorname{BC}} > F_{(b-1)(c-1),abc(n-1),\alpha}\).
- Reject \(H_0\): no A, B, and C interaction if \(F_{\operatorname{ABC}} > F_{(a-1)(b-1)(c-1),abc(n-1),\alpha}\).
Shrimp
Analysis of Variance Table
Response: wg
Df Sum Sq Mean Sq F value Pr(>F)
temp 1 15376 15376 5.30 0.030 *
dens 1 21219 21219 7.31 0.012 *
salt 2 96762 48381 16.66 2.9e-05 ***
temp:dens 1 8711 8711 3.00 0.096 .
temp:salt 2 300855 150428 51.80 2.0e-09 ***
dens:salt 2 674 337 0.12 0.891
temp:dens:salt 2 24038 12019 4.14 0.029 *
Residuals 24 69691 2904
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1References
Kuehl, R. O. 2000. Design of Experiments: Statistical Principles of Research Design and Analysis. Duxbury/Thomson Learning.
Mohr, Donna L, William J Wilson, and Rudolf J Freund. 2021. Statistical Methods. Academic Press.