Activity 11 — Linear Regression

Fit a line, read the output, and turn a tracing mass into a leaf area

tidyverse
descriptive-stats

Hands-on companion to the Regression lecture: fit a calibration curve with lm(), read the summary() output, check the residuals, and use predict() to convert tracing masses into leaf areas. Students hand in the finished script with their written answers.

Author

Bill Perry

Published

October 5, 2026

In-class Activity 11: Linear Regression

Today’s Objectives

  1. Fit a linear regression with lm(Y ~ X, data)
  2. Read the slope, intercept, p-value, and R² off the summary() output
  3. Draw the fitted line with geom_smooth(method = "lm")
  4. Check the residuals with plot(model)
  5. Use predict() to turn a tracing mass into a leaf area — one value, then a whole column
  6. Write a standard regression results sentence
Tip📂 Files for today
ImportantWhat you hand in
  • The code is already written in the skeleton — run it chunk by chunk as you work through Parts 1–9 below.
  • What gets graded is your answers, typed into the Q boxes in the script.
  • Hand in scripts/11_regression.R. Run it top to bottom (Ctrl/Cmd + Shift + Enter) first — it must run with no errors.

Answer boxes look like this — type after ANSWER: and start every line with #:

# # # # # # # # # # # # # # # # # # # # # # # # # # # # #
# Q3. What is the slope, with units?
# ANSWER: 130.52 cm2 per gram of paper
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # #

Keep the plots plain today — no labels, colors, or themes. The goal is to read them.


🧩 Chunk 2 — Fit it and read it (after lecture Chunk 2)

Parts 1–4: load the paper data, plot it, fit the line, read the output.

Part 1 · Load libraries and data

▶ Run this (already at the top of your skeleton):

# Load libraries ---------------------------------------
library(readxl)      # read Excel files
library(tidyverse)   # dplyr + ggplot2
library(janitor)     # clean_names()

# Load the paper calibration data ----------------------
paper_df <- read_excel("data/paper_area_weights.xlsx")

glimpse(paper_df)

118 paper squares. area_cm2 is the known area of each square; mass_g is what it weighed.

✏️ Q1: Which column is the explanatory variable (X) and which is the response (Y)? Why that direction, given that we want leaf areas?


Part 2 · Plot it before you fit it

▶ Run this:

# Scatter plot of area against mass --------------------
paper_df %>%
  ggplot(aes(x = mass_g, y = area_cm2)) +
  geom_point()

Look for three things: is it straight, is the scatter the same all the way along, and is any point off on its own?

✏️ Q2: Straight or curved? Any obvious outlier? Is a straight line the right tool here?


Part 3 · Fit the regression

▶ Run this:

# Fit the regression line ------------------------------
paper_lm_model <- lm(area_cm2 ~ mass_g, data = paper_df)

Nothing prints — lm() stores the answer in a model object. The formula reads Y ~ X, “area is predicted by mass.”

▶ Run this:

# Read the model ---------------------------------------
summary(paper_lm_model)

⚠️ Watch out! Flip the formula to lm(mass_g ~ area_cm2) and R runs happily — it just answers a question you didn’t ask. The response always goes on the left of ~.


Part 4 · Read the summary() output

Everything you need is in four places. Find each one in your output.

What you want Where it is
slope (b) Coefficients table, mass_g row, Estimate column
intercept (a) Coefficients table, (Intercept) row, Estimate column
p-value for the slope Coefficients table, mass_g row, Pr(>|t|) column
R² Multiple R-squared near the bottom
F and df the F-statistic: line at the very bottom

✏️ Q4a: Write down the slope and say what it means with units — finish this sentence: “every extra 1 gram of this paper is another ___ cm² of area.”

✏️ Q4b: Write down the intercept. Its p-value is 0.126 — not significant. Why is that the right answer physically? (What is the area of a piece of paper that weighs nothing?)

✏️ Q4c: Write the fitted equation: area = ____ + ____ × mass

✏️ Q4d: Write down the slope’s p-value and R². What does each one tell you — and why are they answering different questions?

💡 Key idea: R² says how tight the fit is. The p-value says whether the slope is real. A relationship can be real and still scattered — that is most of biology.


🧩 Chunk 3 — Plot it, check it, use it (after lecture Chunk 3)

Parts 5–9: draw the line, check the residuals, predict, report.

Part 5 · Draw the fitted line

▶ Run this:

# Scatter plot with the fitted line --------------------
paper_df %>%
  ggplot(aes(x = mass_g, y = area_cm2)) +
  geom_point() +
  geom_smooth(method = "lm")

geom_smooth(method = "lm") adds the regression line plus a grey 95% confidence band. Without method = "lm" you get a wiggly curve instead of a straight line.

✏️ Q5: Where is the grey band widest, and why? (Think about where the data are thickest.)


Part 6 · Check the residuals

A residual is one point’s miss: observed − predicted. plot() on a fitted model draws four residual plots at once — no extra code needed.

▶ Run this:

# All four residual plots in one call ------------------
par(mfrow = c(2, 2))
plot(paper_lm_model)
Panel Name What it checks
1 Residuals vs Fitted straight line? equal scatter?
2 Normal Q-Q are the misses bell-shaped?
3 Scale-Location equal scatter, another view
4 Residuals vs Leverage is one point running the show?

Panel 1 first. A flat shapeless cloud is good. A curve means the relationship isn’t straight. A funnel means the scatter grows as Y grows.

✏️ Q6a: Describe panel 1. Flat cloud, curve, or funnel?

▶ Run this to see the same thing as numbers:

# How big are the misses at each square size? ---------
paper_df %>%
  mutate(residual = residuals(paper_lm_model)) %>%
  group_by(area_cm2) %>%
  summarize(biggest_miss = round(max(abs(residual)), 2))

✏️ Q6b: How big is the worst miss for a 1 cm² square? For a 567 cm² square? Which assumption does that pattern break?

✏️ Q6c: R² is 0.9999 — the best fit you will ever see — and the residual plot still has a problem. In your own words, how can both be true at once?

💡 Key idea: Our leaf tracings weigh 0.1–0.4 g, down at the tight end of the line where a typical miss is a fraction of a cm². The calibration is fine for the job we are using it for. Fixing the funnel properly is a graduate topic — noticing it and saying so is the undergraduate job.


Part 7 · Predict one tracing

▶ Run this:

# Predict the area of one leaf tracing ----------------
one_tracing <- tibble(mass_g = 0.138)

predict(paper_lm_model, newdata = one_tracing)

▶ Run this to add a 95% interval:

# Two tracings, with 95% prediction intervals ---------
two_tracings <- tibble(mass_g = c(0.092, 0.138))

predict(paper_lm_model,
        newdata = two_tracings,
        interval = "prediction")
  • fit is the prediction; lwr and upr are the 95% range
  • interval = "confidence" would give the range for the average tracing of that mass — narrower. For one real leaf, use "prediction".

⚠️ Watch out! newdata must use the same column name as the model. The model was built on mass_g, so the column must be called mass_g. Call it mass and R stops with “object ‘mass_g’ not found.”

✏️ Q7: Record the predicted area for the 0.138 g tracing and its 95% prediction interval. Which tracing is predicted to be larger, and does that make sense?


Part 8 · Convert a whole column of tracings

paper_mass_g in the leaf data is the mass of each leaf’s paper tracing. Now we turn all of them into areas at once.

▶ Run this:

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

# Turn every tracing mass into an area ----------------
leaf_df <- leaf_df %>%
  mutate(leaf_area_cm2 =
           predict(paper_lm_model,
                   newdata = tibble(mass_g = paper_mass_g)))

leaf_df %>%
  select(teams, shade, mass_g, paper_mass_g, leaf_area_cm2)

Same rule as Part 7: predict() wants a column called mass_g, so we hand it tibble(mass_g = paper_mass_g) — our tracing masses, under the name the model knows.

▶ Run this:

# Predicted leaf area by side -------------------------
leaf_df %>%
  ggplot(aes(x = shade, y = leaf_area_cm2)) +
  geom_boxplot() +
  geom_jitter(width = 0.15)
# Mean predicted area per side ------------------------
leaf_df %>%
  group_by(shade) %>%
  summarize(
    n = sum(!is.na(leaf_area_cm2)),
    mean_area = mean(leaf_area_cm2, na.rm = TRUE)
  )

✏️ Q8a: Record the mean predicted area for each side. Does that match what the t-test told you about leaf mass back in Activities 9–10?

✏️ Q8b: Some leaves have NA for leaf_area_cm2. Why? (Look at paper_mass_g for that team.) Is R right to leave them missing?

✏️ Q8c: Why is “shady leaves are about 1 cm² bigger” a more useful sentence for a reader than “shady leaves weigh 0.02 g more”?


Part 9 · Write the results sentence

Fill in the blanks from your summary() output — do not retype numbers from these instructions, read them off your own screen.

____ was a strong predictor of ____ (linear regression:
F(___, ___) = ____, p ____, R2 = ____). The fitted
equation was area = ____ + ____ x mass.

Then add one sentence about the assumptions — what the residual plot showed and what you are being careful about because of it.

✏️ Q9: Write the full results paragraph in your answer box.

Save your script, then run it top to bottom one last time before you hand it in.

Tip✍️ Optional — the same analysis as a report

If you’d rather hand this in as a Word document than a script, regression_analysis_report.qmd is the same analysis written up in Quarto. Its Discussion section is waiting for the paragraph you just wrote; paste it in and press Render. The .R script is still what’s graded.


Extension — out of class (~20–30 min)

Turn this in with your script, as # comments at the bottom of the file. paper_df is the same for everyone — but the object you calibrate is one only you have.

E1 · Calibrate something only you have

  1. Cut any shape you like out of ordinary printer or notebook paper — trace your hand, your initials in bubble letters, a rough outline of your home state.
  2. Weigh it in grams on a kitchen or postal scale. No scale? Estimate it against a US penny (2.5 g) or a sheet of printer paper (≈4.5 g), and say that’s what you did.
  3. In R, use paper_lm_model and predict(..., interval = "prediction") to estimate its area from its mass.
  4. Measure the real area the low-tech way: trace it onto graph paper and count squares.

Report, as comments: what the shape was, its mass, the predicted area and interval, your counted area, and whether the counted area landed inside the interval.

E2 · Explain your own numbers

  1. Using your R² from Part 4, say in plain language what fraction of the variation in area the line accounts for — and one thing R² does not tell you.
  2. Your prediction interval in E1 is wider than a confidence interval would be at the same mass. Why is predicting one specific object riskier than estimating an average?
  3. Run max(paper_df$mass_g). Suppose you had cut your shape out of cardboard and it weighed 9 g. Why would predicting its area be a different kind of mistake than the prediction you actually made?

📖 Whitlock & Schluter §17.2 p. 550: “Extrapolation is the prediction of Y at values of X beyond the range of X-values in the data. Extrapolation is problematic because there is no way to ensure the relationship between X and Y continues to be linear.”


End of Activity 11. Next: real climate data with GSODR — and the same lm(), where the slope becomes a rate of change per year.