Correlation and Linear Regression

Association vs prediction — r, a fitted line, R², and Cook’s distance on Adelie penguins

Everything on this page runs real R in your browser. You will practice Coding Assignment 4, Part D plus lecture extras on correlation, R², and Cook’s distance.

NoteAssignment vs this activity

On DataHub, Assignment 4 regresses biomass_g on leaf_area_cm2 in the brittlebush dataset. Here we ask whether flipper length predicts body mass in Adelie penguins only — so species clustering does not inflate R². The workflow — scatterplot, lm(), residual plot, R² — matches the assignment; cor.test() and Cook’s distance come from lecture A010.

NoteCorrelation vs regression
  • Correlation (r): strength and direction of a linear association in the data you measured. Use cor.test() when both variables are numeric and you want r and a p-value.
  • Regression: fit a line to predict Y from X. Check assumptions after fitting with residual plots. = proportion of variance in Y explained by the model (r² when there is one predictor).

Unlike a t-test, you fit the model first, then diagnose it.

Research question: Among Adelie penguins, does flipper length predict body mass?

1. Load packages

NoteHint

Load tidyverse and palmerpenguins. You will load broom and ggrepel later for the Cook’s distance plot.

TipSolution
library(tidyverse)
library(palmerpenguins)

2. Build the Adelie table

Keep Adelie penguins with non-missing body_mass_g and flipper_length_mm — use !is.na() on both columns in filter(). Save as penguin_data.

NoteHint

Filter penguins to Adelie rows with non-missing body_mass_g and flipper_length_mm.

Example:

penguins |>
  filter(species == "Gentoo", !is.na(body_mass_g))
TipSolution
penguin_data <- penguins |>
  filter(species == "Adelie", !is.na(body_mass_g), !is.na(flipper_length_mm))

151 Adelie penguins — complete measurements only.

3. Scatterplot with trend line

Plot body_mass_g (y) vs flipper_length_mm (x). Add a linear trend with geom_smooth(method = "lm"). Save as scatter_lm and display it.

NoteHint

Use ggplot() with geom_point() and geom_smooth(method = "lm"). Map flipper length to x and body mass to y.

Example:

ggplot(mtcars, aes(x = wt, y = mpg)) +
  geom_point() +
  geom_smooth(method = "lm")
TipSolution
scatter_lm <- ggplot(penguin_data, aes(x = flipper_length_mm, y = body_mass_g)) +
  geom_point() +
  geom_smooth(method = "lm")
scatter_lm

Positive slope — heavier Adelie penguins tend to have longer flippers.

4. Pearson correlation

Run cor.test() on flipper length and body mass. Save as cor_result and print it — you need both r and a p-value.

NoteHint

Use cor.test() on two numeric vectors from the same data frame.

Example:

cor.test(mtcars$mpg, mtcars$wt)
TipSolution
cor_result <- cor.test(penguin_data$flipper_length_mm, penguin_data$body_mass_g)
cor_result

r ≈ 0.47, p ≪ 0.001 — moderate positive association, statistically significant.

5. Fit linear model

Fit lm(body_mass_g ~ flipper_length_mm, data = penguin_data). Save as lm_model and print summary(lm_model).

NoteHint

Fit with lm(response ~ predictor, data = ...) and inspect with summary().

Example:

fit <- lm(mpg ~ wt, data = mtcars)
summary(fit)
TipSolution
lm_model <- lm(body_mass_g ~ flipper_length_mm, data = penguin_data)
summary(lm_model)

Slope for flipper length is positive and highly significant.

6. Residuals vs fitted

Use plot(lm_model, which = 1) to display the residuals-versus-fitted plot.

NoteHint

After fitting a linear model, use plot(model, which = 1) to check residual spread.

Example:

fit <- lm(mpg ~ wt, data = mtcars)
plot(fit, which = 1)

Look for residuals scattered roughly evenly around zero with no obvious curve.

TipSolution
plot(lm_model, which = 1)

7. Cook’s distance

Load broom and ggrepel, then build the Cook’s distance plot from lecture A010. Flag cases above 4 / nobs(lm_model).

NoteHint

Layers to add after ggplot(aes(case, .cooksd)):

  • geom_segment() from each point down to zero
  • geom_hline() at cook_cut (dashed)
  • geom_point() for all cases
  • red geom_point() for rows in influential
  • geom_text_repel() to label influential case numbers
  • labs() for title and axis names

Example segment plot on a small made-up table:

tibble(case = 1:5, value = c(0.01, 0.02, 0.15, 0.01, 0.03)) |>
  ggplot(aes(case, value)) +
  geom_segment(aes(xend = case, yend = 0)) +
  geom_point()
TipSolution
library(broom)
library(ggrepel)

aug <- augment(lm_model) |>
  mutate(case = row_number())

cook_cut <- 4 / nobs(lm_model)
influential <- aug |>
  filter(.cooksd > cook_cut)

p_cooks <- aug |>
  ggplot(aes(case, .cooksd)) +
  geom_segment(aes(xend = case, yend = 0)) +
  geom_hline(yintercept = cook_cut, linetype = "dashed") +
  geom_point() +
  theme_minimal() +
  geom_point(data = influential, color = "red") +
  geom_text_repel(data = influential, aes(label = case), size = 3) +
  labs(
    title = "Cook's Distance by Case",
    x = "Case index",
    y = "Cook's distance"
  )
p_cooks

Six cases exceed the cutoff — they pull the line more than typical points. Flagging ≠ deleting unless you have a biological reason.

NoteOutliers vs influence

A point can be an outlier in Y without strongly influencing the slope. Cook’s distance combines leverage and residual size. The 4/n cutoff is a screening rule — investigate flagged cases, but do not winsorize or trim without justification.

8. R-squared

Extract R² from summary(lm_model). Save as lm_r_squared and print it.

NoteHint

Extract R² from the model summary object with $r.squared.

Example:

summary(lm(mpg ~ wt, data = mtcars))$r.squared
TipSolution
lm_r_squared <- summary(lm_model)$r.squared
lm_r_squared

R² ≈ 0.22 — flipper length explains about 22% of body-mass variation among Adelie penguins. r ≈ 0.47, so r² ≈ 0.22 — same information, different scale.

9. Interpret the model

Assign slope_significant <- TRUE and keep_influential_points <- TRUE (flag but do not delete influential cases without cause).

NoteHint

Read summary(lm_model) for the slope p-value. Inspect p_cooks — are there influential points worth flagging? Assign TRUE or FALSE for each variable.

TipSolution
slope_significant <- TRUE
keep_influential_points <- TRUE

Residual plot looks reasonable; slope p ≪ 0.001; r and R² tell a consistent moderate-positive story.

10. Recap — functions you practiced

Function Role
ggplot() + geom_point() + geom_smooth(method = "lm") Scatter with fitted line
cor.test(x, y) Pearson r and p-value for linear association
lm(y ~ x, data = ...) Fit simple linear regression
summary(lm_model) Coefficients, R², and p-values
plot(model, which = 1) Residuals vs fitted
augment(model) Add residuals, fitted values, Cook’s d (.cooksd)
4 / nobs(model) Common Cook’s distance cutoff
geom_text_repel() Label influential case indices

If you can explain each row, you have covered Coding Assignment 4, Part D plus lecture A010 extras.

Keep playing

Pool all three species and plot body mass vs flipper length colored by species. Why does R² jump when you ignore species?

R² ≈ 0.76 pooled vs ≈ 0.22 within Adelie — much of the “prediction” was really species clustering, not flipper length alone. That is why we regressed within one species.