Two-Factor Between-Subjects ANOVA

Author

Dr. Anisha Singh

We continue working with readtime.csv, but now we add a second independent variable. In addition to Device (Phone, Tablet, Laptop), participants read under different lighting conditions (Blue, White, Yellow, Natural). Each participant was assigned to exactly one device type and one lighting condition, making this a 3 × 4 between-subjects factorial design (12 groups, 5 participants per group, N = 60).

Note

Why use a factorial design?

Instead of running two separate one-way ANOVAs — one for Device, one for Light — a factorial ANOVA tests both factors simultaneously. This is more efficient, uses the full sample for each test, and most importantly lets us test whether the two factors interact: does the effect of device type depend on what lighting is used? A one-way ANOVA can never answer that question.

Packages

# install.packages("psych")
# install.packages("emmeans")
# install.packages("effectsize")

library(psych)
library(emmeans)
library(effectsize)

Load and Prepare the Data

reader <- read.csv("readtime.csv")

Recoding factor levels

The variables Device and Light are stored as numbers. We need to tell R they are categorical and give them meaningful labels. The cleanest way is with factor(), which sets levels and labels in a single line — no extra packages needed.

reader$Device <- factor(reader$Device,
                        levels = 1:3,
                        labels = c("Phone", "Tablet", "Laptop"))

reader$Light  <- factor(reader$Light,
                        levels = 1:4,
                        labels = c("Blue", "White", "Yellow", "Natural"))

Exploring the Data

Marginal means by Device

describeBy(reader$ReadingTime, reader$Device)

group: Phone
   vars  n    mean     sd  median  min     max   range
X1    1 20 1252.83 305.62 1233.12  572.08 1797.21 1225.13

group: Tablet
   vars  n    mean     sd  median   min     max   range
X1    1 20 1032.20 307.44 1085.92  551.34 1588.64 1037.30

group: Laptop
   vars  n    mean     sd  median   min     max  range
X1    1 20 1014.00 288.58  969.52  543.75 1530.75  987.00

Marginal means by Light

describeBy(reader$ReadingTime, reader$Light)

group: Blue
   vars  n   mean     sd  median   min     max   range
X1    1 15 1306.90 276.55 1295.44  862.69 1797.21  934.52

group: White
   vars  n    mean     sd  median   min     max   range
X1    1 15 1194.90 276.74 1233.46  741.65 1667.19  925.54

group: Yellow
   vars  n   mean     sd  median   min     max   range
X1    1 15  970.70 279.57  991.76  543.75 1494.78  951.03

group: Natural
   vars  n   mean     sd  median   min     max   range
X1    1 15  926.20 285.55  894.61  551.34 1304.63  753.29

Cell means (Device × Light)

tapply(reader$ReadingTime, list(reader$Device, reader$Light), mean)

           Blue    White   Yellow  Natural
Phone  1462.158 1385.960 1093.694 1069.490
Tablet 1250.188 1096.594  914.006  868.010
Laptop 1208.350 1102.146  904.400  841.100

Reading time is highest under blue light for all three devices, and lowest under natural light. Phones are consistently slower than tablets and laptops across all lighting conditions.

Visualizing the Interaction

Before fitting the model, plot the cell means. If the lines are roughly parallel, it suggests no interaction — the effect of one factor is the same regardless of the level of the other.

Quick check: `interaction.plot()`

interaction.plot() is a base R function that requires no additional packages and is useful for a quick look at your data during analysis.

par(mar = c(5, 4, 4, 8))   # widen right margin so the legend is not cut off

interaction.plot(
  x.factor     = reader$Light,
  trace.factor = reader$Device,
  response     = reader$ReadingTime,
  fun          = mean,
  xlab         = "Lighting Condition",
  ylab         = "Mean Reading Time (ms)",
  trace.label  = "Device",
  col          = c("steelblue", "darkorchid", "tomato"),
  lwd          = 2
)

par(mar = c(5, 4, 4, 2))   # reset margins to default

Note

Limitation of interaction.plot()

Because this is a base R function, the legend is drawn outside the plot region and its placement depends on your window size — the device names can be cut off if the plot is too narrow. The overall appearance is also not publication-ready. Use interaction.plot() to get a fast visual check during your analysis, then switch to the ggplot2 version below when you need a clean figure for a paper or presentation.

Publication-ready plot: `ggplot2`

# install.packages("ggplot2")   # run once if needed
# install.packages("dplyr")

