Lecture 04 — Testing Our Hypothesis

The two-sample Welch’s t-test from assumptions to results

tidyverse
descriptive-stats

What are T-Tests and how to run them and test assumptions and interpret.

Author

Bill Perry

Published

July 5, 2026

Where we left off (Lecture 03)

  • Wranglingfilter(), select(), mutate(), arrange() as the core tidyverse verbs
  • Descriptive stats — mean, median, SD, SE with group_by() + summarize()
  • NA handlingsum(!is.na()) is the safe way to count n
  • skimr — fast full-dataset overview in one line
  • Boxplot + jittered points; mean ± SE with stat_summary()
  • Finding: shady leaves appear heavier, taller, and wider than sunny leaves
Note

✅ Key idea from Lecture 03

We described the pattern in our data. Today we formally test whether that pattern is statistically significant.

Goals for Today

  • Understand why comparing means is not enough — we need a test
  • State null and alternate hypotheses formally
  • Check assumptions before running any test
    • normality: histogram, QQ plot, Shapiro-Wilk
    • variance: Levene’s test
  • Run Welch’s two-sample t-test in R
  • Read every line of the t.test() output
  • Report results in scientific writing format
Tip

🖐 Try it yourself

By the end you will run and interpret a complete t-test on our leaf data.

Tools added today:

  • car — Levene’s variance test

References:

Today’s naming:

  • test objects → _model
  • plots → _plot

Our Hypotheses — A Reminder

Biological prediction:

Shady-side leaves will be larger and heavier — they need more surface area to capture the limited light filtering through the canopy.

Formal statement
H₀ — null \(\mu_{shady} = \mu_{sunny}\) — no difference in mean leaf weight
Hₐ — alternate \(\mu_{shady} \neq \mu_{sunny}\) — mean leaf weight differs by side

We use a two-tailed test — we will detect a difference in either direction.

Why two-tailed?

  • We predicted shady > sunny
  • A one-tailed test would only reject H₀ in that direction
  • Two-tailed is more conservative and more broadly accepted
  • If the result is significant two-tailed, it is also significant one-tailed

α = 0.05 — our threshold for rejecting H₀

How to Use These Slides — Predict · Type · Run

This lecture runs in four short chunks. After each chunk you switch to the activity worksheet and type the code yourself.

For every code block, do three things:

  1. Predict — before it runs, say what you think the output will be
  2. Type it out by hand — do not copy-paste
  3. Run it and compare to your prediction
Note

✅ Why bother? (the evidence)

  • Predicting first forces your brain to retrieve what it knows. The gap between your guess and the real answer is what makes the idea stick — vital for slippery ideas like p-values.
  • Typing by hand builds the finger-memory and error-spotting that copy-paste skips.
  • Chunk → immediate practice keeps a new idea in working memory long enough to form a lasting schema.

Load Libraries and Data

# Load all packages at the top ----------------------------
library(readxl) # reading Excel files
library(tidyverse) # data manipulation + ggplot2
library(skimr) # fast summaries
library(car) # Levene's test for variance equality
# Load the leaf data from the data folder -----------------
tree_df <- read_excel("data/2026_06_25_tree_experiment_raw_data.xlsx")
Note

car is the Companion to Applied Regression package. It contains leveneTest(), which we use to check whether our two groups have similar variance.

Install once: install.packages("car")

What We Know So Far

Note

🔮 Predict first: The two group means differ by about 4 g. With only 10 leaves per side — commit now to a yes/no: is that gap big enough to be statistically significant? You’ll find out in Chunk 3.

# Quick recap from Lecture 03 ----------------------------
recap_df <- tree_df %>%
  group_by(side) %>%
  summarize(
    n = sum(!is.na(weight_g)),
    mean_wt = round(mean(weight_g, na.rm = TRUE), 2),
    sd_wt = round(sd(weight_g, na.rm = TRUE), 2),
    se_wt = round(sd_wt / sqrt(n), 2)
  )
recap_df
# A tibble: 2 × 5
  side      n mean_wt sd_wt se_wt
  <chr> <int>   <dbl> <dbl> <dbl>
