calculate_summary_stats <- function(data,
metric = "inverse_simpson",
group_vars = c("site", "year")) {
summary_data <- data %>%
group_by(across(all_of(group_vars))) %>%
summarise(
mean_coralcover = mean(coralcover, na.rm = TRUE),
se_coralcover = sd(coralcover, na.rm = TRUE) / sqrt(n()),
mean_metric = mean(get(metric), na.rm = TRUE),
se_metric = sd(get(metric), na.rm = TRUE) / sqrt(n()),
.groups = 'drop'
)
return(summary_data)
}
plot_mean_and_se <- function(summary_data,
metric_label = "Inverse Simpson Index",
facet_var = "site",
smooth = FALSE) {
# Calculate scale factor to align the two variables
max_coralcover <- max(summary_data$mean_coralcover, na.rm = TRUE)
max_metric <- max(summary_data$mean_metric, na.rm = TRUE)
scale_factor <- max_coralcover / max_metric
p <- ggplot(summary_data, aes(x = year))
if (smooth) {
# Smoothed lines with confidence intervals
p <- p +
geom_smooth(
aes(
y = mean_coralcover,
color = "Coral Cover",
fill = "Coral Cover"
),
method = "loess",
span = 0.5,
se = TRUE,
alpha = 0.2
) +
geom_smooth(
aes(
y = mean_metric * scale_factor,
color = metric_label,
fill = metric_label
),
method = "loess",
span = 0.5,
se = TRUE,
alpha = 0.2
)
} else {
# Mean lines with error ribbons
p <- p +
geom_line(aes(y = mean_coralcover, color = "Coral Cover"), size = 1) +
geom_ribbon(
aes(
ymin = mean_coralcover - se_coralcover,
ymax = mean_coralcover + se_coralcover,
fill = "Coral Cover"
),
alpha = 0.3
) +
geom_line(aes(y = mean_metric * scale_factor, color = metric_label),
size = 1) +
geom_ribbon(aes(
ymin = (mean_metric - se_metric) * scale_factor,
ymax = (mean_metric + se_metric) * scale_factor,
fill = metric_label
),
alpha = 0.3)
}
p <- p +
scale_y_continuous(name = "Mean Coral Cover (%)",
sec.axis = sec_axis( ~ . / scale_factor, name = paste("Mean", metric_label))) +
scale_color_manual(name = "",
values = c("Coral Cover" = "blue", metric_label = "red")) +
scale_fill_manual(name = "",
values = c("Coral Cover" = "blue", metric_label = "red")) +
facet_wrap(as.formula(paste("~", facet_var))) +
labs(title = paste("Mean Coral Cover and", metric_label, "Over Time"),
x = "Year") +
theme_minimal() +
theme(legend.position = "bottom")
return(p)
}