library(ggplot2)
library(dplyr)

reader %>%
  group_by(Device, Light) %>%
  summarise(mean_RT = mean(ReadingTime), .groups = "drop") %>%
  ggplot(aes(x = Light, y = mean_RT, color = Device, group = Device)) +
  geom_line(linewidth = 1.2) +
  geom_point(size = 3) +
  labs(
    x     = "Lighting Condition",
    y     = "Mean Reading Time (ms)",
    color = "Device"
  ) +
  theme_minimal()

ggplot2 reads the factor levels directly, so the legend will correctly label each line “Phone”, “Tablet”, and “Laptop” regardless of window size. The plot can be exported at any resolution with ggsave().

The three lines track closely together and run roughly parallel across all four lighting conditions — a visual cue that Device and Light are acting independently rather than interacting.

Fitting the Model

The * operator fits both main effects and their interaction in one step.

fit2.way <- aov(ReadingTime ~ Device * Light, data = reader)
summary(fit2.way)

             Df  Sum Sq Mean Sq F value   Pr(>F)
Device        2  706968  353484   4.648  0.01428 *
Light         3 1481064  493688   6.492  0.00089 ***
Device:Light  6   21543    3591   0.047  0.99953
Residuals    48 3650203   76046

Note

Reading the factorial ANOVA table

There are three tests:

Device main effect: Does reading time differ by device type, averaging across all lighting conditions? F(2, 48) = 4.65, p = .014 — yes.
Light main effect: Does reading time differ by lighting condition, averaging across all devices? F(3, 48) = 6.49, p < .001 — yes.
Device × Light interaction: Does the effect of device change depending on which light is used (or vice versa)? F(6, 48) = 0.047, p = .9995 — no.

The near-zero F for the interaction (0.047) with a p-value of .9995 tells us the two factors act completely independently. The advantage of phones over tablets is essentially the same under blue light, white light, yellow light, and natural light. This is what statisticians call an additive relationship.

Overall Effect Sizes: Eta Squared

Eta squared (η²) is the proportion of total variance in reading time explained by each effect.

eta_squared(fit2.way, partial = FALSE)

# Effect Size for ANOVA (Type I)

Parameter    |  Eta2 |       95% CI
------------------------------------
Device       |  0.12 | [0.00, 1.00]
Light        |  0.25 | [0.07, 1.00]
Device:Light |  0.004| [0.00, 1.00]

Device: η² = .12 — device type explains 12% of the total variance in reading time.
Light: η² = .25 — lighting condition explains 25% of the total variance. This is the stronger of the two main effects.
Interaction: η² < .01 — the Device × Light combination explains essentially none of the variance.

Note

η² vs. partial η²

In a one-way ANOVA, η² and partial η² are identical. In a factorial design they differ: partial η² estimates each effect’s share of variance after removing the other effects, which inflates the values (e.g., Device partial η² = .16 vs. full η² = .12). Full η² is more conservative and directly interpretable as a proportion of total variance. Partial η² is commonly reported in published articles when comparing across studies with different designs.

When the Interaction Is Not Significant

Tip

The interaction guides your next step

In a factorial ANOVA, always look at the interaction first:

If the interaction IS significant: the main effects are misleading on their own. You must decompose the interaction using simple effects — for example, testing Device differences separately within each level of Light. The main effect averages across conditions that behave differently, which obscures what is really happening.
If the interaction is NOT significant (as here): the two factors act independently and additively. It is safe — and correct — to interpret the main effects directly. The effect of Device is the same across all lighting levels, so averaging across them gives an accurate picture.

Because our interaction is not significant, we proceed to main-effect post-hoc tests for Device and Light separately.

If the Interaction Had Been Significant: Simple Effects

When an interaction is significant, interpreting the main effects alone is misleading — the effect of Device depends on which lighting level you are in, so averaging across lighting levels obscures what is actually happening. The correct next step is simple effects analysis: test Device differences separately within each level of Light (and/or Light differences within each level of Device).

# Estimated means for each Device × Light cell
simple <- emmeans(fit2.way, ~ Device | Light)
simple

Light = Blue:
 Device emmean  SE df lower.CL upper.CL
 Phone    1462 123 48     1214     1710
 Tablet   1250 123 48     1002     1498
 Laptop   1208 123 48      960     1456

Light = White:
 Device emmean  SE df lower.CL upper.CL
 Phone    1386 123 48     1138     1634
 Tablet   1097 123 48      849     1345
 Laptop   1102 123 48      854     1350