1 shady    10     7.8  1.03  0.33
2 sunny    10     3.8  0.79  0.25
  • RESULTS:
    • Shady mean: 7.80 g
    • Sunny mean: 3.80 g
    • Difference: 4.0 g — a big gap!
  • But is that gap real, or could it just be sampling noise?
  • With only n = 10 per group, random chance could produce apparent differences.
  • We need a formal test to quantify the probability.

🧩 Chunk 1 of 4 · Why We Test & the t-Statistic

We will cover: why comparing means isn’t enough, what a two-sample t-test is, why we default to Welch’s, and the formula in action.

Tip

🖐 After this chunk: Activity Parts 1–3 (load data, state hypotheses, recap the descriptive stats).

Why Comparing Means Is Not Enough

# Remind ourselves what the data look like ---------------
recap_plot <- tree_df %>%
  ggplot(aes(x = side, y = weight_g, color = side)) +
  geom_point(
    position = position_jitter(width = 0.15, seed = 42),
    alpha = 0.5,
    size = 2
  ) +
  stat_summary(fun = mean, geom = "point", size = 4) +
  stat_summary(
    fun.data = mean_se,
    geom = "errorbar",
    width = 0.15,
    linewidth = 0.9
  ) +
  labs(x = "Side", y = "Leaf Weight (g)") +
  theme_minimal() +
  theme(legend.position = "none")
recap_plot

The means look different — but:

  • We only have 10 leaves per group
  • There is natural variability within each group
  • The groups’ ranges overlap somewhat

A p-value tells us:

“If H₀ were true, how often would we see a gap this large just by chance?”

Small p → the gap is unlikely to be chance.

What Is a Two-Sample t-Test?

A two-sample t-test compares the means of two independent groups.

It asks: Could both groups plausibly be drawn from populations with the same mean?

The t-statistic is:

\[t = \frac{\bar{x}_1 - \bar{x}_2}{\text{SE}_{\text{diff}}}\]

  • Numerator: the observed difference in means
  • Denominator: the uncertainty in that difference (SE of the difference)

Large |t| → the difference is large relative to the noise → small p-value

Component Our context
Group 1 Shady leaves
Group 2 Sunny leaves
Measurement weight_g
H₀ \(\mu_{shady} = \mu_{sunny}\)
Hₐ \(\mu_{shady} \neq \mu_{sunny}\)

Why Welch’s? — Don’t Assume Equal Variance

Standard t-test — pools the two variances into one estimate. Requires both groups to have the same population variance.

Welch’s t-test — uses each group’s variance separately. Does not require equal variance.

\[t_W = \frac{\bar{x}_1 - \bar{x}_2}{\sqrt{\dfrac{s_1^2}{n_1} + \dfrac{s_2^2}{n_2}}}\]

The denominator is the true SE of the difference — no pooling.

Why always use Welch’s?

  • When variances are equal: almost identical result — minimal loss of power
  • When variances are not equal: standard t-test gives inflated Type I error (false positives)
  • Welch’s protects you either way
Tip

Best practice: use Welch’s as your default for all two-group comparisons.

In R: var.equal = FALSE

The Welch’s Formula in Action

Note

🔮 Predict first: t = difference ÷ SE-of-the-difference, and the difference is 4 g. Do you expect |t| to land nearer 1, 5, or 10? Predict before the numbers appear.

# Pull stats and compute t by hand ----------------------
stats_df <- tree_df %>%
  group_by(side) %>%
  summarize(
    n = sum(!is.na(weight_g)),
    m = mean(weight_g, na.rm = TRUE),
    s = sd(weight_g, na.rm = TRUE)
  )

m_sha <- stats_df$m[stats_df$side == "shady"]
m_sun <- stats_df$m[stats_df$side == "sunny"]
s_sha <- stats_df$s[stats_df$side == "shady"]
s_sun <- stats_df$s[stats_df$side == "sunny"]
n_sha <- stats_df$n[stats_df$side == "shady"]
n_sun <- stats_df$n[stats_df$side == "sunny"]

se_diff <- sqrt((s_sha^2 / n_sha) + (s_sun^2 / n_sun))
t_manual <- (m_sha - m_sun) / se_diff

