Unit 7 Project
Donatella
3/9/2026
Scenario: Inference for Penguin Data using Randomization and Bootstrapping
We are interested in the Palmer Penguins data set, and we will use two inferential tools to study group differences. In this project, we will work with a randomization test to ask whether an observed difference between two groups is real or plausibly due to chance, and a bootstrap confidence interval to estimate how large that difference plausibly is.
We will do this twice. First, we will work through a complete example together: studying the difference in body mass between male and female Gentoo penguins. Then, you will carry out the same analysis yourself for a second comparison: bill depth in Adelie versus Chinstrap penguins.
Part 1: Body Mass in Male vs. Female Gentoo Penguins
Setting up the data
We begin by filtering the penguins data to Gentoo penguins only, keeping only rows with complete data for body mass and sex.
gentoo <- penguins %>%
filter(species == "Gentoo", !is.na(body_mass_g), !is.na(sex))
nrow(gentoo)## [1] 119table(gentoo$sex)##
## female male
## 58 61We can visualize the raw data to get a sense of whether a difference exists before doing any inference.
ggplot(gentoo, aes(x = sex, y = body_mass_g, fill = sex)) +
geom_boxplot(alpha = 0.7) +
labs(title = "Body mass of Gentoo penguins by sex",
x = "Sex", y = "Body mass (g)") +
theme_minimal() +
theme(legend.position = "none")
There appears to be a visible difference. But could this have arisen by chance? We use a randomization test to find out.
The Observed Difference
Our test statistic is the difference in mean body mass: male mean minus female mean.
gentoo %>%
group_by(sex) %>%
summarise(mean_body_mass = mean(body_mass_g))## # A tibble: 2 × 2
## sex mean_body_mass
## <fct> <dbl>
## 1 female 4680.
## 2 male 5485.observed_diff <- mean((gentoo %>% filter(sex == "male"))$body_mass_g) -
mean((gentoo %>% filter(sex == "female"))$body_mass_g)
observed_diff## [1] 805.0947The observed difference in mean body mass is approximately 805 grams, with males heavier on average.
The Randomization Test
To ask whether this difference is real, we simulate a world in which sex has no relationship to body mass. We do this by shuffling the sex labels and recomputing the difference in means, repeating the process 1,000 times to build a null distribution.
one_shuffle_diff <- function() {
gentoo %>%
mutate(sex_shuffle = sample(sex)) %>%
group_by(sex_shuffle) %>%
summarise(mean_mass = mean(body_mass_g)) %>%
summarise(diff = mean_mass[sex_shuffle == "male"] -
mean_mass[sex_shuffle == "female"]) %>%
pull(diff)
}
null_distribution <- replicate(1000, one_shuffle_diff())We can visualize the null distribution and mark where our observed difference falls.
null_df <- data.frame(sim_diff = null_distribution)
ggplot(null_df, aes(x = sim_diff)) +
geom_histogram(bins = 40, fill = "steelblue", color = "white", alpha = 0.8) +
geom_vline(xintercept = observed_diff, color = "firebrick", linewidth = 1.2,
linetype = "dashed") +
annotate("text", x = observed_diff - 200, y = 60,
label = paste0("Observed\nDifference\n", round(observed_diff, 0), "g"),
color = "firebrick", hjust = 0, size = 3.5) +
labs(title = "Null distribution: 1,000 shuffles of sex labels",
x = "Simulated difference in mean body mass (g)",
y = "Count") +
theme_minimal()
The observed difference falls far out in the tail of the null distribution — none of the 1,000 shuffles came close to producing a difference as large as what we observed.
We can quantify this with an informal p-value: the proportion of simulated differences that were as large or larger than what we observed.
p_value <- sum(null_distribution >= observed_diff) / length(null_distribution)
p_value## [1] 0The p-value is 0. Out of 1,000 random shuffles, zero produced a difference as large as the observed difference of 805 grams. This is strong evidence that the difference in body mass between male and female Gentoo penguins is not due to chance.
The Bootstrap Confidence Interval
The randomization test tells us the difference is real. Now we ask: how large is it? We use bootstrapping to estimate a plausible range for the true difference. Unlike the randomization test, bootstrapping resamples from the real data with replacement.
gentoo_males <- gentoo %>% filter(sex == "male")
gentoo_females <- gentoo %>% filter(sex == "female")
one_boot_diff <- function() {
boot_males <- sample(gentoo_males$body_mass_g, size = nrow(gentoo_males), replace = TRUE)
boot_females <- sample(gentoo_females$body_mass_g, size = nrow(gentoo_females), replace = TRUE)
mean(boot_males) - mean(boot_females)
}
boot_diffs <- replicate(1000, one_boot_diff())Once we have our bootstrapped distribution, we can look at it, and calculate the 95% confidence interval by chopping off the most extreme values (2.5% on each side).
data.frame(boot_diff = boot_diffs) %>%
ggplot(aes(x = boot_diff)) +
geom_histogram(bins = 40, fill = "steelblue", color = "white") +
geom_vline(xintercept = observed_diff, color = "firebrick", linewidth = 1.2) +
labs(title = "Bootstrap distribution: difference in mean body mass",
subtitle = "Red = observed difference",
x = "Bootstrap difference in mean body mass (male - female)",
y = "Count") +
theme_minimal()
ci_diff <- quantile(boot_diffs, probs = c(0.025, 0.975))
ci_diff## 2.5% 97.5%
## 703.0469 915.8022The 95% confidence interval runs from approximately 710g to 910g (although exact values will change depending on the randomness in our specific simulation). We are 95% confident that male Gentoo penguins are, on average, somewhere between 710 and 910 grams heavier than female Gentoo penguins. Notice that this interval does not include zero, which is consistent with our randomization test conclusion that the difference is real.
Part 2: Do a Similar Analysis for another Penguin Species!
Now it is your turn. You will carry out the same two-part analysis for a different comparison. Pick a different species (Adelie or Chinstrap) and a quantitative variable (body mass is fine to do, but you may find more interesting results with bill depth, bill length, or flipper length) and address the following questions:
library(dplyr)
library(ggplot2)
library(palmerpenguins)
data("penguins")
adelie <- penguins %>%
filter(species == "Adelie", !is.na(bill_depth_mm), !is.na(sex))
nrow(adelie)## [1] 146table(adelie$sex)##
## female male
## 73 73Question A: Consider the difference in mean between the male and female penguins of your chosen species. Is that difference real, or is it plausibly due to chance?
adelie %>%
group_by(sex) %>%
summarise(mean_bill_depth = mean(bill_depth_mm))## # A tibble: 2 × 2
## sex mean_bill_depth
## <fct> <dbl>
## 1 female 17.6
## 2 male 19.1observed_diff <- mean((adelie %>% filter(sex=="male"))$bill_depth_mm) -
mean((adelie %>% filter(sex=="female"))$bill_depth_mm)
observed_diff## [1] 1.450685null_df <- data.frame(sim_diff = null_distribution)
ggplot(null_df, aes(x = sim_diff)) +
geom_histogram(bins = 40, fill = "steelblue", color = "white") +
geom_vline(xintercept = observed_diff, color = "red", linewidth = 1.2) +
labs(title = "Null distribution: shuffled bill depth differences",
x = "Simulated difference in mean bill depth (male - female)",
y = "Count") +
theme_minimal()
p_value <- sum(null_distribution >= observed_diff) / length(null_distribution)
p_value## [1] 0.491adelie_males <- adelie %>% filter(sex == "male")
adelie_females <- adelie %>% filter(sex == "female")one_boot_diff <- function() {
boot_males <- sample(adelie_males$bill_depth_mm,
size = nrow(adelie_males),
replace = TRUE)
boot_females <- sample(adelie_females$bill_depth_mm,
size = nrow(adelie_females),
replace = TRUE)
mean(boot_males) - mean(boot_females)
}data.frame(boot_diff = boot_diffs) %>%
ggplot(aes(x = boot_diff)) +
geom_histogram(bins = 40, fill = "purple", color = "white") +
geom_vline(xintercept = observed_diff, color = "red", linewidth = 1.2) +
labs(title = "Bootstrap distribution of bill depth difference",
x = "Bootstrapped difference in mean bill depth",
y = "Count") +
theme_minimal()
ci_diff <- quantile(boot_diffs, probs = c(0.025, 0.975))
ci_diff## 2.5% 97.5%
## 703.0469 915.8022The histogram shows there is a difference unlikely to happen by random chance. There is also a difference between male and female Adelie penguins.
Bootstrap distribution shows a range of possibilities in difference for bill depth.