Lecture 06 — Summary Statistics

Mean, median, spread, and describing groups the tidy way

tidyverse

descriptive-stats

Mean, median, variance, SD, and SE by hand; counting correctly with NAs; group_by() + summarize(), across(), and skimr for describing our leaf data.

Author

Bill Perry

Published

September 10, 2026

Where we left off (Lecture 05 — Wrangling)

filter(), select(), mutate(), arrange() — the core wrangling verbs
Chaining verbs into one pipeline with %>%
Our leaf data: leaf_df — shade, mass_g, petiole_mm, thickness_mm, paper_mass_g

Note

✅ Key idea from Lecture 05

You can already reshape the data into exactly the rows and columns you need. Today we learn to describe what’s in it — one number at a time, then all at once.

Our Scientific Question

Biological prediction: Leaves on the shady side of a tree will be larger and heavier to capture more of the limited light.

	Hypothesis
H₀ — null	No difference in leaf size between sunny and shady sides
Hₐ — alternate	Leaf size differs between sunny and shady sides

Plan for this unit:

Today — describe the data
Lecture 09 — formally test the hypotheses with a t-test

Statistics helps us decide which hypothesis is better supported by the data.

References:

📖 Whitlock & Schluter, Ch. 3 — Describing Data
📖 Whitlock & Schluter, Ch. 4 — Estimating with Uncertainty
📖 R4DS Ch 13 — Numbers
📖 R4DS Ch 18 — Missing Values

Both W&S PDFs are in the course readings/ folder.

How to Use These Slides — Predict · Type · Run

This lecture runs in two chunks. After each chunk you switch to the activity and type the code yourself into your R script.

For every code block, do three things:

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

Note

✅ Why bother?

A single mean or SD looks obvious once it’s on the slide — predicting it first is what tells you whether you actually understand where the number comes from.
Typing sum(!is.na(x)) yourself, instead of pasting it, is what makes you notice why it’s there instead of length(x).

Load Libraries and Data

# Load all packages at the top of every script ----------
library(readxl) # reading Excel files
library(tidyverse) # data wrangling + ggplot2
library(janitor) # clean_names()
library(skimr) # fast descriptive summaries

# Read the leaf data from the data folder ---------------
leaf_df <- read_excel("data/2026_09_03_data_sci_leaf_area.xlsx") %>%
  clean_names()

Note

Install once — load every session

Run in the Console one time only: install.packages("skimr")
Then every session: library(skimr) activates it.

🧩 Chunk 1 of 2 · Centre, Spread, and Counting Correctly

We will cover: the mean, median, variance, SD, SE, and the length() trap.

🛑 After this chunk: Activity Parts 1–3 (load, the NA problem, and base-R stats).

What Is the Mean?

\[\bar{x} = \frac{\sum_{i=1}^{n} x_i}{n}\]

Add up all values, divide by n
The balance point of the distribution
Sensitive to outliers — one very large value pulls the mean up

Example:

Values: 5, 6, 7, 8, 50 → mean = 15.2

Is 15.2 a fair summary? Not really — the outlier dominates!

Why it matters:

The most commonly reported statistic
Used in every t-test, ANOVA, and regression
Always pair it with a measure of spread (SD or SE)

📖 Whitlock & Schluter, Ch. 3 — Describing Data; R4DS §13 — Numbers

Calculating the Mean — filter() and pull()

# Use filter() + pull() to extract one group's values --
sunny_mass <- leaf_df %>%
  filter(shade == "sunny") %>%
  pull(mass_g)

sunny_mass # look at the raw values

 [1] 0.4000 0.4800 0.3400 0.6500 0.2700 0.4300 0.7366 0.4263 0.5682 0.7600
[11] 0.5777 0.1933 0.7912 0.4990 0.6882 0.4162 0.6185 1.0842 0.6649 0.3571
[21] 0.3484 0.5107 0.2303 0.3041 0.2817

# Calculate the mean -----------------------------------
mean_sunny <- mean(sunny_mass, na.rm = TRUE)
cat("Mean sunny mass:", round(mean_sunny, 3), "g\n")

Mean sunny mass: 0.505 g

filter(shade == "sunny") — keeps only sunny rows
pull(mass_g) — extracts that column as a vector
na.rm = TRUE — removes any NA before computing

Key pattern:

filter() → picks rows

pull() → extracts a column as a vector

📖 R4DS §3 — Data Transformation

What Is the Median?

The middle value when data are sorted

Sorted data	Median
5, 6, 7, 8, 50	7
5, 6, 7, 8, 50, 51	(7 + 8) / 2 = 7.5

