Activity 06 - Summary statistics

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

tidyverse

descriptive-stats

Hands-on companion to the summary statistics lecture. Handle NAs correctly with sum(!is.na()), source a reusable summary_stats() helper, use it by group with group_by(), explore with skimr, then summarize many columns at once with across().

Author

Bill Perry

Published

October 1, 2026

Describing Our Leaf Data

Recap from Activity 05

Used filter(), select(), mutate(), and arrange() to wrangle the leaf data
Chained all four verbs into a single pipeline
Used the pipe %>% to chain steps

Today’s Objectives

Find missing values with filter(is.na()) and filter(!is.na())
Understand why NA values require special handling with sum(!is.na())
Write the whole summary by hand — group_by() + summarize() with mean, variance, sd, n, and se
Source a reusable summary_stats() function that packages up that same code
Explore data quickly with skimr
Summarize many columns at once with across()

How this activity works

You are building one R script this whole class, and you turn it in. Create it now: scripts/06_summary_stats.R.

Every line you run goes in that script — not in the Console, not typed into this page. Type it, don’t paste it.

Start every code chunk in your script with a short # comment that says what it does. The comments are part of the grade.

The top of your script, in this order: a title comment, then your library() calls, then the line that loads the data into leaf_df.

Code marked ▶ Run this is typed into your script exactly as shown. Code marked ✏️ Your turn is a change you make in that same script and run.

🔮 Predict before you run

Before you run any ▶ Run this block, predict the answer — a row count, a mean, which group comes out larger. Predicting is what turns “I saw it on a slide” into “I can write it.”

Part 1 · Load libraries and data

📂 Get the data

Same leaf file you have used since Activity 02 — it should already be in your project’s data/ folder. If not: 2026_09_03_data_sci_leaf_area.xlsx → put it in data/.

📂 Get the skeleton script

06_summary_stats_skeleton.R → save it as scripts/06_summary_stats.R. It already has the library calls, the line that loads leaf_df, and the two source() lines below — so you don’t retype the plumbing, just the analysis.

Also grab r_themes_for_3_sizes.R and summary_stats_function.R if they aren’t already in your project’s themes/ folder (you made that folder back in Activity 04).

Open the skeleton and read the top of it — this part is done for you:

# ---- Activity 06: Summary statistics ---------------------
# your name, today's date

# ---- Libraries -------------------------------------------
library(readxl)      # read Excel files
library(tidyverse)   # dplyr + ggplot2
library(janitor)     # clean_names()
library(skimr)       # fast descriptive summaries

# ---- Load data ------------------------------------------
leaf_df <-
  read_excel("data/2026_09_03_data_sci_leaf_area.xlsx") %>%
  clean_names()

# ---- Load our theme file and summary_stats() helper ------
source("themes/r_themes_for_3_sizes.R")
source("themes/summary_stats_function.R")

✏️ Your turn: Run this whole top section once, top to bottom, so leaf_df, theme_regular(), and summary_stats() are all ready before Part 2.

Part 2 · The NA problem — counting observations correctly

What is an NA?

NA stands for “Not Available” — it marks a missing value. Missing data are common in real ecological studies. Our leaf data has some: not every leaf got a paper_mass_g measurement.

Step 1 — look at the missing values directly

▶ Run this:

# which rows are MISSING a paper_mass_g value?
leaf_df %>%
  filter(is.na(paper_mass_g))

How many rows came back?

▶ Run this:

# flip it with ! — keep only the rows that HAVE a value
leaf_df %>%
  filter(!is.na(paper_mass_g))

How many rows came back?
Do the two row counts add up to nrow(leaf_df)?

💡 Key idea: is.na() asks “is this missing?” and ! flips the answer. Look at the missing rows before you compute a single statistic — you cannot trust a mean until you know what is missing from it.

Step 2 — the problem with `length()`

▶ Run this:

# a small vector with two missing values, for demonstration
x <- c(0.4, 0.3, NA, 0.7, NA)

# length() counts ALL positions — including the NAs
length(x)

# sum(!is.na()) counts only the REAL values
sum(!is.na(x))

length(x) returned:
sum(!is.na(x)) returned:

✏️ Your turn — in your script: Now do the same on the real data. Compare nrow(leaf_df) to sum(!is.na(leaf_df$paper_mass_g)).

nrow(leaf_df):
sum(!is.na(leaf_df$paper_mass_g)):
How many paper_mass_g values are missing?

⚠️ Watch out! Any stats function on a vector containing NA returns NA — not an error. R will not tell you something went wrong. Always include na.rm = TRUE.

Part 3 · Do it by hand — `summarize()`

Before we use any helper function, write every statistic out yourself. This is the code the function will replace later — you need to be able to read it.

One group at a time

▶ Run this:

# every statistic, written out, for the sunny leaves
leaf_df %>%
  filter(shade == "sunny") %>%
  summarize(
    mean     = mean(mass_g, na.rm = TRUE),
    variance = var(mass_g, na.rm = TRUE),
    sd       = sd(mass_g, na.rm = TRUE),
    n        = sum(!is.na(mass_g)),
    se       = sd / sqrt(n)
  )

Every function gets na.rm = TRUE — otherwise one NA makes the whole answer NA
n = sum(!is.na(mass_g)) — the correct count, not n()
se = sd / sqrt(n) — you can reuse sd and n later in the same summarize(), because by then they are already columns

sunny: mean, variance, sd, n, se

✏️ Your turn: Do the shady side the same way — change "sunny" to "shady" and run it again.

shady: mean, variance, sd, n, se
Which side has the larger mean? The larger sd?

Both groups at once — `group_by()`

🔮 Predict first: How many rows will this return? (How many values does shade take?)

▶ Run this:

# group_by() splits the data; summarize() runs once per group
leaf_df %>%
  group_by(shade) %>%
  summarize(
    mean     = mean(mass_g, na.rm = TRUE),
    variance = var(mass_g, na.rm = TRUE),
    sd       = sd(mass_g, na.rm = TRUE),
    n        = sum(!is.na(mass_g)),
    se       = sd / sqrt(n)
  )

💡 Key idea: the summarize() did not change at all — the only new line is group_by(shade). Running the filter version twice, once per side, is exactly the repetition group_by() exists to remove.

Do these two rows match what you got with filter() above?  Y / N

Part 4 · Package it — `summary_stats()`

You just typed that summarize() twice. Rather than retype it for every variable for the rest of the term, we wrap it in a function once. It is the same code — plus a 95% confidence interval for the mean.

summary_stats <- function(data, variable) {
  data %>%
    summarize(
      n        = sum(!is.na({{ variable }})),
      mean     = mean({{ variable }}, na.rm = TRUE),
      variance = var({{ variable }}, na.rm = TRUE),
      sd       = sd({{ variable }}, na.rm = TRUE),
      se       = sd / sqrt(n),
      ci_lower = mean - qt(0.975, df = n - 1) * se,
      ci_upper = mean + qt(0.975, df = n - 1) * se,
      .groups  = "drop"
    )
}

n = sum(!is.na({{ variable }})) — the correct count, not n()
{ variable } lets you call summary_stats(leaf_df, mass_g) with a bare column name
se and the confidence interval both reuse mean and sd computed earlier in the same summarize()

📂 Get the file

summary_stats_function.R — already sourced for you at the top of the skeleton script from Part 1, in the same themes/ folder you made in Activity 04.

✏️ Your turn — in your script: Filter leaf_df to shade == "sunny" and pipe the result straight into summary_stats(mass_g).

Same n, mean, variance, sd, se as your
hand-written version in Part 3?  Y / N

▶ Run this:

# group first, then call the SAME function — no rewriting it
stats_df <- leaf_df %>%
  group_by(shade) %>%
  summary_stats(mass_g)

stats_df

💡 Key idea: summary_stats() did not have to change to work on groups. Grouped data goes in, one row per group comes out — because summarize() inside the function behaves that way.

✏️ Your turn: Run summary_stats() on a different column — petiole_mm or thickness_mm — grouped by shade. That one-word change is the whole point of writing the function.

Column you chose, and its two means:

Part 5 · Several variables at once

summary_stats() handles one variable at a time. For several measurements at once, write summarize() directly:

▶ Run this:

# mean of all three leaf measurements, by shade
size_stats_df <- leaf_df %>%
  group_by(shade) %>%
  summarize(
    n          = sum(!is.na(mass_g)),
    mean_mass  = mean(mass_g, na.rm = TRUE),
    mean_pet   = mean(petiole_mm, na.rm = TRUE),
    mean_thick = mean(thickness_mm, na.rm = TRUE)
  )

size_stats_df

✏️ Your turn: Looking at size_stats_df, fill in the table:

Measurement    | Sunny mean | Shady mean | Shady larger? (Y/N)
---------------|------------|------------|---------------------
mass_g         |            |            |
petiole_mm     |            |            |
thickness_mm   |            |            |

Part 6 · Fast overview with `skimr`

▶ Run this:

# skim() grouped by shade — n_missing,
# mean, SD, and percentiles at once
leaf_df %>%
  group_by(shade) %>%
  skim()

✏️ Your turn: Look at the n_missing row for paper_mass_g, and the mean vs p50 (median) for mass_g.

n_missing paper_mass_g (sunny / shady):
Does it match the row count from your filter(is.na()) in Part 2?
Are mean and median of mass_g close for each
group? What does that suggest about skew?

💡 Key idea: one line of skim() gets you n_missing, mean, SD, median, and percentiles for every column, grouped by shade — a fast first look, not a replacement for the summary you write yourself.

Part 7 · The fast way — `across()`

▶ Run this:

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

mean_all_df

✏️ Your turn — in your script: Change mean inside across() to sd and run it.

Same numbers as size_stats_df?  Y / N
What changed:

Part 8 · Review and checkpoint

At this point you should be able to:

Use filter(is.na()) and filter(!is.na()) to see the missing data before summarizing
Explain why length() / n() is unreliable with NAs and use sum(!is.na()) instead
Write a group_by() + summarize() by hand with mean, variance, sd, n, and se
Source summary_stats() from themes/ and use it on a filtered group and on grouped data
Use skim() for a fast full-dataset overview
Use across() to apply one function to many columns at once

✏️ Your turn — before you move on: Run your entire script top to bottom with Ctrl/Cmd + Shift + Enter (Source). Does it complete without errors?

Ran cleanly?  Y / N
If not, what error appeared:

End of the Summary Statistics activity. Next: Quarto — writing reproducible reports that mix prose, code, and results in one document.

Other Formats

Describing Our Leaf Data

Recap from Activity 05

Today’s Objectives

How this activity works

🔮 Predict before you run

Part 1 · Load libraries and data

Part 2 · The NA problem — counting observations correctly

What is an NA?

Step 1 — look at the missing values directly

Step 2 — the problem with length()

Part 3 · Do it by hand — summarize()

One group at a time

Both groups at once — group_by()

Part 4 · Package it — summary_stats()

Part 5 · Several variables at once

Part 6 · Fast overview with skimr

Part 7 · The fast way — across()

Part 8 · Review and checkpoint

Step 2 — the problem with `length()`

Part 3 · Do it by hand — `summarize()`

Both groups at once — `group_by()`

Part 4 · Package it — `summary_stats()`

Part 6 · Fast overview with `skimr`

Part 7 · The fast way — `across()`