cat("Difference in means:", round(m_sha - m_sun, 2), "g\n")
Difference in means: 4 g
cat("SE of difference:   ", round(se_diff, 4), "\n")
SE of difference:    0.411 
cat("t-statistic:        ", round(t_manual, 3), "\n")
t-statistic:         9.733 

Step by step:

  1. Difference in means: 7.80 − 3.80 = 4.0 g
  2. SE of the difference: \(\sqrt{\frac{s^2_{shady}}{10} + \frac{s^2_{sunny}}{10}}\)
  3. t = difference / SE

The observed difference is ~10 standard errors above zero — very far from what we’d expect if H₀ were true.

🛑 Pause — Do Activity Parts 1–3 Now

Load the data, write your hypotheses before looking at any results, and recompute the group means / SD / SE. Predict each output, type it, run it.

🧩 Chunk 2 of 4 · Framework, Hypotheses & Assumptions

We will cover: the five-step framework, stating H₀ / Hₐ, and checking normality (histogram, QQ, Shapiro-Wilk) and equal variance (Levene’s).

Tip

🖐 After this chunk: Activity Parts 4–7 (histogram, QQ, Shapiro, Levene).

Five Steps to a Hypothesis Test

We will always follow this sequence:

Step Action
1 State H₀ and Hₐ
2 Check normality — histogram/box plot
3 Check normality — QQ plot + Shapiro-Wilk
4 Check variance — Levene’s test
5 Run the test; interpret and report
Important

Check assumptions BEFORE running the test.

If your data seriously violate assumptions, the p-value from the t-test is unreliable.

With n = 10 per group, our main concern is gross non-normality. The t-test is quite robust to mild departures.

Step 1 — Hypotheses (Formal)

Null hypothesis:

\[H_0: \mu_{shady} = \mu_{sunny}\]

Shady and sunny leaves have the same mean weight in the population.

Alternate hypothesis (two-tailed):

\[H_A: \mu_{shady} \neq \mu_{sunny}\]

The mean leaf weights differ between sides — in either direction.

Significance level: α = 0.05

We will reject H₀ if p < 0.05.

Why state hypotheses first?

  • Keeps the analysis honest
  • Prevents “HARKing” (Hypothesizing After Results are Known)
  • Forces you to commit to a direction (one-tailed) or not (two-tailed) before seeing the data

Our commitment: two-tailed, α = 0.05, before running any test.

Step 2 — Check Normality: Histogram

# Histogram per group to visually check normality -------
hist_norm_plot <- tree_df %>%
  ggplot(aes(x = weight_g, fill = side)) +
  geom_histogram(binwidth = 1, color = "white", alpha = 0.8) +
  facet_wrap(~side, ncol = 2) +
  labs(
    title = "Leaf Weight Distribution by Side",
    x = "Leaf Weight (g)",
    y = "Count"
  ) +
  theme_minimal() +
  theme(legend.position = "none")
hist_norm_plot

What to look for:

  • Roughly bell-shaped — not heavily skewed
  • No obvious bimodal peaks
  • No extreme outliers

With n = 10 per group, histograms look lumpy — that is normal for small samples.

Histograms are a visual check only. We also use QQ plots and Shapiro-Wilk.

Step 3 — Check Normality: QQ Plot

# Q-Q plot — points should fall along the diagonal line --
qq_norm_plot <- tree_df %>%
  ggplot(aes(sample = weight_g, color = side)) +
  stat_qq() +
  stat_qq_line(color = "black", linewidth = 0.8) +
  facet_wrap(~side, scales = "free") +
  labs(
    title = "Normal Q-Q Plots by Side",
    x = "Theoretical Quantiles",
    y = "Sample Quantiles"
  ) +
  theme_minimal() +
  theme(legend.position = "none")
qq_norm_plot

How to read a QQ plot:

  • Points on the line → data follow a normal distribution
  • Points curving away at ends → heavy tails or skew
  • With only n = 10, some deviation is expected and acceptable

If points follow the line roughly, we are comfortable proceeding with the t-test.

Step 3 — Check Normality: Shapiro-Wilk Test