Not pulled by outliers → robust
If mean ≫ median → right-skewed
If mean ≈ median → roughly symmetric

# Median for each side ---------------------------------
shady_mass <- leaf_df %>%
  filter(shade == "shady") %>%
  pull(mass_g)

med_sunny <- median(sunny_mass, na.rm = TRUE)
med_shady <- median(shady_mass, na.rm = TRUE)

cat("Median sunny:", round(med_sunny, 3), "\n")

Median sunny: 0.48

cat("Median shady:", round(med_shady, 3), "\n")

Median shady: 0.521

What Is Variance?

\[s^2 = \frac{\sum_{i=1}^{n}(x_i - \bar{x})^2}{n - 1}\]

For each value, compute distance from mean: \((x_i - \bar{x})\)
Square those distances → all positive
Sum them up, divide by \(n - 1\)

Dividing by \(n - 1\) gives an unbiased estimate of population variance.

Larger variance = more spread = more uncertainty

# Variance for each side ------------------------------
var_sunny <- var(sunny_mass, na.rm = TRUE)
var_shady <- var(shady_mass, na.rm = TRUE)

cat("Variance sunny:", round(var_sunny, 4), "\n")

Variance sunny: 0.0442

cat("Variance shady:", round(var_shady, 4), "\n")

Variance shady: 0.0219

Note

Variance is in squared units (g²) — hard to interpret. That is why we use SD.

📖 R4DS §13.3

What Is Standard Deviation?

\[s = \sqrt{s^2} = \sqrt{\frac{\sum(x_i - \bar{x})^2}{n-1}}\]

The square root of variance → back in original units (g) ✓
For a normal distribution:
- ~68% of values fall within mean ± 1 SD
- ~95% fall within mean ± 2 SD
Reports how spread out individual data points are

# Standard deviation for each side -------------------
sd_sunny <- sd(sunny_mass, na.rm = TRUE)
sd_shady <- sd(shady_mass, na.rm = TRUE)

cat("SD sunny:", round(sd_sunny, 3), "g\n")

SD sunny: 0.21 g

cat("SD shady:", round(sd_shady, 3), "g\n")

SD shady: 0.148 g

What Is Standard Error?

\[SE = \frac{s}{\sqrt{n}}\]

Measures how precisely we know the mean
Not the spread of raw data — the uncertainty of the mean itself
Gets smaller as n increases → more data = more precise estimate

	Describes
SD	spread of individual data points
SE	precision of the sample mean

We report SE when making claims about group means.

📖 Whitlock & Schluter, Ch. 4 — Estimating with Uncertainty

# SE = SD / sqrt(n) — need correct n -----------------
n_sunny <- sum(!is.na(sunny_mass))
n_shady <- sum(!is.na(shady_mass))

se_sunny <- sd_sunny / sqrt(n_sunny)
se_shady <- sd_shady / sqrt(n_shady)

cat("n sunny  =", n_sunny, "\n")

n sunny  = 25

cat("SE sunny =", round(se_sunny, 3), "g\n\n")

SE sunny = 0.042 g

cat("n shady  =", n_shady, "\n")

n shady  = 28

cat("SE shady =", round(se_shady, 3), "g\n")

SE shady = 0.028 g

The Problem with `length()`

# length() counts ALL positions — including NAs -------
x <- c(4, 3, NA, 7, NA) # a vector with missing values

length(x) # returns 5 — but 2 are NA!

[1] 5

sum(is.na(x)) # how many NAs are there?

[1] 2

# Using length() in SE gives the WRONG answer ---------
wrong_n <- length(x)
wrong_se <- sd(x, na.rm = TRUE) / sqrt(wrong_n)

cat("Wrong n  =", wrong_n, "\n")

Wrong n  = 5

cat("Wrong SE =", round(wrong_se, 3), "\n")

Wrong SE = 0.931

length(x) counts all positions in the vector — including NA
This is a silent error — R will not warn you!
Our leaf data has real NAs too: not every leaf got a paper_mass_g measurement
Using n = 5 in SE = SD / √n when n should be 3 gives the wrong answer

📖 R4DS §18 — Missing Values

Fix — `sum(!is.na())`

# Break it down step by step -------------------------
x <- c(4, 3, NA, 7, NA)

is.na(x) # TRUE where NA, FALSE where real

[1] FALSE FALSE  TRUE FALSE  TRUE

!is.na(x) # FLIP: TRUE where real, FALSE where NA

