One-Way ANOVA

Introduction: Analysis of Variance

• We have previously discussed testing the difference between two groups.

  • What about when there are three or more groups?

• We will use a method called analysis of variance (ANOVA).

  • This method partitions the variance of the outcome into variance due to the groups and variance due to “other” factors.

• Fun fact: the two-sample t-test is a special case of ANOVA.

  • If you perform ANOVA when comparing two means, you will obtain the same results as the two-sample t-test.

R Setup

library(tidyverse)
library(ssstats)
apple_harvest <- read_csv("https://raw.githubusercontent.com/samanthaseals/STA4173Fa26/refs/heads/main/lectures/data/apple_harvest.csv")
baltimare_harbor <- read_csv("https://raw.githubusercontent.com/samanthaseals/STA4173Fa26/refs/heads/main/lectures/data/baltimare_harbor.csv")
boutique_gems <- read_csv("https://raw.githubusercontent.com/samanthaseals/STA4173Fa26/refs/heads/main/lectures/data/boutique_gems.csv")
klugetown_bazaar <- read_csv("https://raw.githubusercontent.com/samanthaseals/STA4173Fa26/refs/heads/main/lectures/data/klugetown_bazaar.csv")
ponyville_diner_tips <- read_csv("https://raw.githubusercontent.com/samanthaseals/STA4173Fa26/refs/heads/main/lectures/data/ponyville_diner_tips.csv")
wonderbolts_laptimes <- read_csv("https://raw.githubusercontent.com/samanthaseals/STA4173Fa26/refs/heads/main/lectures/data/wonderbolts_laptimes.csv")

Analysis of Variance

• Consider comparing k groups.

ANOVA Tables

• We will see our ANOVA results displayed in a table.

• Here is what we are working towards:

ANOVA: Sums of Squares

• We are partitioning the total variability into different components.

• The total sum of squares measures the overall variability in the data:

\text{SS}_{\text{Tot}} = \sum_{i=1}^k \sum_{j=1}^n (y_{ij} - \bar{y}_{..})^2

• and is partitioned into SS_Trt and SS_E:

\text{SS}_\text{Trt}=n \sum_{i=1}^k \left( \bar{y}_{i.}-\bar{y}_{..} \right)^2 \ \ \ \ \ \ \ \ \ \ \text{SS}_\text{E} = \sum_{i=1}^k \sum_{j=1}^n (y_{ij} - \bar{y}_{i.})^2

ANOVA: Degrees of Freedom

• Next, we have the degrees of freedom. Each sum of square will have its own degrees of freedom.

• The total df:

\text{df}_{\text{Tot}} = n-1

• and is partitioned into df_Trt and df_E:

\text{df}_\text{Trt}=k - 1 \ \ \ \ \ \ \ \ \ \ \text{df}_\text{E} = n - k

ANOVA: Mean Squares and Test Statistic

• To find the mean squares, we simply divide the sum of squares by the appropriate degrees of freedom,

\text{MS}_\text{Trt}=\frac{\text{SS}_{\text{Trt}}}{\text{df}_{\text{Trt}}} \ \ \ \ \ \ \ \ \ \ \text{MS}_\text{E}=\frac{\text{SS}_{\text{E}}}{\text{df}_{\text{E}}}

• Then, we can find our test statistic,

F_{\text{Trt}} = \frac{\text{MS}_{\text{Trt}}}{\text{MS}_{\text{E}}}

ANOVA: Mean Squares and Test Statistic

• Note for the future: the MS and F can be generalized –

\text{MS}_\text{X}=\frac{\text{SS}_{\text{X}}}{\text{df}_{\text{X}}} \ \ \ \ \ \ \ \ \ \ F_{\text{X}} = \frac{\text{MS}_{\text{X}}}{\text{MS}_{\text{E}}}

• When constructing a test statistic in ANOVA, we always divide by MS_E.

ANOVA Tables

• Revisiting the ANOVA table,

• Note the relationships:

  • SS_Trt + SS_E = SS_Tot
  • df_Trt + df_E = df_Tot

Hypothesis Test for \mu_1 = \mu_2 = ... = \mu_k

• Hypotheses

  • H_0: \mu_1 = \mu_2 = ... = \mu_k
  • H_1: at least one mean is different.

• Test Statistic:

  • F_0 = \text{MS}_{\text{Trt}} / \text{MS}_{\text{E}}

Statistical Inference on k Means (R)

• We will use the anovas() function from library(ssstats) to construct the ANOVA table.

1dataset_name |> anovas(outcome = continuous_variable,
2                       grouping = grouping_variable,
3                       alpha = alpha_level)

Example 1: Sweet Apple Acres