Note

🔮 Predict first: Shapiro-Wilk’s H₀ is “the data are normal.” From the histogram and QQ plot you just saw, predict for each side: will p be above or below 0.05?

# Shapiro-Wilk test for the sunny side -----------------
tree_df %>%
  filter(side == "sunny") %>%
  pull(weight_g) %>%
  shapiro.test()

    Shapiro-Wilk normality test

data:  .
W = 0.8197, p-value = 0.02513
# Shapiro-Wilk test for the shady side -----------------
tree_df %>%
  filter(side == "shady") %>%
  pull(weight_g) %>%
  shapiro.test()

    Shapiro-Wilk normality test

data:  .
W = 0.89461, p-value = 0.191

Shapiro-Wilk H₀: the data are normally distributed.

  • p > 0.05 → fail to reject normality → proceed with t-test
  • p < 0.05 → evidence of non-normality → consider a non-parametric alternative
Important

Small-sample caution:

With n = 10, this test has low power — it rarely detects non-normality even when present. The QQ plot is equally important.

Step 4 — Check Variance: Levene’s Test

# Levene's test — are the variances similar? -----------
leveneTest(weight_g ~ side, data = tree_df)
Levene's Test for Homogeneity of Variance (center = median)
      Df F value Pr(>F)
group  1     0.6 0.4486
      18               

Levene’s test H₀: the two groups have equal variance.

  • p > 0.05 → variances are not significantly different
  • p < 0.05 → variances are significantly different

Our response regardless of the result:

We will use Welch’s t-test (var.equal = FALSE).

Welch’s is valid whether or not the variances are equal. We run Levene’s to understand our data, not to choose our test.

Live Demo — Watch It Break (on purpose)

Suppose I forgot library(car) and run Levene’s:

leveneTest(weight_g ~ side, data = tree_df)

R stops with:

Error in leveneTest(weight_g ~ side, data = tree_df) :
  could not find function "leveneTest"

The fix — load the package that contains the function:

library(car)
leveneTest(weight_g ~ side, data = tree_df)
Tip

✅ Why show a broken run?

“could not find function” almost always means a missing library() — not a typo. Watching me diagnose it here means you’ll fix it in seconds when it happens to you.

🛑 Pause — Do Activity Parts 4–7 Now

Check every assumption before you test: histogram, QQ plot, Shapiro-Wilk, and Levene’s. Predict each p-value before you run it.

🧩 Chunk 3 of 4 · Run & Interpret the Test

We will cover: running Welch’s t.test(), reading every line of the output, and making a formal reject / fail-to-reject decision.

Tip

🖐 After this chunk: Activity Parts 8–11 (run, decode, calculate t by hand, decide).

Step 5 — Run the Welch’s t-Test

Note

🔮 Predict first: Write down a guess for the p-value before you run t.test() — above or below 0.05? Then see how close you were.

# Welch's two-sample t-test — unequal variance ---------
leaf_ttest_model <- t.test(
  weight_g ~ side, # formula: response ~ grouping variable
  data = tree_df,
  var.equal = FALSE, # Welch's — no pooled variance
  alternative = "two.sided" # two-tailed test
)

leaf_ttest_model

    Welch Two Sample t-test

data:  weight_g by side
t = 9.7333, df = 16.834, p-value = 2.514e-08
alternative hypothesis: true difference in means between group shady and group sunny is not equal to 0
95 percent confidence interval:
 3.132297 4.867703
sample estimates:
mean in group shady mean in group sunny 
                7.8                 3.8 

Key arguments:

argument what it does
weight_g ~ side compare weight_g between levels of side
var.equal = FALSE use Welch’s correction
alternative = "two.sided" two-tailed test

The result is a model object. We decode each line of the output on the next slide.

Reading the t.test() Output