[1]  TRUE  TRUE FALSE  TRUE FALSE

sum(!is.na(x)) # count the TRUEs = number of real values

[1] 3

# Apply to our real leaf data ------------------------
n_sun_check <- sum(!is.na(sunny_mass))
cat("Non-missing sunny masses =", n_sun_check, "\n")

Non-missing sunny masses = 25

How it works step by step:

is.na(x) → TRUE for each NA
!is.na(x) → flips it: TRUE for real values
sum(...) → counts TRUEs

Tip

Best practice: always use sum(!is.na(x)) when you need n — even when you think there are no NAs.

All Stats — Base R, One Group at a Time

# All stats for the sunny side -----------------------
n_s <- sum(!is.na(sunny_mass))
mean_s <- mean(sunny_mass, na.rm = TRUE)
med_s <- median(sunny_mass, na.rm = TRUE)
sd_s <- sd(sunny_mass, na.rm = TRUE)
se_s <- sd_s / sqrt(n_s)

cat("--- Sunny Leaves ---\n")

--- Sunny Leaves ---

cat("n      =", n_s, "\n")

n      = 25

cat("Mean   =", round(mean_s, 3), "\n")

Mean   = 0.505

cat("Median =", round(med_s, 3), "\n")

Median = 0.48

cat("SD     =", round(sd_s, 3), "\n")

SD     = 0.21

cat("SE     =", round(se_s, 3), "\n")

SE     = 0.042

Every function takes na.rm = TRUE
cat() prints clean labelled output
This is clear but repetitive for many groups → that is why we use group_by() next

🛑 Pause — Do Activity Parts 1–3 Now

🛑 Do Activity Parts 1–3 now

Load the data, work the NA problem, and compute mean / median / spread in base R. Predict each output, type it, then run it.

🧩 Chunk 2 of 2 · Tidy Stats at Scale

We will cover: group_by() + summarize(), summarizing many columns at once with across(), and skimr.

🛑 After this chunk: Activity Parts 4–6.

Stats the Tidy Way — `group_by()` + `summarize()`

Note

🔮 Predict first: How many rows will stats_df have? (Hint: how many values does shade take?) Predict the number before you run it.

# Calculate all stats for both sides at once ----------
stats_df <- leaf_df %>%
  group_by(shade) %>%
  summarize(
    n = sum(!is.na(mass_g)),
    mean_mass = round(mean(mass_g, na.rm = TRUE), 3),
    med_mass = round(median(mass_g, na.rm = TRUE), 3),
    sd_mass = round(sd(mass_g, na.rm = TRUE), 3),
    se_mass = round(sd_mass / sqrt(n), 3)
  )

stats_df

# A tibble: 2 × 6
  shade     n mean_mass med_mass sd_mass se_mass
  <chr> <int>     <dbl>    <dbl>   <dbl>   <dbl>
1 shady    28     0.523    0.521   0.148   0.028
2 sunny    25     0.505    0.48    0.21    0.042

group_by(shade) — splits the data by the shade column
summarize() — collapses each group to one row of stats
You get both groups in one clean table
You can reference a column just created (sd_mass → se_mass)

📖 R4DS §3.5 — summarize()

Multiple Variables at Once

# One group_by() pipe summarizes all measurements -----
size_stats_df <- leaf_df %>%
  group_by(shade) %>%
  summarize(
    n = sum(!is.na(mass_g)),
    mean_mass = round(mean(mass_g, na.rm = TRUE), 3),
    mean_thick = round(mean(thickness_mm, na.rm = TRUE), 3),
    mean_pet = round(mean(petiole_mm, na.rm = TRUE), 1)
  )

size_stats_df

# A tibble: 2 × 5
  shade     n mean_mass mean_thick mean_pet
  <chr> <int>     <dbl>      <dbl>    <dbl>
1 shady    28     0.523      0.356     56.7
2 sunny    25     0.505      0.29      54.4

With one group_by() pipe you summarize all three measurements across both sides.

Do the numbers support our hypothesis?

Shady heavier? → check mean_mass
Shady leaves thicker? → check mean_thick
Shady petioles longer? → check mean_pet

The Fast Way — `across()`

# Apply the SAME function to every numeric column -----
mean_all_df <- leaf_df %>%
  group_by(shade) %>%
  summarize(
    across(
      c(mass_g, petiole_mm, thickness_mm),
      ~ round(mean(.x, na.rm = TRUE), 3)
    )
  )

mean_all_df

# A tibble: 2 × 4
  shade mass_g petiole_mm thickness_mm
  <chr>  <dbl>      <dbl>        <dbl>