• Applejack has come to us: “Between the Red Delicious row, the Honeycrisp row, and that new Zap Apple row Granny Smith insisted on plantin’, I reckon one of ‘em oughta be doin’ better than the others. Can you tell me if they’re really any different, or am I just imaginin’ things?”

• We will use the apple_harvest dataset to answer the following question:

  • Does apple harvest (in bushels) differ across the three orchard rows?

Example 1

• What are the levels in our factor (row)?

apple_harvest |> n_pct(row)

Example 1

• What are the summary statistics by level of row?

apple_harvest |> group_by(row) |> mean_median(bushels)

# A tibble: 3 × 4
  row               variable mean_sd    median_iqr
  <chr>             <chr>    <chr>      <chr>     
1 Honeycrisp Row    bushels  41.6 (8.0) 41.0 (6.9)
2 Red Delicious Row bushels  40.6 (3.8) 42.0 (5.1)
3 Zap Apple Row     bushels  37.5 (4.4) 38.5 (5.1)

Example 1

• Constructing the ANOVA table,

apple_harvest |> anovas(outcome = bushels,
                        grouping = row)

One-way ANOVA

Source     df         SS       MS      F  p-value 
──────────────────────────────────────────── 
row         2   160.9693  80.4846  2.468    0.095 
Residuals  51  1663.4383  32.6164                 
──────────────────────────────────────────── 


Hypotheses (row):
H₀: All group means are equal (no effect of row)
H₁: At least one group mean differs (effect of row exists)
Test statistic: F(2, 51) = 2.468
p-value: 0.095
Conclusion: Fail to reject the null hypothesis (p = 0.095 ≥ α = 0.050)

Example 1

• Hypotheses

  • H_0: \mu_{\text{HC}} = \mu_{\text{RD}} = \mu_{\text{ZA}}
  • H_1: at least one mean is different.

• Test Statistic and p-Value:

  • F_0 = 2.468
  • p = 0.095

• Rejection Region

  • Reject H_0 if p < \alpha; \alpha = 0.05

• Conclusion / Interpretation

  • Fail to reject H_0. There is not sufficient evidence to suggest that the production of the rows are different.

Example 2: Rarity’s Boutique

• Rarity flits between racks of fabric: “Darling, I need your professional opinion. Between my Canterlot Couture line, my everyday Ponyville Casual line, and my absolutely dazzling Gala Elegance line, I use wildly different amounts of gemstones. Is the difference really as dramatic as I think?”

• We will use the boutique_gems dataset to answer the following question:

  • Does the number of gems used per dress differ across the three collections?

Example 2

• What are the levels in our factor (collection)?

boutique_gems |> n_pct(collection)

Example 2

• What are the summary statistics by level of collection?

boutique_gems |> group_by(collection) |> mean_median(gem_count)

# A tibble: 3 × 4
  collection        variable  mean_sd     median_iqr 
  <chr>             <chr>     <chr>       <chr>      
1 Canterlot Couture gem_count 55.9 (8.4)  57.0 (8.0) 
2 Gala Elegance     gem_count 124.3 (7.6) 124.0 (8.5)
3 Ponyville Casual  gem_count 87.9 (6.4)  87.0 (9.0) 

Example 2

• Constructing the ANOVA table,

boutique_gems |> anovas(outcome = gem_count,
                       grouping = collection)

One-way ANOVA

Source      df          SS          MS        F  p-value 
─────────────────────────────────────────────────── 
collection   2  35140.5778  17570.2889  311.110  < 0.001 
Residuals   42   2372.0000     56.4762                   
─────────────────────────────────────────────────── 


Hypotheses (collection):
H₀: All group means are equal (no effect of collection)
H₁: At least one group mean differs (effect of collection exists)
Test statistic: F(2, 42) = 311.11
p-value: < 0.001
Conclusion: Reject the null hypothesis (p = < 0.001 < α = 0.050)

Example 2

• Hypotheses

  • H_0: \mu_{\text{CC}} = \mu_{\text{PC}} = \mu_{\text{GE}}
  • H_1: at least one mean is different.

• Test Statistic and p-Value:

  • F_0 = 311.11
  • p < 0.001

• Rejection Region

  • Reject H_0 if p < \alpha; \alpha = 0.05

• Conclusion / Interpretation

  • Reject H_0. There is sufficient evidence to suggest that the number of gems are different between the fashion lines.

Example 3: Ponyville Diner

• Sunny Skies wipes down the counter after the lunch rush: “Feels like Caramel always walks away with better tips than me and Lyra put together. Before I start complainin’ to the manager, can you check if there’s actually a difference in how much we’re bringing in?”

