Gulf of America/Mexico Loggerhead Sea Turtle Age & Growth
Olivia Trahan Dept. Honor’s Project
Created: 2026-02-09
Population and Food Web Ecology Lab – Dr. Matthew Ramirez – UNCW
R Markdown Information
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. These files are much easier to share and follow the code and outputs.
Load in Packages & Setup Your Environment
Change the directory to where your files are located. This establishes where your R markdown file and project are saved. In the bottom right panel, you can view the files within your directory. In the top right, you can view your environment where your data frames and their variables are located (helpful to view what data type the variables are).
## Import packages and setup directory
## Only run these two lines once
#install.packages("patchwork")
#install.packages("ggplot2")
library(librarian)
shelf(dplyr, readxl, writexl, openxlsx, here, pls, ggpubr, mdatools, FSA, purrr, data.table, tibble, ggplot2, patchwork, mgcv, quiet = TRUE)
## Set working directory
setwd("~/Desktop/UNCW_PhD Docs/Ramirez_Skeletochronology/Olivia Trahan_Proj")Load in Metadata
Saved the ‘age and growth’ tab of the spreadsheet ‘Cc GoA samples & consensus reads’ as a .csv. This reduces the file to a simple plain-text format where each row is a new line and fields (like columns in a spreadsheet) are separated by commas. Great for import into R as a data.frame, but it cannot save more than one Excel tabs.
## Read in metadata as a .csv file
meta.cc <- read.csv(file="Cc GoA_Age and Growth.csv", header=T)
## Display the column headers
head(meta.cc)## Region Species STSSN Other.ID.. SCL Stranding.Date Sex State
## 1 GoMx Cc AMO20220430-01 71 4/30/22 F TX
## 2 GoMx Cc AMO20220430-01 71 4/30/22 F TX
## 3 GoMx Cc AMO20220430-01 71 4/30/22 F TX
## 4 GoMx Cc AMO20220430-01 71 4/30/22 F TX
## 5 GoMx Cc AMO20220430-01 71 4/30/22 F TX
## 6 GoMx Cc AMO20220430-01 71 4/30/22 F TX
## Nutritional_Status Age_at_stranding SCL_initial SCL_final Mean_SCL
## 1 1 16.75 50.3 53.8 52.0
## 2 1 16.75 53.8 59.4 56.6
## 3 1 16.75 59.4 63.0 61.2
## 4 1 16.75 63.0 64.9 64.0
## 5 1 16.75 64.9 66.9 65.9
## 6 1 16.75 66.9 68.0 67.4
## Growth_Rate Age Year Year_Before_Death
## 1 3.5 7.75 2013 -17
## 2 5.6 8.75 2014 -16
## 3 3.7 9.75 2015 -15
## 4 1.9 10.75 2016 -14
## 5 2.0 11.75 2017 -13
## 6 1.1 12.75 2018 -12## Convert Year to numeric
## May be stored as a factor/character, in which case min() and max() operate on factor levels instead of true numeric values, leading to incorrect year ranges.
meta.cc <- meta.cc %>%
mutate(Year = as.numeric(as.character(Year)))
## Round ages to nearest whole number, creating 'Age' as a new column with the rounded Ages and adding 'Age_group' as a column
meta.cc <- meta.cc %>%
mutate(
Age_round = round(Age),
Age_group = case_when(
Age_round >= 3 & Age_round <= 7 ~ "Age 3-7",
Age_round >= 8 & Age_round <= 11 ~ "Age 8-11",
Age_round >= 12 & Age_round <= 15 ~ "Age 12-15",
Age_round >= 16 & Age_round <= 19 ~ "Age 16-19",
Age_round >= 20 & Age_round <= 24 ~ "Age 20-24"))
meta.cc <- meta.cc %>%
mutate(Nutritional_Status_code = Nutritional_Status) %>% # keep numeric code
mutate(Nutritional_Status = case_when(
Nutritional_Status_code %in% c(1,2) ~ "Emaciated",
Nutritional_Status_code == 3 ~ "Good",
Nutritional_Status_code == 4 ~ "Robust"
)) %>%
mutate(Nutritional_Status = factor(Nutritional_Status,
levels = c("Emaciated", "Good", "Robust")))
## Filter to remove NAs in Growth Rate, Year Before Death and Age_Round
meta.cc_clean <- meta.cc %>% filter(!is.na(Growth_Rate) & !is.na(Year_Before_Death) & !is.na(Age_round))
## Filter for Year_Before_Death between -10 and 0
meta.cc_clean_10y <- meta.cc %>%
filter(Year_Before_Death >= -10 & Year_Before_Death <= 0)
## Compute mean and SE by Year_Before_Death and Nutritional_Status
summary_meta <- meta.cc %>%
group_by(Year_Before_Death, Nutritional_Status) %>%
summarise(
Mean = mean(Growth_Rate, na.rm = TRUE),
SD = sd(Growth_Rate, na.rm = TRUE),
N = sum(!is.na(Growth_Rate)),
.groups = 'drop'
) %>%
mutate(SE = SD / sqrt(N),
CI_lower = Mean - 1.96*SE,
CI_upper = Mean + 1.96*SE)
## Export meta.cc_clean
write_xlsx(meta.cc_clean, "meta.cc_clean.xlsx")Create Summary Table - by Age_at_stranding
Group by Age_at_stranding and changed sig. figs. to 3 yet excluded the Year to provide correct ranges (min and max).
## Uses dplyr (group_by, mutate, summarise)
summary_table_age <- meta.cc %>%
group_by(Age_at_stranding) %>%
summarise(
N = sum(!is.na(Growth_Rate)),
Mean = mean(Growth_Rate, na.rm = TRUE),
SD = sd(Growth_Rate, na.rm = TRUE),
Year_min = min(Year, na.rm = TRUE),
Year_max = max(Year, na.rm = TRUE)
) %>%
ungroup() %>%
mutate(across(c(N, Mean, SD), ~ signif(., 4)))
## Display the created table
summary_table_age## # A tibble: 11 × 6
## Age_at_stranding N Mean SD Year_min Year_max
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 11.8 8 2.82 2.37 2014 2021
## 2 12.8 6 2.4 2.22 2013 2021
## 3 13.8 11 1.96 1.49 2012 2022
## 4 14.8 15 2.83 2.63 2012 2021
## 5 15.8 36 2.08 2.01 2009 2021
## 6 16.8 43 2.65 1.82 2010 2022
## 7 17.8 22 2.06 1.75 2009 2021
## 8 18.8 13 2.46 1.51 2009 2021
## 9 20.8 13 1.85 1.75 2007 2021
## 10 22.8 17 1.24 1.34 2005 2021
## 11 23.8 14 0.757 1.37 2009 2022Create Summary Table - by ‘SCL_final’ grouped per 10 cm SCL for Florida & Texas
Group by SCL_final and changed sig. figs. to 3 yet excluded the Year to provide correct ranges (min and max).
## Creates a new variable called 'SCL_bin' from SCL_final per size bin.
meta.cc <- meta.cc %>%
mutate(
SCL_bin = cut(
SCL_final,
breaks = seq(10, ceiling(max(SCL_final, na.rm = TRUE) / 10) * 10, by = 10),
right = FALSE,
labels = paste(
seq(10, ceiling(max(SCL_final, na.rm = TRUE) / 10) * 10 - 10, by = 10),
seq(19, ceiling(max(SCL_final, na.rm = TRUE) / 10) * 10 - 1, by = 10),
sep = "-")))
## FLORIDA
## Uses dplyr (filter, group_by, mutate, summarise)
summary_table_size_FL <- meta.cc %>%
filter(State == "FL") %>%
group_by(SCL_bin) %>%
summarise(
N = sum(!is.na(Growth_Rate)),
Mean = mean(Growth_Rate, na.rm = TRUE),
SD = sd(Growth_Rate, na.rm = TRUE),
Year_min = min(Year, na.rm = TRUE),
Year_max = max(Year, na.rm = TRUE)
) %>%
ungroup() %>%
mutate(across(c(N, Mean, SD), ~ signif(., 4))) %>%
filter(!is.na(SCL_bin))
## Display the created table
summary_table_size_FL## # A tibble: 6 × 6
## SCL_bin N Mean SD Year_min Year_max
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 30-39 4 1.85 2.23 2012 2015
## 2 40-49 6 3.27 1.31 2010 2019
## 3 50-59 8 3.35 2.07 2012 2022
## 4 60-69 13 1.4 1.44 2014 2021
## 5 70-79 15 1.07 1.78 2009 2022
## 6 80-89 1 5.3 NA 2022 2022## TEXAS
## Uses dplyr (filter, group_by, mutate, summarise)
summary_table_size_TX <- meta.cc %>%
filter(State == "TX") %>%
group_by(SCL_bin) %>%
summarise(
N = sum(!is.na(Growth_Rate)),
Mean = mean(Growth_Rate, na.rm = TRUE),
SD = sd(Growth_Rate, na.rm = TRUE),
Year_min = min(Year, na.rm = TRUE),
Year_max = max(Year, na.rm = TRUE)
) %>%
ungroup() %>%
mutate(across(c(N, Mean, SD), ~ signif(., 4))) %>%
filter(!is.na(SCL_bin))
## Display the created table
summary_table_size_TX## # A tibble: 5 × 6
## SCL_bin N Mean SD Year_min Year_max
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 30-39 7 1.27 0.531 2009 2015
## 2 40-49 30 2.26 1.76 2005 2017
## 3 50-59 40 2.84 2.09 2010 2021
## 4 60-69 57 1.90 1.77 2013 2021
## 5 70-79 17 1.82 2.05 2016 2021Plot mean growth rates per year & age class
Using base plot to build a plot of mean growth rate with SDs over Age, with Age group in diffierent colors
## Plot Means by Year & Age Class
## Round ages to nearest whole number, creating 'Age' as a new column with the rounded Ages
summary_table_age_sub <- summary_table_age %>%
mutate(
Age_round = round(Age_at_stranding),
Age_group = case_when(
Age_round >= 12 & Age_round <= 15 ~ "Age 12-15",
Age_round >= 16 & Age_round <= 19 ~ "Age 16-19",
Age_round >= 20 & Age_round <= 24 ~ "Age 20-24"
),
SE = SD / sqrt(N)
) %>%
filter(!is.na(Age_group))
## Assign colors
age_colors <- c("Age 12-15" = "red",
"Age 16-19" = "blue",
"Age 20-24" = "green")
## Initialize empty plot with fixed y-axis range 0-6
par(mar = c(5,5,2,1))
plot(NA, NA, xlim = c(12, 24),
ylim = c(0,6),
xlab = "Age (years)",
ylab = expression(paste("Mean Growth Rate (cm SCL/yr)")),
las = 1)
## Loop over age groups
for(group in unique(summary_table_age_sub$Age_group)) {
dat <- summary_table_age_sub %>% filter(Age_group == group)
if(nrow(dat) == 0) next # skip empty groups
I <- order(dat$Age_round)
## Shaded SD area
polygon(c(dat$Age_round[I], rev(dat$Age_round[I])),
c(dat$Mean[I]-dat$SD[I], rev(dat$Mean[I]+dat$SD[I])),
col = adjustcolor(age_colors[group], alpha.f=0.2),
border = NA)
## Mean line
lines(dat$Age_round[I], dat$Mean[I], col = age_colors[group], lwd = 2)
## Points
points(dat$Age_round[I], dat$Mean[I], pch = 16, col = age_colors[group])
}
## Add legend
legend("topright", legend = names(age_colors), col = age_colors, lwd = 2, pch = 16, bty = "n")
## TO SAVE: COPY & PASTE THIS CODE INTO THE CONSOLE BELOW AND RUN. EXPORT THE PLOT FROM THE RIGHT-HAND PANEL.## Plot Average Growth Rates for the last 17 years of life per Nutritional Status
## Assign colors to match Nutritional Status labels
status_colors <- c("Emaciated" = "red", "Good" = "green", "Robust" = "blue")
Growth.Yrsbefdeath <- ggplot(meta.cc_clean, aes(x = Year_Before_Death, y = Growth_Rate,
color = Nutritional_Status, group = Nutritional_Status)) +
stat_summary(fun = mean, geom = "line", size = 1) +
stat_summary(fun = mean, geom = "point", size = 2) +
stat_summary(fun.data = mean_cl_normal, geom = "ribbon",
aes(fill = Nutritional_Status), alpha = 0.2, color = NA) +
scale_color_manual(values = status_colors) +
scale_fill_manual(values = status_colors) +
labs(x = "Years Before Death", y = "Mean Growth Rate (cm SCL/yr)",
color = "Nutritional Status", fill = "Nutritional Status") +
theme_minimal(base_size = 14) +
theme(legend.position = "top")## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## View Plot
Growth.Yrsbefdeath## Warning: Removed 14 rows containing missing values or values outside the scale range
## (`geom_ribbon()`).
## Save plot to desired location
ggsave("Plots/Growth.Rates.over.Yrs.Before.Death.png", plot = Growth.Yrsbefdeath,
width = 7, height = 5, dpi = 300)## Warning: Removed 14 rows containing missing values or values outside the scale range
## (`geom_ribbon()`).## Plot Average Growth Rates for the last 10 years of life per Nutritional Status
Growth.Yrsbefdeath_10 <- ggplot(meta.cc_clean_10y, aes(x = Year_Before_Death, y = Growth_Rate,
color = Nutritional_Status, group = Nutritional_Status)) +
stat_summary(fun = mean, geom = "line", size = 1) +
stat_summary(fun = mean, geom = "point", size = 2) +
stat_summary(fun.data = mean_cl_normal, geom = "ribbon",
aes(fill = Nutritional_Status), alpha = 0.2, color = NA) +
scale_color_manual(values = status_colors) +
scale_fill_manual(values = status_colors) +
labs(x = "Years Before Death", y = "Mean Growth Rate (cm SCL/yr)",
color = "Nutritional Status", fill = "Nutritional Status") +
theme_minimal(base_size = 14) +
theme(legend.position = "top")
## View Plot
Growth.Yrsbefdeath_10## Warning: Removed 7 rows containing non-finite outside the scale range
## (`stat_summary()`).
## Removed 7 rows containing non-finite outside the scale range
## (`stat_summary()`).
## Removed 7 rows containing non-finite outside the scale range
## (`stat_summary()`).## Warning: Removed 10 rows containing missing values or values outside the scale range
## (`geom_ribbon()`).
## Save plot to desired location
ggsave("Plots/Growth.Rates.over.Yrs.Before.Death_10.png", plot = Growth.Yrsbefdeath_10,
width = 7, height = 5, dpi = 300)## Warning: Removed 7 rows containing non-finite outside the scale range
## (`stat_summary()`).## Warning: Removed 7 rows containing non-finite outside the scale range
## (`stat_summary()`).
## Removed 7 rows containing non-finite outside the scale range
## (`stat_summary()`).## Warning: Removed 10 rows containing missing values or values outside the scale range
## (`geom_ribbon()`).## Read in metadata as a .csv file
LIA.TX <- read.csv(file="TX.csv", header=T)
LIA.FL <- read.csv(file="FL.csv", header=T)
LIA.ALL <- read.csv(file="All_States.csv", header=T)
## Display the column headers
head(LIA.ALL)## Age N Mean SD Year_min Year_max
## 1 0.00 9 9.95 3.41 1996 2008
## 2 0.75 1 NA NA NA NA
## 3 1.75 9 3.56 2.02 1988 2010
## 4 2.75 24 2.73 2.07 1989 2011
## 5 3.75 32 2.86 1.64 1984 2011
## 6 4.75 40 3.29 2.09 1985 2012LIA.ALL <- LIA.ALL %>%
mutate(
Age_round = round(Age),
Age_group = case_when(
Age_round >= 3 & Age_round <= 7 ~ "Age 3-7",
Age_round >= 8 & Age_round <= 11 ~ "Age 8-11",
Age_round >= 12 & Age_round <= 15 ~ "Age 12-15",
Age_round >= 16 & Age_round <= 19 ~ "Age 16-19",
Age_round >= 20 & Age_round <= 24 ~ "Age 20-24"))
## Plot
# Ensure Age_round is numeric
LIA.ALL$Age_round <- as.numeric(LIA.ALL$Age_round)
meta.cc_clean$Age_round <- as.numeric(meta.cc_clean$Age_round)
growth_comparison_plot <- ggplot() +
# Bottom: box & whiskers (this study)
geom_boxplot(
data = meta.cc_clean,
aes(
x = Age_round,
y = Growth_Rate,
group = Age_round,
color = "This study"
),
fill = "grey80",
alpha = 0.6,
width = 0.6,
outlier.alpha = 0.3
) +
# Top: spline with 95% CI (Avens et al. 2025)
geom_smooth(
data = LIA.ALL,
aes(x = Age_round, y = Mean),
method = "gam",
formula = y ~ s(x, k = 6),
se = TRUE,
color = "blue",
fill = "lightblue",
linewidth = 1.2,
show.legend = FALSE # <-- hides CI from legend
) +
# Mean points (Avens et al. 2025)
geom_point(
data = LIA.ALL,
aes(
x = Age_round,
y = Mean,
color = "Avens et al. 2025"
),
size = 2
) +
# Manual legend colors
scale_color_manual(
name = "Data Source",
values = c(
"Avens et al. 2025" = "blue",
"This study" = "black"
)
) +
scale_x_continuous(
limits = c(0, 45),
breaks = seq(0, 45, 5)
) +
scale_y_continuous(
limits = c(0, 10)
) +
labs(
x = "Rounded Age (yr)",
y = "Mean Growth Rate (cm SCL / yr)",
title = "Growth rate comparison: literature vs this study"
) +
theme_minimal(base_size = 14) +
theme(
legend.position = "top"
)
## View Plot
growth_comparison_plot## Warning: Removed 2 rows containing non-finite outside the scale range
## (`stat_smooth()`).## Warning: Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).
## Save plot to desired location
ggsave("Plots/growth_comparison_plot_LIA.data.All.states.png", plot = growth_comparison_plot,
width = 7, height = 5, dpi = 300)## Warning: Removed 2 rows containing non-finite outside the scale range (`stat_smooth()`).
## Removed 2 rows containing missing values or values outside the scale range
## (`geom_point()`).## Plot
## View Plot
## Save plot to desired location
## ggsave("Plots/PLOT NAME.png", plot = PLOT_NAME,
## width = 7, height = 5, dpi = 300)## Plot
## View Plot
## Save plot to desired location
## ggsave("Plots/PLOT NAME.png", plot = PLOT_NAME,
## width = 7, height = 5, dpi = 300)