Light = Yellow:
 Device emmean  SE df lower.CL upper.CL
 Phone    1094 123 48      846     1342
 Tablet    914 123 48      666     1162
 Laptop    904 123 48      656     1152

Light = Natural:
 Device emmean  SE df lower.CL upper.CL
 Phone    1069 123 48      822     1317
 Tablet    868 123 48      620     1116
 Laptop    841 123 48      593     1089

# Pairwise Device comparisons within each Light level
pairs(simple, adjust = "tukey")

Light = Blue:
 contrast        estimate  SE df t.ratio p.value
 Phone - Tablet    211.97 174 48   1.215  0.4500
 Phone - Laptop    253.81 174 48   1.455  0.3213
 Tablet - Laptop    41.84 174 48   0.240  0.9688

Light = White:
 contrast        estimate  SE df t.ratio p.value
 Phone - Tablet    289.37 174 48   1.659  0.2314
 Phone - Laptop    283.81 174 48   1.627  0.2441
 Tablet - Laptop    -5.55 174 48  -0.032  0.9994

Light = Yellow:
 contrast        estimate  SE df t.ratio p.value
 Phone - Tablet    179.69 174 48   1.030  0.5616
 Phone - Laptop    189.29 174 48   1.085  0.5277
 Tablet - Laptop     9.61 174 48   0.055  0.9983

Light = Natural:
 contrast        estimate  SE df t.ratio p.value
 Phone - Tablet    201.48 174 48   1.155  0.4854
 Phone - Laptop    228.39 174 48   1.310  0.3969
 Tablet - Laptop    26.91 174 48   0.154  0.9870

P value adjustment: tukey method for comparing a family of 3 estimates

# Confidence intervals for simple effects
confint(pairs(simple, adjust = "tukey"))

Light = Blue:
 contrast        estimate  SE df lower.CL upper.CL
 Phone - Tablet    211.97 174 48     -210      634
 Phone - Laptop    253.81 174 48     -168      676
 Tablet - Laptop    41.84 174 48     -380      464

Light = White:
 contrast        estimate  SE df lower.CL upper.CL
 Phone - Tablet    289.37 174 48     -132      711
 Phone - Laptop    283.81 174 48     -138      706
 Tablet - Laptop    -5.55 174 48     -427      416

Light = Yellow:
 contrast        estimate  SE df lower.CL upper.CL
 Phone - Tablet    179.69 174 48     -242      601
 Phone - Laptop    189.29 174 48     -233      611
 Tablet - Laptop     9.61 174 48     -412      431

Light = Natural:
 contrast        estimate  SE df lower.CL upper.CL
 Phone - Tablet    201.48 174 48     -220      623
 Phone - Laptop    228.39 174 48     -193      650
 Tablet - Laptop    26.91 174 48     -395      449

Note

Why are the simple effects CIs wider than the main-effect CIs?

When testing Device at each level of Light separately, each comparison uses only the 5 participants per cell rather than the 20 per device collapsed across all lighting conditions. Smaller cell n means larger standard errors and wider confidence intervals. This is the cost of decomposing an interaction — you gain specificity (knowing where the effect occurs) but lose precision.

Post-hoc Tests

Device: Pairwise Comparisons

em_device <- emmeans(fit2.way, ~ Device)
device_pairs <- pairs(em_device, adjust = "tukey")
device_pairs

 contrast        estimate   SE df t.ratio p.value
 Phone - Tablet     220.6 87.2 48   2.530  0.0385
 Phone - Laptop     238.8 87.2 48   2.739  0.0231
 Tablet - Laptop     18.2 87.2 48   0.209  0.9763

Results are averaged over the levels of: Light
P value adjustment: tukey method for comparing a family of 3 estimates

Phone users read significantly slower than both Tablet (p = .039) and Laptop users (p = .023). Tablet and Laptop do not differ from each other (p = .976).

Confidence Intervals: Device

confint(device_pairs)

 contrast        estimate   SE df lower.CL upper.CL
 Phone - Tablet     220.6 87.2 48      9.7      432
 Phone - Laptop     238.8 87.2 48     27.9      450
 Tablet - Laptop     18.2 87.2 48   -192.7      229

Confidence level used: 0.95
Conf-level adjustment: tukey method for comparing a family of 3 estimates

Both Phone − Tablet [9.7, 432] and Phone − Laptop [27.9, 450] have confidence intervals that exclude 0, confirming significant differences. The Tablet − Laptop interval [−192.7, 229] includes 0.

