03 — Two-Way ANOVA

What changes when there are two factors?

Scientific question

So far we have treated all 60 replicates per method as exchangeable.

But they are not.

They come from three different environmental matrices: Fill, Soil, and Sediment.

Ignoring Matrix leaves potentially important structure unmodelled.

Residual variance that could be explained by Matrix remains in the error term, reducing the precision of our method comparisons.

More importantly, if Matrix matters, any conclusion about methods that ignores Matrix may be incomplete or misleading.

The natural next question is therefore:

Does Matrix influence the relative bias, and does it do so independently of Method?

This is the domain of two-way ANOVA.

The dataset

We continue with the same balanced dataset.

The response is RelativeBias_pct.

We now have two factors:

  • Method (4 levels: Reference, A, B, C)
  • Matrix (3 levels: Fill, Soil, Sediment)

The design remains balanced: 20 observations per cell, 240 observations overall.

Load the data

library(here)
library(data.table)
library(ggplot2)

dat <- fread(
  here("data", "dataset_v01_clean_balanced.csv")
)

dat[, Method := factor(
  Method,
  levels = c("Reference", "Method_A", "Method_B", "Method_C")
)]

dat[, Matrix := factor(
  Matrix,
  levels = c("Fill", "Soil", "Sediment")
)]

dat[, .N, by = .(Method, Matrix)]
       Method   Matrix     N
       <fctr>   <fctr> <int>
 1: Reference     Fill    20
 2:  Method_A     Fill    20
 3:  Method_B     Fill    20
 4:  Method_C     Fill    20
 5: Reference     Soil    20
 6:  Method_A     Soil    20
 7:  Method_B     Soil    20
 8:  Method_C     Soil    20
 9: Reference Sediment    20
10:  Method_A Sediment    20
11:  Method_B Sediment    20
12:  Method_C Sediment    20

All 12 cells have 20 observations.

The design is balanced.

This matters.

In a balanced design, the sums of squares for Method and Matrix are orthogonal.

Each factor can be tested after the other, and the order does not matter.

This clean separation will disappear when we introduce unbalanced designs in Chapter 05.

For now, we exploit it to focus on the concepts.

Why add Matrix?

Matrix effects are the daily bread of analytical chemistry.

A method that works well in clean Fill may fail in high-organic Sediment.

Extraction efficiency, interference, and recovery all depend on the sample matrix.

If we ignore Matrix, we are implicitly assuming that the method effect is the same in all matrices.

That assumption may be wrong.

And if it is wrong, the one-way model is incomplete.

Look at the data: cell means

Before fitting any model, we examine the cell means.

cell_means <- dat[, .(
  Mean = mean(RelativeBias_pct),
  SD   = sd(RelativeBias_pct),
  SE   = sd(RelativeBias_pct) / sqrt(.N),
  N    = .N
), by = .(Method, Matrix)]

cell_means
       Method   Matrix      Mean       SD        SE     N
       <fctr>   <fctr>     <num>    <num>     <num> <int>
 1: Reference     Fill  -7.44550 1.719187 0.3844220    20
 2:  Method_A     Fill -14.30650 2.371942 0.5303823    20
 3:  Method_B     Fill   3.42045 3.332940 0.7452681    20
 4:  Method_C     Fill  -0.80340 2.878025 0.6435459    20
 5: Reference     Soil -10.77400 2.797611 0.6255648    20
 6:  Method_A     Soil -17.56500 3.198729 0.7152576    20
 7:  Method_B     Soil   2.48005 3.495915 0.7817105    20
 8:  Method_C     Soil  -1.35020 2.974563 0.6651325    20
 9: Reference Sediment -10.34650 3.826971 0.8557368    20
10:  Method_A Sediment -21.01000 3.417124 0.7640922    20
11:  Method_B Sediment  -0.04235 3.765225 0.8419299    20
12:  Method_C Sediment  -5.64900 3.154135 0.7052860    20

The means range from -21% (Method A in Sediment) to 3.4% (Method B in Fill).

Even within a single method, the bias varies across matrices.

Method A, for example, is approximately -14.3% in Fill but drops to -21% in Sediment.

That is a difference of -6.7 percentage points within a single method.

This pattern suggests that Matrix is not merely a nuisance variable.

It may interact with Method.

But before testing the interaction, we fit the additive model.

The additive model

The additive two-way model is:

\[ Y_{ijk} = \mu + \alpha_i + \beta_j + \varepsilon_{ijk} \]

where:

  • \(\alpha_i\) is the effect of method \(i\);
  • \(\beta_j\) is the effect of matrix \(j\);
  • \(\varepsilon_{ijk} \sim N(0, \sigma^2)\).

This model assumes that the effect of Method is the same in every Matrix.

It estimates an “average” method effect across matrices, and an “average” matrix effect across methods.

The assumption is strong.

We will test it.

What does the model test?

Two separate null hypotheses:

  1. \(H_0: \alpha_1 = \alpha_2 = \alpha_3 = \alpha_4 = 0\) (no method effect)
  2. \(H_0: \beta_1 = \beta_2 = \beta_3 = 0\) (no matrix effect)

Because the design is balanced, these tests are independent of each other.

The sum of squares for Method is the same whether Matrix is in the model or not.

This orthogonality is a property of the design, not of the model.

Fitting the additive model

fit_add <- lm(RelativeBias_pct ~ Method + Matrix, data = dat)
summary(fit_add)

Call:
lm(formula = RelativeBias_pct ~ Method + Matrix, data = dat)

Residuals:
    Min      1Q  Median      3Q     Max 