• We will use the ponyville_diner_tips dataset to answer the following question:

  • Does tip amount differ across the three servers?

Example 3

• What are the levels in our factor (server)?

ponyville_diner_tips |> n_pct(server)

Example 3

• What are the summary statistics by level of server?

ponyville_diner_tips |> group_by(server) |> mean_median(tip_bits)

# A tibble: 3 × 4
  server            variable mean_sd   median_iqr
  <chr>             <chr>    <chr>     <chr>     
1 Caramel           tip_bits 6.7 (4.4) 6.0 (5.9) 
2 Lyra Heartstrings tip_bits 4.4 (4.7) 2.9 (2.9) 
3 Sunny Skies       tip_bits 7.4 (8.9) 4.3 (8.5) 

Example 3

• Constructing the ANOVA table,

ponyville_diner_tips |> anovas(outcome = tip_bits,
                       grouping = server)

One-way ANOVA

Source     df         SS       MS      F  p-value 
──────────────────────────────────────────── 
server      2   109.2923  54.6461  1.356    0.265 
Residuals  63  2539.3768  40.3076                 
──────────────────────────────────────────── 


Hypotheses (server):
H₀: All group means are equal (no effect of server)
H₁: At least one group mean differs (effect of server exists)
Test statistic: F(2, 63) = 1.356
p-value: 0.265
Conclusion: Fail to reject the null hypothesis (p = 0.265 ≥ α = 0.050)

Example 3

• Hypotheses

  • H_0: \mu_{\text{C}} = \mu_{\text{LH}} = \mu_{\text{SS}}
  • H_1: at least one mean is different.

• Test Statistic and p-Value:

  • F_0 = 1.356
  • p = 0.265

• Rejection Region

  • Reject H_0 if p < \alpha; \alpha = 0.05

• Conclusion / Interpretation

  • Fail to reject H_0. There is not sufficient evidence to suggest that the servers are earning different tips.

ANOVA Assumptions

• We can represent ANOVA with the following model:

y_{ij} = \mu + \tau_i + \varepsilon_{ij}

• where:

  • y_{ij} is the j^{\text{th}} observation in the i^{\text{th}} group,
  • \mu is the overall (grand) mean,
  • \tau_i is the treatment effect for group i, and
  • \varepsilon_{ij} is the error term for the j^{\text{th}} observation in the i^{\text{th}} group.

ANOVA Assumptions

• We assume that the error term follows a normal distribution with mean 0 and a constant variance, \sigma^2. i.e.,

\varepsilon_{ij} \overset{\text{iid}}{\sim} N(0, \sigma^2)

• Very important note: the assumption is on the error term (\varepsilon_{ij}) and NOT on the outcome (y_{ij})!

ANOVA Assumptions

• We will use the residual (the difference between the observed value and the predicted value) to assess assumptions:

e_{ij} = y_{ij} - \hat{y}_{ij}

• Note that in the case of one-way ANOVA, the predicted value of y is the group mean.

ANOVA Assumptions: Graphical Assessment

• Normality: quantile-quantile plot

  • Should have points close to the 45^\circ line
  • We will focus on the “center” portion of the plot

• Variance: scatterplot of the residuals against the predicted values

  • Should be “equal spread” between the groups
  • No “pattern”

Graphical Assessment of Assumptions (R)

• We will use the anova_check() function from library(ssstats) to assess normality.

1dataset_name |> anova_check(outcome = continuous_variable,
2                       grouping = grouping_variable)

• This will provide the the Q-Q plots and variance scatterplot necessary to assess the ANOVA assumptions.

Example 1

• Checking the normality of our Sweet Apple Acres dataset,

apple_harvest |> anova_check(outcome = bushels,
                             grouping = row)  

Example 2

• Checking the normality of Rarity’s dataset,

boutique_gems |> anova_check(outcome = gem_count,
                             grouping = collection)

Example 3

• Checking the normality of the Ponyville Diner dataset,

ponyville_diner_tips |> anova_check(outcome = tip_bits,
                                    grouping = server)

Kruskal-Wallis

• If we break either ANOVA assumption, we should implement the nonparametric version, the Kruskal-Wallis.

  • The Kruskal-Wallis is an extension of the Wilcoxon rank sum (as ANOVA is an extension of the two-sample t-test).

• The Kruskal-Wallis test determines if k independent samples come from populations with the same distribution.

• Note: we can transform the data to make the variance “behave”

  • Box-Cox, or power, transformation

Hypothesis Test for M_1 = M_2 = ... = M_k

• Hypotheses

  • H_0: M_1 = ... = M_k
  • H_1: at least one M_i is different

