BIO320 Final Project
Carly Blundi and Sofia Chavarria
2026-04-15
Introduction to Golem Grad Tortoises
Female Hermann’s tortoises on Golem Grad (“Snake Island”) in North Macedonia are facing an ongoing extinction event due to aggressive mating behaviors by male tortoises (Arsovski et al., 2026). Sexually aggressive males on the island greatly outnumber the females, and inflict severe injuries on the females during copulation; these injuries put the females at a high risk of fatal falls from the island’s plateau (Golubović et al., 2025). The harassed females are emaciated, reproduce less frequently, and have lower survival rates than females from the neighboring mainland population. This study accumulates sixteen years of recapture data and predicts that the last island female will die in 2083 (Arsovski et al., 2026). Despite the high population, a highly skewed sex ratio can lead to a collapse in population over time. By comparing population data from both the mainland and the Golem Grad island population, we aim to assess if the island tortoises have smaller carapaces and smaller body mass compared to the mainland tortoises due to the ongoing aggression associated with mating.
Our Research Question
Do the female Golem Grad Island tortoises have a significantly lower body mass and carapace length compared to their mainland tortoise counterparts?
Female Golem Grad island tortoise caraspace injury after falling
Loading the data:
knitr::opts_chunk$set(echo = TRUE)
tortoise_body_condition_cleaned <- readr::read_csv('https://raw.githubusercontent.com/rfordatascience/tidytuesday/main/data/2026/2026-03-03/tortoise_body_condition_cleaned.csv')## Rows: 10174 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (4): individual, season, locality, sex
## dbl (5): year, year_recode, body_mass_grams, body_condition_index, straight_...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.library(dplyr)##
## Attaching package: 'dplyr'
##
## The following objects are masked from 'package:stats':
##
## filter, lag
##
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(ggplot2)
library(car)## Loading required package: carData
##
## Attaching package: 'car'
##
## The following object is masked from 'package:dplyr':
##
## recodelibrary(knitr)# Filtering out males, combining beach and plateu into one locality "Mainland", grouping by individual and locality to find average body mass and caraspace length for each individual
# Filter females, recode locality, and average recapture data per individual
female_tortoises_avg <- tortoise_body_condition_cleaned %>%
filter(sex != "m") %>%
mutate(locality = case_when(
locality %in% c("Beach", "Plateau") ~ "Mainland",
locality == "Konjsko" ~ "Island",
TRUE ~ locality
)) %>%
group_by(individual, locality) %>%
summarise(
body_mass_grams = mean(body_mass_grams, na.rm = TRUE),
straight_carapace_length_mm = mean(straight_carapace_length_mm, na.rm = TRUE),
.groups = "drop"
)
#Histograms to visualizes the body mass of female tortoises based on their locality
ggplot(female_tortoises_avg, aes(x = body_mass_grams, fill = locality)) +
geom_histogram(color = "black", bins = 15) +
facet_wrap(~locality, ncol = 1) +
scale_fill_manual(values = c("Mainland" = "steelblue", "Island" = "darkolivegreen3")) +
labs(
x = "Body Mass (grams)",
y = "Count",
title = "Body Mass Distribution of Female Tortoises by Locality",
fill = "Locality"
) +
theme_minimal()
summary(female_tortoises_avg)## individual locality body_mass_grams
## Length:646 Length:646 Min. : 676.3
## Class :character Class :character 1st Qu.:1200.8
## Mode :character Mode :character Median :1482.6
## Mean :1482.5
## 3rd Qu.:1771.0
## Max. :2748.0
## straight_carapace_length_mm
## Min. :151.0
## 1st Qu.:178.5
## Median :193.0
## Mean :191.3
## 3rd Qu.:204.8
## Max. :243.5# Normality testing for Island tortoises
ggplot(female_tortoises_avg %>% filter(locality == "Island"),
aes(x = body_mass_grams)) +
geom_histogram(color = "black", fill = "steelblue", bins = 15) +
labs(
x = "Body Mass (grams)",
y = "Count",
title = "Body Mass Distribution of Island Female Tortoises"
) +
theme_minimal()
# Base R qqnorm for Island
qqnorm(female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Island"],
main = "Q-Q Plot - Island Female Tortoises",
col = "steelblue");qqline(female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Island"])
shapiro.test(female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Island"]) #pvalue =0.01721##
## Shapiro-Wilk normality test
##
## data: female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Island"]
## W = 0.99304, p-value = 0.01721# Normality testing for Mainland tortoises
ggplot(female_tortoises_avg %>% filter(locality == "Mainland"),
aes(x = body_mass_grams)) +
geom_histogram(color = "black", fill = "darkolivegreen3", bins = 15) +
labs(
x = "Body Mass (grams)",
y = "Count",
title = "Body Mass Distribution of Mainland Female Tortoises"
) +
theme_minimal()
qqnorm(female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Mainland"],
main = "Q-Q Plot - Mainland Female Tortoises",
col = "darkolivegreen3");qqline(female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Mainland"])
shapiro.test(female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Mainland"]) #pvalue 0.01193##
## Shapiro-Wilk normality test
##
## data: female_tortoises_avg$body_mass_grams[female_tortoises_avg$locality == "Mainland"]
## W = 0.97349, p-value = 0.01193# Add log transformed body mass column
female_tortoises_avg <- female_tortoises_avg %>%
mutate(log_body_mass_grams = log(body_mass_grams))
# Island tortoises only - log transformed
ggplot(female_tortoises_avg %>% filter(locality == "Island"),
aes(x = log_body_mass_grams)) +
geom_histogram(color = "black", fill = "steelblue", bins = 15) +
labs(
x = "Log Body Mass (grams)",
y = "Count",
title = "Log Transformed Body Mass Distribution - Island Female Tortoises"
) +
theme_classic()
ggplot(female_tortoises_avg %>% filter(locality == "Island"),
aes(sample = log_body_mass_grams)) +
stat_qq(color = "steelblue") +
stat_qq_line(color = "black") +
labs(
x = "Theoretical Quantiles",
y = "Sample Quantiles",
title = "Q-Q Plot Log Body Mass - Island Female Tortoises"
) +
theme_classic()
# Shapiro-Wilk test for Island tortoises - log transformed
shapiro.test(female_tortoises_avg$log_body_mass_grams[female_tortoises_avg$locality == "Island"])##
## Shapiro-Wilk normality test
##
## data: female_tortoises_avg$log_body_mass_grams[female_tortoises_avg$locality == "Island"]
## W = 0.9879, p-value = 0.0002838# Mainland tortoises only - log transformed
ggplot(female_tortoises_avg %>% filter(locality == "Mainland"),
aes(x = log_body_mass_grams)) +
geom_histogram(color = "black", fill = "darkolivegreen3", bins = 15) +
labs(
x = "Log Body Mass (grams)",
y = "Count",
title = "Log Transformed Body Mass Distribution - Mainland Female Tortoises"
) +
theme_classic()
ggplot(female_tortoises_avg %>% filter(locality == "Mainland"),
aes(sample = log_body_mass_grams)) +
stat_qq(color = "darkolivegreen3") +
stat_qq_line(color = "black") +
labs(
x = "Theoretical Quantiles",
y = "Sample Quantiles",
title = "Q-Q Plot Log Body Mass - Mainland Female Tortoises"
) +
theme_classic()
shapiro.test(female_tortoises_avg$log_body_mass_grams[female_tortoises_avg$locality == "Mainland"])##
## Shapiro-Wilk normality test
##
## data: female_tortoises_avg$log_body_mass_grams[female_tortoises_avg$locality == "Mainland"]
## W = 0.99102, p-value = 0.5701wilcox.test(body_mass_grams ~ locality, data = female_tortoises_avg)##
## Wilcoxon rank sum test with continuity correction
##
## data: body_mass_grams by locality
## W = 61216, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0# Female Tortoise Straight Carapace Length Analysis
localities <- c("Island", "Mainland")
hist_colors <- c("Island" = "steelblue", "Mainland" = "darkolivegreen3")
for (loc in localities) {
df_loc <- female_tortoises_avg %>% filter(locality == loc)
# Histogram
print(
ggplot(df_loc, aes(x = straight_carapace_length_mm)) +
geom_histogram(color = "black", fill = hist_colors[loc], bins = 15) +
labs(
x = "Straight Carapace Length (mm)",
y = "Count",
title = paste("Straight Carapace Length Distribution -", loc, "Female Tortoises")
) +
theme_classic()
)
# Q-Q Plot
qqnorm(df_loc$straight_carapace_length_mm,
main = paste("Q-Q Plot -", loc, "Female Tortoises"),
col = hist_colors[loc])
qqline(df_loc$straight_carapace_length_mm)
# Shapiro-Wilk Test
cat("\nShapiro-Wilk Test —", loc, "\n")
print(shapiro.test(df_loc$straight_carapace_length_mm))
}
##
## Shapiro-Wilk Test — Island
##
## Shapiro-Wilk normality test
##
## data: df_loc$straight_carapace_length_mm
## W = 0.99316, p-value = 0.01895
##
## Shapiro-Wilk Test — Mainland
##
## Shapiro-Wilk normality test
##
## data: df_loc$straight_carapace_length_mm
## W = 0.97178, p-value = 0.008188# --- 1. Levene's Test: Body Mass (grams) ---
cat("Levene's Test — Body Mass (grams)\n")## Levene's Test — Body Mass (grams)leveneTest(body_mass_grams ~ locality, data = female_tortoises_avg)## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 1 34.437 7.048e-09 ***
## 644
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# --- 2. Levene's Test: Straight Carapace Length (mm) ---
cat("\nLevene's Test — Straight Carapace Length (mm)\n")##
## Levene's Test — Straight Carapace Length (mm)leveneTest(straight_carapace_length_mm ~ locality, data = female_tortoises_avg)## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 1 3.4463 0.06385 .
## 644
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1female_tortoises_avg <- female_tortoises_avg %>%
mutate(log_straight_carapace_length_mm = log(straight_carapace_length_mm))
for (loc in localities) {
df_loc <- female_tortoises_avg %>% filter(locality == loc)
# Histogram — log scale
print(
ggplot(df_loc, aes(x = log_straight_carapace_length_mm)) +
geom_histogram(color = "black", fill = hist_colors[loc], bins = 15) +
labs(
x = "Log Straight Carapace Length (mm)",
y = "Count",
title = paste("Log Straight Carapace Length Distribution -", loc, "Female Tortoises")
) +
theme_classic()
)
# Q-Q Plot — log scale
print(
ggplot(df_loc, aes(sample = log_straight_carapace_length_mm)) +
stat_qq(color = hist_colors[loc]) +
stat_qq_line(color = "black") +
labs(
x = "Theoretical Quantiles",
y = "Sample Quantiles",
title = paste("Q-Q Plot Log Straight Carapace Length -", loc, "Female Tortoises")
) +
theme_classic()
)
# Shapiro-Wilk Test — log scale
cat("\nShapiro-Wilk Test (log) —", loc, "\n")
print(shapiro.test(df_loc$log_straight_carapace_length_mm))
}
##
## Shapiro-Wilk Test (log) — Island
##
## Shapiro-Wilk normality test
##
## data: df_loc$log_straight_carapace_length_mm
## W = 0.98392, p-value = 1.795e-05
##
## Shapiro-Wilk Test (log) — Mainland
##
## Shapiro-Wilk normality test
##
## data: df_loc$log_straight_carapace_length_mm
## W = 0.97701, p-value = 0.02622wilcox.test(straight_carapace_length_mm ~ locality, data = female_tortoises_avg)##
## Wilcoxon rank sum test with continuity correction
##
## data: straight_carapace_length_mm by locality
## W = 57312, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0Discussion
After analyzing the data, we found that the female island tortoises have a significantly lower body mass when compared to their mainland tortoise counterparts (Mann-Whitney test: W= 61,216, p=2.2e-16). Similarly, the female island tortoises also have a significantly shorter carapace length in comparison to the mainland tortoises(Mann-Whitney test: W= 57,312, p=2.2e-16). Based on these results, we can reject the null hypothesis. These findings provide insight into the ecological context of the female Hermann’s tortoises based on location. The statistical analysis of the recapture data demonstrates that, on average, the island female tortoises have a lower body mass as well as a shorter carapace length compared to the mainland female tortoises. These findings suggest that the sexually aggressive behavior of the male tortoises impacts the physical health of the female tortoises that inhabit the island. Since the data were not normally distributed and did not have equal variance, we used the nonparametric alternative to the two-sample t-test.
Conclusion
In conclusion, the data analysis shows that male aggression significantly alters the average body mass and carapace length of females on Golem Grad Island. The Mann-Whitney showed significant results (2.2e-16, W=9767.5), which displays that the aggressive mating behavior of the males on the island has led to a decrease in the overall well-being of the female tortoises on the island, in comparison to the mainland female tortoises. These findings are meaningful because it further supports the claim that the last female tortoise on the island will likely die out by 2083 (Arsovski et al., 2026) and that it is critical to protect the female Golem Grad tortoises to preserve the population.
References
Arsovski D, Bonnet X, Golubović A, Tomović L. Sex Ratio Bias Triggers Demographic Suicide in a Dense Tortoise Population. Ecol Lett. 2026 Jan;29(1):e70296. doi: 10.1111/ele.70296. PMID: 41587391.
Golubović, Ana & Arsovski, Dragan & Bjelica, Vukašin & Bonnet, Xavier & Tomović, Ljiljana. (2025). Coercive mating in
Arsovski Dragon, Hermann’s tortoise (Testudo hermanni) on Golem Grad island affects personality traits of harassed females. A female Hermann’s tortoise from the island of Golem Grad in North Macedonia, with an injury from a fall on her carapace. The New York Times. https://www.nytimes.com/2026/02/14/science/tortoises-island-sex-cliff.html
Anthropic. (2026). Claude (Claude Sonnet 4.6) [Large language model]. https://www.claude.ai