1 shady  0.523       56.7        0.356
2 sunny  0.505       54.4        0.29

across(columns, function) applies one function to many columns at once
.x is a stand-in for “whichever column we’re on right now”
Same idea as mean_mass/mean_thick/mean_pet above — just without typing each one out

Tip

Once you have more than 3–4 measurement columns, across() saves real typing — and real typos.

Fast Summaries with `skimr`

# skim() grouped by shade — full column summaries -----
leaf_df %>%
  group_by(shade) %>%
  skim()

Data summary
Name	Piped data
Number of rows	53
Number of columns	8
_______________________
Column type frequency:
character	3
numeric	4
________________________
Group variables	shade

Variable type: character

skim_variable	shade	n_missing	complete_rate	min	max	n_unique
twig_id	shady	22	0.21	6	6	3
twig_id	sunny	19	0.24	6	6	3
leaf_id	shady	16	0.43	1	3	12
leaf_id	sunny	13	0.48	1	3	12
teams	shady	0	1.00	5	18	5
teams	sunny	0	1.00	5	18	5

Variable type: numeric

skim_variable	shade	n_missing	complete_rate	mean	sd	p0	p25	p50	p75	p100	hist
mass_g	shady	0	1.00	0.52	0.15	0.34	0.42	0.52	0.59	0.92	▇▅▃▁▁
mass_g	sunny	0	1.00	0.51	0.21	0.19	0.35	0.48	0.65	1.08	▇▇▆▃▁
petiole_mm	shady	0	1.00	56.66	17.53	0.38	50.28	58.69	66.85	91.00	▁▁▅▇▂
petiole_mm	sunny	0	1.00	54.39	16.92	29.00	44.85	49.00	64.00	85.00	▆▇▂▃▅
thickness_mm	shady	1	0.96	0.36	0.16	0.12	0.23	0.36	0.48	0.67	▇▆▃▇▂
thickness_mm	sunny	0	1.00	0.29	0.19	0.10	0.14	0.23	0.40	0.69	▇▂▃▁▂
paper_mass_g	shady	7	0.75	0.28	0.04	0.20	0.25	0.27	0.30	0.36	▁▇▃▂▃
paper_mass_g	sunny	6	0.76	0.27	0.09	0.11	0.21	0.25	0.34	0.43	▂▇▃▅▅

Per column and per group: n_missing, mean, sd, percentiles (p0 to p100), and a mini histogram.

Tip

skim() is your fastest first look at any dataset. Pair it with group_by() to compare groups instantly.

→ ACTIVITY 6 starts now

🛑 Go to Activity 6 — Summary Statistics

Close the slides. Finish Activity Parts 4–6: group_by() + summarize(), across(), and skimr.

What We Learned Today

Statistics concepts:
- Mean — balance point; sensitive to outliers
- Median — middle value; robust
- Variance / SD — spread of individual data
- SE — precision of the sample mean
R skills:
- sum(!is.na()) — the right way to count n
- group_by() + summarize() — all stats for all groups at once
- across() — one function, many columns
- skim() — instant full dataset overview

References:

📖 Whitlock & Schluter, Ch. 3 — Describing Data
📖 Whitlock & Schluter, Ch. 4 — Estimating with Uncertainty
📖 R4DS §13 — Numbers
📖 R4DS §18 — Missing Values

Both W&S PDFs are in the course readings/ folder.

Up next — Lecture 07, Quarto:

Writing reproducible reports
Mixing prose, code, and results in one document

Other Formats

Where we left off (Lecture 05 — Wrangling)

Our Scientific Question

How to Use These Slides — Predict · Type · Run

Load Libraries and Data

🧩 Chunk 1 of 2 · Centre, Spread, and Counting Correctly

What Is the Mean?

Calculating the Mean — filter() and pull()

What Is the Median?

What Is Variance?

What Is Standard Deviation?

What Is Standard Error?

The Problem with length()

Fix — sum(!is.na())

All Stats — Base R, One Group at a Time

🛑 Pause — Do Activity Parts 1–3 Now

🧩 Chunk 2 of 2 · Tidy Stats at Scale

Stats the Tidy Way — group_by() + summarize()

Multiple Variables at Once

The Fast Way — across()

Fast Summaries with skimr

→ ACTIVITY 6 starts now

What We Learned Today

The Problem with `length()`

Fix — `sum(!is.na())`

Stats the Tidy Way — `group_by()` + `summarize()`

The Fast Way — `across()`

Fast Summaries with `skimr`