# Extract individual values from the model object ------
cat("t-statistic  :", round(leaf_ttest_model$statistic, 3), "\n")
t-statistic  : 9.733 
cat("df (Welch)   :", round(leaf_ttest_model$parameter, 2), "\n")
df (Welch)   : 16.83 
cat("p-value      :", signif(leaf_ttest_model$p.value, 3), "\n")
p-value      : 2.51e-08 
cat(
  "95% CI       :",
  round(leaf_ttest_model$conf.int[1], 2),
  "to",
  round(leaf_ttest_model$conf.int[2], 2),
  "g\n"
)
95% CI       : 3.13 to 4.87 g
cat("Mean shady   :", round(leaf_ttest_model$estimate[1], 2), "g\n")
Mean shady   : 7.8 g
cat("Mean sunny   :", round(leaf_ttest_model$estimate[2], 2), "g\n")
Mean sunny   : 3.8 g

Line by line:

output meaning
t our calculated t-statistic
df Welch-Satterthwaite df (non-integer is normal!)
p-value probability of this result if H₀ is true
95% CI plausible range for the true difference in means
mean of shady/sunny the two group means

The 95% CI is the range we could expect to see if we repeated this study many times.

Step 5 — Make a Decision

# State the decision formally --------------------------
# if/else: R checks the condition and runs one of two blocks
# We will learn these formally in the functions lecture
p_val <- leaf_ttest_model$p.value
alpha <- 0.05

if (p_val < alpha) {
  cat("p =", signif(p_val, 3), "< α =", alpha, "\n")
  cat("Decision: REJECT H₀\n")
  cat("Conclusion: leaf weight differs between sides.\n")
} else {
  cat("p =", signif(p_val, 3), ">= α =", alpha, "\n")
  cat("Decision: FAIL TO REJECT H₀\n")
  cat("Conclusion: insufficient evidence of a difference.\n")
}
p = 2.51e-08 < α = 0.05 
Decision: REJECT H₀
Conclusion: leaf weight differs between sides.
Important

What a small p-value means:

If H₀ were true (no real difference), we would get a t-statistic this large or larger less than p × 100% of the time just by chance.

We set α = 0.05 before the test. If p < α, we reject H₀.

What a small p-value does NOT mean:

  • It is not the probability that H₀ is true
  • Statistical significance ≠ biological importance
  • Report the effect size (difference in means) alongside p

🛑 Pause — Do Activity Parts 8–11 Now

Run the test, decode t / df / p / CI, reproduce t by hand, and state your decision in formal language. Predict the p-value first.

🧩 Chunk 4 of 4 · Visualize & Report

We will cover: embedding the test result in a boxplot and a mean ± SE plot, writing a scientific results sentence, and saving publication-quality figures.

Tip

🖐 After this chunk: Activity Parts 12–13 (plots and the results sentence).

Visualize the Result — Boxplot

# Final boxplot with all points -----------------------
leaf_box_04_plot <- tree_df %>%
  ggplot(aes(x = side, y = weight_g, fill = side)) +
  geom_boxplot(alpha = 0.5, outlier.shape = NA) +
  geom_point(
    position = position_jitter(width = 0.15, seed = 42),
    alpha = 0.6,
    size = 2.5
  ) +
  labs(
    title = "Leaf Weight by Tree Side",
    subtitle = paste0(
      "Welch's t-test: t = ",
      round(leaf_ttest_model$statistic, 2),
      ", p = ",
      signif(leaf_ttest_model$p.value, 2)
    ),
    x = "Side of Tree",
    y = "Leaf Weight (g)"
  ) +
  theme_minimal() +
  theme(legend.position = "none")
leaf_box_04_plot

Including the test result in the figure:

  • subtitle pulls t and p directly from the model object
  • R fills in the actual numbers automatically
  • If you re-run with new data, the subtitle updates
Tip

Embedding test statistics in the plot subtitle saves you from copying numbers by hand — and avoids typos in your report.

Visualize the Result — Mean ± SE

# Mean ± SE plot with individual points ---------------
leaf_mean_se_04_plot <- tree_df %>%
  ggplot(aes(x = side, y = weight_g, color = side)) +
  geom_point(
    position = position_jitter(width = 0.15, seed = 42),
    alpha = 0.35,
    size = 2
  ) +
  stat_summary(fun = mean, geom = "point", size = 4) +
  stat_summary(
    fun.data = mean_se,
    geom = "errorbar",
    width = 0.15,
    linewidth = 0.9
  ) +
  labs(
    title = "Mean ± SE Leaf Weight by Tree Side",
    subtitle = paste0(
      "Welch's t: t(",
      round(leaf_ttest_model$parameter, 1),
      ") = ",
      round(leaf_ttest_model$statistic, 2),
      ", p < 0.001"
    ),
    x = "Side of Tree",
    y = "Leaf Weight (g)"
  ) +
  theme_minimal() +
  theme(legend.position = "none")