Light: Pairwise Comparisons

em_light <- emmeans(fit2.way, ~ Light)
light_pairs <- pairs(em_light, adjust = "tukey")
light_pairs

 contrast         estimate  SE df t.ratio p.value
 Blue - White        112.0 101 48   1.112  0.6838
 Blue - Yellow       336.2 101 48   3.339  0.0086
 Blue - Natural      380.7 101 48   3.781  0.0024
 White - Yellow      224.2 101 48   2.227  0.1307
 White - Natural     268.7 101 48   2.668  0.0492
 Yellow - Natural     44.5 101 48   0.442  0.9709

Results are averaged over the levels of: Device
P value adjustment: tukey method for comparing a family of 4 estimates

Blue light produced significantly longer reading times than Yellow (p = .009) and Natural (p = .002). White light was significantly slower than Natural (p = .049). Blue vs. White, White vs. Yellow, and Yellow vs. Natural did not reach significance after Tukey correction.

Confidence Intervals: Light

confint(light_pairs)

 contrast         estimate  SE df lower.CL upper.CL
 Blue - White        112.0 101 48   -156.0      380
 Blue - Yellow       336.2 101 48     68.2      604
 Blue - Natural      380.7 101 48    112.7      649
 White - Yellow      224.2 101 48    -43.8      492
 White - Natural     268.7 101 48      0.7      537
 Yellow - Natural     44.5 101 48   -223.5      312

Confidence level used: 0.95
Conf-level adjustment: tukey method for comparing a family of 4 estimates

The Blue − Yellow [68.2, 604] and Blue − Natural [112.7, 649] intervals exclude 0 comfortably. The White − Natural interval [0.7, 537] just clears 0, consistent with its marginal p = .049.

Effect Sizes for Pairwise Comparisons

Cohen’s d quantifies how large each mean difference is in standard deviation units, independent of sample size. We use the pooled model error (sigma(fit2.way)) as the standardizer.

Device pairwise d

eff_size(em_device, sigma = sigma(fit2.way), edf = df.residual(fit2.way))

 contrast          effect.size    SE df lower.CL upper.CL
 (Phone - Tablet)        0.800 0.327 48    0.143    1.457
 (Phone - Laptop)        0.866 0.328 48    0.206    1.526
 (Tablet - Laptop)       0.066 0.316 48   -0.570    0.702

Results are averaged over the levels of: Light
sigma used for effect sizes: 275.8

Phone vs. Tablet: d = 0.80 (large). Phone vs. Laptop: d = 0.87 (large). Tablet vs. Laptop: d = 0.07 (negligible).

Light pairwise d

eff_size(em_light, sigma = sigma(fit2.way), edf = df.residual(fit2.way))

 contrast           effect.size    SE df lower.CL upper.CL
 (Blue - White)           0.406 0.367 48   -0.333    1.145
 (Blue - Yellow)          1.219 0.386 48    0.444    1.995
 (Blue - Natural)         1.381 0.391 48    0.594    2.167
 (White - Yellow)         0.813 0.374 48    0.060    1.566
 (White - Natural)        0.974 0.378 48    0.214    1.735
 (Yellow - Natural)       0.161 0.366 48   -0.574    0.896

Results are averaged over the levels of: Device
sigma used for effect sizes: 275.8

The largest effects involve blue light: Blue vs. Natural d = 1.38 (very large) and Blue vs. Yellow d = 1.22 (very large). White vs. Natural d = 0.97 (large) and White vs. Yellow d = 0.81 (large). Comparisons involving only adjacent conditions (Blue vs. White; Yellow vs. Natural) show small, non-significant effects.

Packages

Load and Prepare the Data

Recoding factor levels

Exploring the Data

Marginal means by Device

Marginal means by Light

Cell means (Device × Light)

Visualizing the Interaction

Quick check: interaction.plot()

Publication-ready plot: ggplot2

Fitting the Model

Overall Effect Sizes: Eta Squared

When the Interaction Is Not Significant

If the Interaction Had Been Significant: Simple Effects

Post-hoc Tests

Device: Pairwise Comparisons

Confidence Intervals: Device

Light: Pairwise Comparisons

Confidence Intervals: Light

Effect Sizes for Pairwise Comparisons

Device pairwise d

Light pairwise d

Quick check: `interaction.plot()`

Publication-ready plot: `ggplot2`