-8.2998 -2.0300  0.2267  2.0741 10.9446 

Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)     -7.3564     0.5084 -14.470  < 2e-16 ***
MethodMethod_A  -8.1052     0.5871 -13.807  < 2e-16 ***
MethodMethod_B  11.4747     0.5871  19.546  < 2e-16 ***
MethodMethod_C   6.9211     0.5871  11.790  < 2e-16 ***
MatrixSoil      -2.0186     0.5084  -3.970 9.55e-05 ***
MatrixSediment  -4.4782     0.5084  -8.808 2.90e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3.215 on 234 degrees of freedom
Multiple R-squared:  0.852, Adjusted R-squared:  0.8489 
F-statistic: 269.5 on 5 and 234 DF,  p-value: < 2.2e-16

The coefficients are parameterized with Reference and Fill as baselines.

The method coefficients are differences from Reference, averaged across matrices.

The matrix coefficients are differences from Fill, averaged across methods.

This is the additive assumption in action.

The ANOVA table

anova(fit_add)
Analysis of Variance Table

Response: RelativeBias_pct
           Df  Sum Sq Mean Sq F value    Pr(>F)    
Method      3 13127.4  4375.8  423.24 < 2.2e-16 ***
Matrix      2   804.8   402.4   38.92  2.56e-15 ***
Residuals 234  2419.3    10.3                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The output shows:

  • Method: \(F_{3, 234} = 423.2\), \(p < 0.001\).
  • Matrix: \(F_{2, 234} = 38.9\), \(p < 0.001\).

Both effects are highly significant.

Method remains the dominant source of variation.

But Matrix also contributes meaningfully.

The residual degrees of freedom are 234, down from 236 in the one-way model because we have spent 2 additional parameters on Matrix.

How much variation does the model explain?

summary(fit_add)$r.squared
[1] 0.852044

The additive model explains approximately 85.2% of the variance.

This is an improvement over the one-way model, which explained 80.3%.

Adding Matrix has captured real structure.

The interaction model: a preview

The additive model assumes that the method effect is constant across matrices.

We can test that assumption by fitting the full interaction model:

\[ Y_{ijk} = \mu + \alpha_i + \beta_j + (\alpha\beta)_{ij} + \varepsilon_{ijk} \]

fit_int <- lm(RelativeBias_pct ~ Method * Matrix, data = dat)
anova(fit_int)
Analysis of Variance Table

Response: RelativeBias_pct
               Df  Sum Sq Mean Sq  F value    Pr(>F)    
Method          3 13127.4  4375.8 446.7043 < 2.2e-16 ***
Matrix          2   804.8   402.4  41.0777 5.815e-16 ***
Method:Matrix   6   185.9    31.0   3.1624   0.00534 ** 
Residuals     228  2233.4     9.8                       
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The interaction term is significant: \(F_{6, 228} = 3.16\), \(p = 5.3e-03\).

This tells us that the additive model is not fully adequate.

The effect of Method is not the same in every Matrix.

But we defer full interpretation of the interaction to the next chapter.

For now, we note that the additive model is a useful stepping stone.

It tells us that both Method and Matrix matter.

The interaction model tells us that the way they matter is more complex.

Comparing the two models

anova(fit_add, fit_int)
Analysis of Variance Table

Model 1: RelativeBias_pct ~ Method + Matrix
Model 2: RelativeBias_pct ~ Method * Matrix
  Res.Df    RSS Df Sum of Sq      F  Pr(>F)   
1    234 2419.3                               
2    228 2233.4  6    185.87 3.1624 0.00534 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The nested F-test compares the additive and interaction models.

The F-statistic is 3.16 on 6 and 228 degrees of freedom, with \(p = 5.3e-03\).

The interaction model fits significantly better.

This is not surprising given what we saw in the cell means.

Method A deteriorates markedly in Sediment, while the other methods are more stable.

That pattern is exactly what an interaction captures.

What the additive model tells us — and what it does not

What it tells us

  • Method has a strong effect on relative bias, even after accounting for Matrix.
  • Matrix also has a significant effect.
  • Together, Method and Matrix explain approximately 85.2% of the variance.

What it does not tell us

  • Which matrix differs from which? The F-test for Matrix is global. We would need contrasts on Matrix levels to pinpoint differences.
  • Why does Matrix matter? The model detects the effect but does not explain the mechanism. Organic matter content varies across matrices and may be the underlying driver. We will explore this with ANCOVA.
  • Is the interaction large enough to matter? The F-test is significant, but the practical impact depends on the application. If you never analyse Sediment, Method A’s poor performance there is irrelevant.

Chemical relevance

Knowing that Matrix matters is not merely a statistical curiosity.

It is a warning.

Any method validation that tests only one matrix is incomplete.

A method certified on Fill may fail on Sediment.

However, the effect size is modest compared to Method.

In practical terms, choosing the right method matters more than choosing the right matrix.

But you still need to know how your method behaves in your matrix.

Take-home message

Adding a second factor is not just “more ANOVA”.

It changes the question from:

“Do methods differ?”

to:

“Do methods differ after accounting for matrix?”

In a balanced design, the answers are clean and orthogonal.

The real world is rarely balanced.

But understanding the balanced case first is essential before tackling the messier ones.

The additive model is a useful approximation.

It tells us that both factors matter.

But the significant interaction tells us that the approximation is incomplete.

The full story requires understanding how Method and Matrix combine.

Next: Interactions

We now know that both Method and Matrix influence the bias.

But the additive model assumes that the method effect is the same in every matrix.

The data tell us otherwise.

Method A deteriorates in Sediment.

The other methods are more stable.

This is an interaction.

And an interaction changes the question.

Next, we ask:

Does the difference between methods change depending on the matrix?