leaf_mean_se_04_plot

What the error bars tell you:

  • Bars show ± 1 SE — precision of the mean estimate
  • Non-overlapping bars suggest a significant difference

Rule of thumb:

If error bars overlap by less than half the bar length, the difference is likely significant at α = 0.05. The p-value from the t-test is the definitive answer.

How to Report Results in Science

# Compute reporting values from the model object -------
t_val <- round(leaf_ttest_model$statistic, 2)
df_val <- round(leaf_ttest_model$parameter, 1)
p_val2 <- signif(leaf_ttest_model$p.value, 2)

cat("--- Results sentence ---\n")
--- Results sentence ---
cat("Shady-side leaves were significantly heavier than\n")
Shady-side leaves were significantly heavier than
cat(
  "sunny-side leaves (Welch's t-test: t(",
  df_val,
  ") =",
  t_val,
  ", p =",
  p_val2,
  ").\n"
)
sunny-side leaves (Welch's t-test: t( 16.8 ) = 9.73 , p = 2.5e-08 ).
cat(
  "Mean ± SE: shady =",
  round(leaf_ttest_model$estimate[1], 2),
  "±",
  round(recap_df$se_wt[recap_df$side == "shady"], 2),
  "g,\n"
)
Mean ± SE: shady = 7.8 ± 0.33 g,
cat(
  "           sunny =",
  round(leaf_ttest_model$estimate[2], 2),
  "±",
  round(recap_df$se_wt[recap_df$side == "sunny"], 2),
  "g.\n"
)
           sunny = 3.8 ± 0.25 g.

Standard format:

“[Finding] (Welch’s t-test: t(df) = X.XX, p = X.XXX; mean ± SE: group1 = X.X ± X.X units, group2 = X.X ± X.X units).”

Always include:

  • Test type (Welch’s t-test)
  • t-statistic and df
  • Exact p-value (or < 0.001)
  • Group means ± SE
  • Units

Never report “p = 0.000” — write “p < 0.001” instead.

Save Your Plots

# Save both plots to the figures folder ----------------
ggsave(
  "figures/leaf_weight_boxplot_ttest.png",
  plot = leaf_box_04_plot,
  width = 5,
  height = 5,
  units = "in",
  dpi = 300
)

ggsave(
  "figures/leaf_weight_mean_se_ttest.png",
  plot = leaf_mean_se_04_plot,
  width = 5,
  height = 5,
  units = "in",
  dpi = 300
)

Both plots are now in your figures/ folder with the test statistics embedded in the subtitle — ready for a lab report or paper.

Tip

dpi = 300 is publication quality. Use it for any figure that will appear in a paper, poster, or formal report.

🛑 Pause — Do Activity Parts 12–13 Now

Build both figures with the statistics embedded in the subtitle, then write your results sentence. Predict what the subtitle will say before you run it.

What We Learned Today

Concepts:

  • Why we need a hypothesis test — means alone are not enough
  • Welch’s t-test — use var.equal = FALSE as the safe default
  • Two-tailed test — detects differences in either direction
  • Five-step framework — hypotheses → normality → variance → test → report
  • Reading output — t, df, p, CI, group means

R skills:

  • library(car) + leveneTest(y ~ group, data) — variance check
  • filter() %>% pull() %>% shapiro.test() — normality per group
  • t.test(y ~ group, data, var.equal = FALSE, alternative = "two.sided")
  • model$statistic, model$p.value, model$conf.int — extract results
  • paste0(...) in subtitle — embed test statistics in figures

References:

Up next — Lecture 05:

  • Linear regression — predicting leaf area from paper tracing mass
  • Building a calibration curve with lm()
  • R², residuals, and assumption checking