Hypothesis Testing: Kruskal-Wallis

• Test Statistic

\chi^2_0 = \frac{12}{n(n+1)} \sum_{i=1}^k \frac{R_i^2}{n_i} - 3(n+1) \sim \chi^2_{\text{df}}

• where

  • R_i is the sum of the ranks for group i,
  • n_i is the sample size for group i,
  • n = \sum_{i=1}^k n_i = total sample size,
  • k is the number of groups, and
  • \text{df} = k-1

Statistical Inference on k Medians (R)

• We will use the kruskal_wallis() function from library(ssstats) to perform the Kruskal-Wallis test.

1dataset_name |> kruskal_wallis(outcome = continuous_variable,
2                               grouping = grouping_variable,
3                               alpha = specified_alpha)

Example 1

• Applying the Kruskal-Wallis to this example,

apple_harvest |> kruskal_wallis(outcome = bushels,
                                grouping = row)

Kruskal-Wallis test

M[Honeycrisp Row] = 41 (n = 18)
M[Red Delicious Row] = 42 (n = 18)
M[Zap Apple Row] = 38.45 (n = 18)

Hypotheses:
H₀: All group distributions are equal (no effect of row)
H₁: At least one group distribution differs in location (effect of row exists)
Test statistic: H(2) = 5.953
p-value: 0.051
Conclusion: Fail to reject the null hypothesis (p = 0.051 ≥ α = 0.050)

Example 1

• Hypotheses

  • H_0: \mu_{\text{HC}} = \mu_{\text{RD}} = \mu_{\text{ZA}}
  • H_1: at least one mean is different.

• Test Statistic and p-Value:

  • F_0 = 5.953
  • p = 0.051

• Rejection Region

  • Reject H_0 if p < \alpha; \alpha = 0.05

• Conclusion / Interpretation

  • Fail to reject H_0. There is not sufficient evidence to suggest that the production of the rows are different.

Example 3

• Applying the Kruskal-Wallis to this example,

ponyville_diner_tips |> kruskal_wallis(outcome = tip_bits,
                       grouping = server)

Kruskal-Wallis test

M[Caramel] = 6.03 (n = 22)
M[Lyra Heartstrings] = 2.94 (n = 22)
M[Sunny Skies] = 4.305 (n = 22)

Hypotheses:
H₀: All group distributions are equal (no effect of server)
H₁: At least one group distribution differs in location (effect of server exists)
Test statistic: H(2) = 4.235
p-value: 0.120
Conclusion: Fail to reject the null hypothesis (p = 0.120 ≥ α = 0.050)

Example 3

• Hypotheses

  • H_0: M_{\text{C}} = M_{\text{LH}} = M_{\text{SS}}
  • H_1: at least one median is different.

• Test Statistic and p-Value:

  • \chi^2_0 = 4.235
  • p = 0.120

• Rejection Region

  • Reject H_0 if p < \alpha; \alpha = 0.05

• Conclusion / Interpretation

  • Fail to reject H_0. There is not sufficient evidence to suggest that the servers are earning different tips.

Let’s Practice! Wonderbolts Academy

• Fleetfoot drops a stack of stopwatch printouts on the table: “Command wants to know if lap times actually differ across our four training squads: Thunderhead, Stormwing, Skybolt, and Downburst.”

• We will use the wonderbolts_laptimes dataset to answer the following question:

  • Do lap times differ across the four training squads?

Let’s Practice!: Klugetown Bazaar

• A weary trader waves you over between haggling matches: “Rustbolt, Glimmering, and Junker each swear their stall pulls in the best prices. I’ve got a season’s worth of sale records — do the three stalls actually sell for different amounts, or is everyone just talkin’ up their own wares?”

• We will use the klugetown_bazaar dataset to answer the following question:

  • Does sale price differ across the three market stalls?

Let’s Practice!: Baltimare Harbor

• The harbormaster spreads a stack of logbooks across the table: “Four dock crews unload every ship that comes through: Anchor, Tidewater, Saltmarsh, and Driftwood. Some of the sailors have been grumbling that one crew is dragging everyone’s schedule behind. Before I say anything to the crews, I want to know if unloading time actually differs across them.”

• We will use the baltimare_harbor dataset to answer the following question:

  • Does unloading time differ across the four dock crews?

Wrap Up

• This lecture introduced one-way ANOVA.

  • One-way ANOVA examines the effect of a single factor.

  • (We will learn about multiple factors in the future!)

• One-way ANOVA only allows us to determine if there is an overall, or omnibus, effect.

  • It does not tell us about pairwise differences.

• In the future:

  • Posthoc testing
  • Two-way ANOVA
  • Profile plots