Beyond bleaching: Collapse of net coral reef carbonate budgets in the Eastern Tropical Pacific after the 4th global coral bleaching event.

Description of the data and file structure Title:** Beyond bleaching: Collapse of net coral reef carbonate budgets in the Eastern Tropical Pacific after the 4th global coral bleaching event.** Francisco Medellín-Maldonado1,2*, Rebeca Granja-Fernández3, Tania González-Mendoza4, Rafael Andrés Cabral-Tena5, Lorenzo Álvarez-Filip1, Andrés López-Pérez2 Abstract Marine heatwaves' increasing frequency and severity drive unprecedented ecological changes in coral reefs. Although extensive research has documented mass coral bleaching events and subsequent declines in live coral cover worldwide, the impacts on structural functionality and the disruption of key mechanisms governing reef accretion remain poorly understood. Here, we investigate the effects of the fourth global bleaching event on carbonate production and degradation processes in coral reefs of the southern Eastern Tropical Pacific. Using a combination of field surveys, carbonate production models, and thermal stress analyses, we show that prolonged SST anomalies peaked at 2.24 °C in August 2023, with degree heating weeks reaching a record-breaking 16, far surpassing established bleaching and mortality thresholds. This extreme thermal event led to regional coral mortality exceeding 76%, with some reefs experiencing complete loss of live coral cover. Consequently, these reefs shifted from highly productive carbonate states (4.25 ± 2.13 kg·CaCO₃·m⁻²·yr⁻¹) to net erosional states (-2.67 ± 2.82 kg·CaCO₃·m⁻²·yr⁻¹). Furthermore, we observed significant disruptions in the mechanisms underlying both carbonate production and erosion. Our findings underscore the vulnerability of reef structural functionality to climate change-induced thermal stress and provide crucial insights into the potential status of coral reefs affected by 4th global coral bleaching event. Repository contains the following: R script, and raw data Figures data.txt code.txt Statics code code.txt data.txt Files DATA_RSPB-2025-0482.xlsx CODE_RSPB-2025-0482_reviwers_2.R data_cover_SST_3_c.txt cobertura_huatus_23vs24.txt SST_Huatulco_promedio_2016.txt data_grupos.txt H23vs24.txt Files and variables File: DATA_RSPB-2025-0482.xlsx Description: Data of SST and DHW from 1985 to 2023, Coral cover from 2016 to 2024 **Files: **data_cover_SST_3_c.txt to H23vs24.txt is the data necesar to run de code in R Variables Data base of surveys and historical monthly sea surface temperature (SST, DHW), records for the study area were obtained from the NOAA Coral Reef Watch Coral Bleaching Thermal Stress dataset (NOAA/CRW), which is based on satellite observations with a 5-km resolution and covers the period from January 1, 1985 to December 31, 2023 (accessed on January 12, 2024). Code/software ################ CODE RSPB-2025-0482 ############################# ##### ##################### ############################################33 # FIT MODELS# library(glmmTMB) library(ggeffects) library(DHARMa) library(ggplot2) ########################################################### data_cover_SST_3_c names(data_cover_SST_3_c) view(data_cover_SST_3_c) data_cover_SST_3_c$Year<- as.factor(data_cover_SST_3_c$Year) data_cover_SST_3_c$Reef<- as.factor(data_cover_SST_3_c$Reef) data_cover_SST_3_c$Seson<- as.factor(data_cover_SST_3_c$Seson) data_cover_SST_3_c<-within(data_cover_SST_3_c, anidado <- factor(Reef:Year))## C O N C A T E N A R data_cover_SST_3_c$anidado<- as.factor(data_cover_SST_3_c$anidado) hist(data_cover_SST_3_c$coral.cover, breaks = 20, main = "Distribución de cobertura coralina") boxplot(coral.cover ~ Seson, data = data_cover_SST_3_c) shapiro.test(data_cover_SST_3_c$coral.cover) boxplot(Max_Av._DHW ~ Reef, data = data_cover_SST_3_c) levels(data_cover_SST_3_c$Seson) data_cover_SST_3_c$Seson<- factor(data_cover_SST_3_c$Seson, levels = c ("Before Niño", "After Niño" )) levels(data_cover_SST_3_c$Reef) data_cover_SST_3_c$Reef <- factor(data_cover_SST_3_c$Reef, levels = c ("San Agistín", "Riscalillo" ,"Jicaral" , "Dos Hermanas", "La India","Cacaluta", "Maguey","Organo", "Violín" , "La Entrega", "Chahué" )) ########## #########################################33 data_cover_SST_3_c$prop_cover <- (data_cover_SST_3_c$coral.cover + 0.001) / 100.002 mod_beta <- glmmTMB(prop_cover ~ Reef + Year, family = beta_family(link = "logit"), data = data_cover_SST_3_c) summary(mod_beta) simulateResiduals(mod_beta, n=250, plot= T) # Anova(mod_beta) Reef_unit<- ggpredict(mod_beta, terms = c("Reef", "Year"), type = "fixed") %>% plot() + theme_classic(base_size = 25)+ labs(x = "", y= "% Coral cover")+ guides()+ggtitle("") + theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_color_manual(values = c("#F38D30", "#B2D8EE","#87B6D9","#5F90BB","#3B6895", "navy", "blue","red")) Reef_unit #1200*550 #============================================ # Predicciones pred_beta <- ggpredict(mod_beta, terms = c("Reef", "Year"), type = "fixed") # Pasar a data.frame y escalar a % pred_df <- as.data.frame(pred_beta) pred_df$predicted <- pred_df$predicted * 100 pred_df$conf.low <- pred_df$conf.low * 100 pred_df$conf.high <- pred_df$conf.high * 100 # Gráfico categórico con error bars Reef_unit2 <- ggplot(pred_df, aes(x = x, y = predicted, color = group)) + geom_pointrange(aes(ymin = conf.low, ymax = conf.high), position = position_dodge(width = 0.6), size = 1) + geom_point(data = data_cover_SST_3_c, aes(x = Reef, y = prop_cover * 100, color = Year), position = position_jitter(width = 0.2, height = 0), size = 2, alpha = 0.6) + theme_classic(base_size = 20) + labs(x = "", y = "Coral cover (%)", color = "year") + theme( axis.line = element_line(color = "black"), panel.border = element_rect(color = "black", fill = NA, size = 0.5), axis.text.x = element_text(angle = 40, hjust = 1), legend.position = "top", legend.direction = "horizontal", # leyenda en horizontal legend.justification = "center" # centrada ) + #guides(colour = "none")+ scale_color_manual(values = c("#F38D30", "#B2D8EE","#87B6D9","#5F90BB","#3B6895", "navy", "blue","red")) Reef_unit2 ############################################################################################################# Year_plot<- ggpredict(mod_beta, terms = c("Year"), type = "fixed") %>% plot() + theme_classic(base_size = 25)+ labs(x = "Year", y= "% Coral cover")+ guides()+ggtitle(" ") + geom_line(color = "#8B0000", size = 1) + annotate("rect",xmin=2023,xmax=2024,ymin=0,ymax=Inf, alpha=0.4, fill='#FF4500', )+ annotate("rect",xmin=2015,xmax=2016,ymin=0,ymax=Inf, alpha=0.4, fill='#FF4500', )+ scale_x_continuous(breaks=c(2015,2016,2017,2018,2019,2020,2021,2022,2023,2024))+ scale_y_continuous(breaks=c(0,0.20,0.40,0.60,0.80, 1.10))+ theme(axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5)) Year_plot #1200*550 #============================================ ############################################################################################################# view(data_cover_SST_3_c) mod_beta_DHW <- glmmTMB(prop_cover ~ Max_Es._DHW + (1|Reef), family = beta_family(link = "logit"), data = data_cover_SST_3_c) mod_beta_DHW1 <- glmmTMB(prop_cover ~ Max_Es._DHW + (1|Reef)+(1|Year), family = beta_family(link = "logit"), data = data_cover_SST_3_c) mod_beta_DHW3 <- glmmTMB(prop_cover ~ Max_Es._DHW +Reef+ (1|Year), family = beta_family(link = "logit"), data = data_cover_SST_3_c) mod_beta_DHW2 <- glmmTMB(prop_cover ~ Max_Av._DHW + (1|Reef), family = beta_family(link = "logit"), data = data_cover_SST_3_c) AIC(mod_beta_DHW, mod_beta_DHW1, mod_beta_DHW2, mod_beta_DHW3) summary(mod_beta_DHW3) Anova(mod_beta_DHW1) simulateResiduals(mod_beta_DHW, n=250, plot= T) # simulateResiduals(mod_beta_DHW1, n=250, plot= T) # the best simulateResiduals(mod_beta_DHW3, n=250, plot= T) # ggpredict(mod_beta_DHW1, terms = c("Max_Es._DHW"), type = "fixed") %>% plot(data = T, colors = "#8B0000") + theme_classic(base_size = 25)+ labs(x = "Max annual DHW", y= "% Coral cover")+ guides()+ggtitle(" ") + theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5)) ggpredict(mod_beta_DHW3, terms = c("Max_Es._DHW","Reef"), type = "fixed") %>% plot() + theme_classic(base_size = 25)+ labs(x = "Max annual DHW", y= "% Coral cover")+ guides()+ggtitle(" ") + scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"), name = "Reef")+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5)) ggpredict(mod_beta_DHW2, terms = c("Max_Av._DHW"), type = "fixed") %>% plot(data = T, colors = "#8B0000") + theme_classic(base_size = 25)+ labs(x = "Max annual DHW", y= "% Coral cover")+ guides()+ggtitle(" ") + theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5)) pred_DHW <- ggpredict(mod_beta_DHW1, terms = "Max_Es._DHW") # --- 4. Gráfico final --- ggplot() + # Datos observados (prop_cover en porcentaje) geom_point(data = data_cover_SST_3_c, aes(x = Max_Es._DHW, y = prop_cover*100, color = Reef), alpha = 0.7, size = 3) + # Curva del modelo geom_line(data = pred_DHW, aes(x = x, y = predicted*100), color = "#8B0000", size = 1.2) + # Banda de confianza geom_ribbon(data = pred_DHW, aes(x = x, ymin = conf.low*100, ymax = conf.high*100), fill = "#8B0000", alpha = 0.4) + geom_hline(yintercept =50, color= "red",linetype="dashed", lwd=1)+ ###media geom_vline(xintercept =5, color= "black",linetype="dashed", lwd=1)+ ###media scale_x_continuous(breaks=c(0,2,4,6,8,10,12,14,16,18))+ # Estilo #guides(colour = "none")+ theme_classic(base_size = 30) + theme(axis.line = element_line(color = "black"), panel.border = element_rect(color = "black", fill = NA, size = 1))+ labs(x = "Annual maximum DHW (C°-weeks)", y = "% Coral cover") + scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"), name = "Reef") ########################################################################################### ################################################################################################################################### data_cober2<-cobertura_huatus_23vs24 names(data_cober2) view(data_cober2) data_cober2$Loc<- as.factor(data_cober2$Loc) data_cober2$Año<- as.factor(data_cober2$Año) data_cober2$Mes<- as.factor(data_cober2$Mes) data_cober2$date<-as.Date(data_cober2$date) levels(data_cober2$Loc) data_cober2$Loc <- factor(data_cober2$Loc, levels = c ("San Agustín", "Riscalillo","Jicaral", "Dos Hermanas", "La India","Cacaluta", "Maguey","Organo", "Violín", "La Entrega", "Chahué" )) Cober_hut<- ggplot(data_cober2, aes(x = Año, y = Vivo, group =Loc, color=Loc)) + geom_point(size=7, pch = 21, stroke = 3) + geom_line(size=1.2, linewidth = 1.5,position = position_dodge(width = 0.01)) + scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + xlab('')+ scale_y_continuous(breaks=c(0,20,40,60,80,100))+ ylab('% Coral cover')+ theme_classic(base_size = 25) + theme( axis.line = element_line(color = "black"), # Asegura que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añade un borde cerrado ) + labs(title = "") Cober_hut ## 1150*1650 Cober_hut2<- ggplot(data_cober2, aes(x = Año, y = Vivo, group =Loc, fill=Loc, color=Loc)) + geom_point(size=7, pch = 21, stroke = 3) + geom_smooth(method = "lm")+ theme_classic(base_size = 20) + scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Arrecife") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + xlab('')+ ylab('Cobertura de coral (% diferencia)') ############################################################################################### mea_dat_SST<-SST_Huatulco_promedio_2016 names(mea_dat_SST) summary(mea_dat_SST) str(mea_dat_SST) mea_dat_SST$day <-as.factor(mea_dat_SST$day) mea_dat_SST$Month <-as.factor(mea_dat_SST$Month) mea_dat_SST$Year <-as.factor(mea_dat_SST$Year) mea_dat_SST$dates<-as.Date(mea_dat_SST$dates) mea_dat_SST<-within(mea_dat_SST, Time <- factor(day:Month))## C O N C A T E N A R mea_dat_SST$Month<-as.ts (mea_dat_SST$Month) SST_huatus_mean <- ggplot(mea_dat_SST, aes(x = Month, y = SST, colour=Year)) + #facet_wrap(~year, nrow=2) + # a panel for each mountain range geom_point(alpha = 0.3) + theme_classic(base_size = 25) + theme( axis.line = element_line(color = "black"), # Asegura que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añade un borde cerrado ) + geom_smooth(method = "gam", aes(color = Year))+ scale_x_continuous(breaks=c(0,2,4,6,8,10,12))+ scale_y_continuous(breaks=c(22,24,26,28,30,32,34))+ xlab("Month of each year") + ylab("SST")+ labs(title = " ")+ annotate("text", x = 6.5, y = 31.7, label = "max SST 2023", size = 5, color = "#22292F") + scale_color_manual(values = c("#1E90FF","#1E90FF","#1C86EE","#87CEFA", "#104E8B","#00BFFF","#00BFFF","red","#009ACD","#00688B","#87CEFA","#B0E2FF","#A4D3EE", "#CCD6EB","#D7E8FF","#CEE2FF","#C6DDFD","#BDD7FB","#AACAF5", "#F4F4FD","#F4F3FF","#E9EBF5","#F0EEFF","#EBE9FF","#E6E4FE","#E1DEFC", "#DBD7FB","#FA4E5A","#F53F4D","#F02B3F","#E42237","#D61F33","#C81B2E","#BA182A","#AC1425","#9E1021","#910B1C","#760312","red"))+ scale_fill_manual(values = c("#1E90FF","#1E90FF","#1C86EE","#87CEFA", "#104E8B","#00BFFF","#00BFFF","red","#009ACD","#00688B","#87CEFA","#B0E2FF","#A4D3EE", "#CCD6EB","#D7E8FF","#CEE2FF","#C6DDFD","#BDD7FB","#AACAF5", "#F4F4FD","#F4F3FF","#E9EBF5","#F0EEFF","#EBE9FF","#E6E4FE","#E1DEFC", "#DBD7FB","#FA4E5A","#F53F4D","#F02B3F","#E42237","#D61F33","#C81B2E","#BA182A","#AC1425","#9E1021","#910B1C","#760312","red"))+ geom_hline(yintercept =31.2 , color= "black",linetype="dashed", lwd=1.2)+ ### guides(colour = "none") SST_huatus_mean ### 550*1000 DHW_huatus_mean <- ggplot(mea_dat_SST, aes(x = Month, y = DHW, colour=Year)) + #facet_wrap(~year, nrow=2) + geom_line(alpha = 2.5) + theme_classic(base_size = 25) + theme( axis.line = element_line(color = "black"), # Asegura que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añade un borde cerrado ) + #geom_smooth(method = "gam")+ scale_x_continuous(breaks=c(0,2,4,6,8,10,12))+ scale_y_continuous(breaks=c(0,2,4,6,8,10,12,14,16))+ xlab("Month of each year") + ylab("DHW")+ labs(title = "")+ guides(colour = "none")+ guides(fill = "none")+ geom_hline(yintercept =16.2 , color= "black",linetype="dashed", lwd=1.2)+ #geom_hline(yintercept =4.9 , color= "#760312",linetype="dashed", lwd=1.2)+ #annotate("text", x = 3, y = 5.5, label = "avarege DHW 2023", size = 4, color = "#22292F") + annotate("text", x = 8, y = 17.2, label = "max DHW 2023", size = 5, color = "#22292F") + scale_color_manual(values = c("#1E90FF","#1E90FF","#1C86EE","#87CEFA", "#104E8B","#00BFFF","#00BFFF","red","#009ACD","#00688B","#87CEFA","#B0E2FF","#A4D3EE", "#CCD6EB","#D7E8FF","#CEE2FF","#C6DDFD","#BDD7FB","#AACAF5", "#F4F4FD","#F4F3FF","#E9EBF5","#F0EEFF","#EBE9FF","#E6E4FE","#E1DEFC", "#DBD7FB","#FA4E5A","#F53F4D","#F02B3F","#E42237","#D61F33","#C81B2E","#BA182A","#AC1425","#9E1021","#910B1C","#760312","red"))+ scale_fill_manual(values = c("#1E90FF","#1E90FF","#1C86EE","#87CEFA", "#104E8B","#00BFFF","#00BFFF","red","#009ACD","#00688B","#87CEFA","#B0E2FF","#A4D3EE", "#CCD6EB","#D7E8FF","#CEE2FF","#C6DDFD","#BDD7FB","#AACAF5", "#F4F4FD","#F4F3FF","#E9EBF5","#F0EEFF","#EBE9FF","#E6E4FE","#E1DEFC", "#DBD7FB","#FA4E5A","#F53F4D","#F02B3F","#E42237","#D61F33","#C81B2E","#BA182A","#AC1425","#9E1021","#910B1C","#760312","red")) DHW_huatus_mean ### 550*1000 panel1bc<-SST_huatus_mean/DHW_huatus_mean panel1bc<-panel1bc+plot_layout(guides = "collect") panel1bc str() ###################3333 names (data_grupos) data_grupos <- data_grupos %>% mutate(Year = as.factor(Year)) Hut_grupos<-ggplot(data_grupos, aes(width = 0.9, fill=Group, y=as.numeric(Cover), x=Year ))+ geom_bar(position="stack", stat="identity" )+ #scale_x_continuous(breaks = c(2023, 2024)) + labs(x = "", y = "% Cover groups") + labs(title = "")+ scale_fill_manual(values = c("#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7")) + theme_classic(base_size = 22)+theme( axis.line = element_line(color = "black"), # Asegura que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añade un borde cerrado ) zplor<-Hut_grupos+coord_polar(theta = "x",start=0) + theme_bw() + theme(legend.key.size = unit(1.8, "cm")) + ylim(c(-100,100))+ theme(legend.position = c(0.1, 0.1),axis.text.y = element_text(size = 30,colour="black"),axis.text.x=element_text(size = 22,colour="black"))+ geom_hline(yintercept =0 , color= "black",linetype="dashed", lwd=1)+ geom_hline(yintercept =100 , color= "black",linetype="dashed", lwd=1)+ theme(axis.text.y=element_blank(), axis.ticks.y=element_blank())+ theme(plot.title=element_text(size=22), axis.title=element_text(size=22,face="bold"), legend.text = element_text(size = 18, face="bold")) Cober_hut<- ggplot(data_cober2, aes(x = Año, y = Vivo, group =Loc, color=Loc)) + geom_point(size=7, pch = 21, stroke = 3) + geom_line(size=1.2, linewidth = 1.5,position = position_dodge(width = 0.01)) + theme_classic(base_size = 22) + scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Arrecife") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + xlab('')+ ylab('% Coral cover')+ labs(title = "b)") Hut_grupos/Cober_hut ## 1150*1650 ############################################################################################### ######################### datH<-H23vs24 names(datH) View(datH) datH$Reef<- as.factor(datH$Reef) datH$Year<- as.factor(datH$Year) datH$Status<- as.factor(datH$Status) datH<-within(datH, Status_Year <- factor(Status:Year)) datH<-within(datH, Status_Reef <- factor(Status:Reef)) levels(datH$Reef) datH$Reef <- factor(datH$Reef, levels = c ("Sn Agustin", "Riscalillo","Jicaral", "Dos Hermanas", "La india " ,"Cacaluta", "Maguey","Organo", "Violin" , "La Entrega", "Chahué" )) levels(datH$Status) datH$Status <- factor(datH$Status, levels = c ("Before ENSO", "After ENSO")) levels(datH$Status_Year) datH$Status_Reef <- factor(datH$Status_Reef, levels = c ( "Before ENSO:2023", "After ENSO:2024" )) ######################### Net G ####################### mod.HG<- lmer(Net.G~Year+(1|Reef), data=datH, REML = FALSE) Anova(mod.HG) lsmeans(mod.HG,pairwise ~ Year) Gnet<-ggpredict(mod.HG, terms = c("Year" ), type = "fixed") %>% plot() + theme_classic(base_size = 25)+ xlab("") + labs(y=expression(GkgCaCO[3]~m{-2}~year{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media annotate("rect",xmin=2023,xmax=2024,ymin=0,ymax=-Inf, alpha=0.3, fill="#FF4500", )+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) Gnet ############################################################################ ######################### CORAL ####################### mod.HGC<- lmer(Coral~Status+(1|Reef), data=datH, REML = FALSE) Anova(mod.HG) lsmeans(mod.HG,pairwise ~ Status) simulateResiduals(mod.HGC, n=250, plot= T) # Coral_NCC<-ggpredict(mod.HGC, terms = c("Status" ), type = "fixed") %>% plot() + theme_classic(base_size = 25)+ xlab("") + labs(y=expression(CoralkgCaCO[3]~m{-2}~year{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_y_continuous(breaks=c(0,5,10))+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) Coral_NCC # Predicciones gp <- ggpredict(mod.HGC, terms = "Status") Coral_NCC <- ggplot(gp, aes(x = x, y = predicted)) + # 1. Intervalos de confianza como líneas verticales sin "caps" geom_linerange(aes(ymin = conf.low, ymax = conf.high), color = "#99EDFF", linewidth = 1) + # 2. Conectar las predicciones con línea geom_line(aes(group = 1), color = "#99EDFF", linewidth = 1) + # 3. Puntos del modelo geom_point(color = "#99EDFF", size = 3, alpha = 0.5) + # 4. Agregar los datos crudos (add.data = TRUE) geom_jitter(data = datH, aes(x = Status, y = Coral), width = 0.1, alpha = 0.3, color = "gray30") + # Línea horizontal = 0 geom_hline(yintercept = 0, color = "red", linetype = "dashed", linewidth = 1) + theme_classic(base_size = 25) + #scale_y_continuous(breaks = c(0,2,6,8,10)) + labs(x = "", y = expression(Coral~kg~CaCO[3]~m^{-2}~year^{-1})) + theme(axis.line = element_line(color = "black"), panel.border = element_rect(color = "black", fill = NA, size = 0.5)) Coral_NCC ############################################################################ ######################### CCA ####################### str(datH) mod.HGCCA<- lmer(CCA~Status+(1|Reef), data=datH, REML = FALSE) simulateResiduals(mod.HGCCA, n=250, plot= T) # Anova(mod.HGCCA) lsmeans(mod.HGCCA,pairwise ~ Status) mod_beta <- glmmTMB(CCA ~ Status + (1|Reef), family = beta_family(link = "logit"), data = datH) summary(mod_beta) simulateResiduals(mod_beta, n=250, plot= T) # CCA_NCC<-ggpredict(mod_beta, terms = c("Status"), type = "fixed") %>% plot() + theme_classic(base_size = 25)+ xlab("") + #scale_y_continuous(breaks=c(0,5,10))+ labs(y=expression(CCAandotherscalcifierskg~CaCO[3]~m{-2}~year{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media annotate("rect",xmin=0.85,xmax=Inf,ymin=0,ymax=-Inf, alpha=0.3, fill='#7a5d7e', )+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) CCA_NCC #============================================ ############################################################################################################# ##### ######## ### library(ggeffects) library(ggplot2) library(dplyr) # Predicciones gp <- ggpredict(mod_beta, terms = "Status") CCA_NCC <- ggplot(gp, aes(x = x, y = predicted)) + # 1. Intervalos de confianza como líneas verticales sin "caps" geom_linerange(aes(ymin = conf.low, ymax = conf.high), color = "#FFFF00", linewidth = 1) + # 2. Conectar las predicciones con línea geom_line(aes(group = 1), color = "#FFFF00", linewidth = 1) + # 3. Puntos del modelo geom_point(color = "#FFFF00", size = 3, alpha = 0.5) + # 4. Agregar los datos crudos (add.data = TRUE) geom_jitter(data = datH, aes(x = Status, y = CCA), width = 0.1, alpha = 0.3, color = "gray30") + # Línea horizontal = 0 geom_hline(yintercept = 0, color = "red", linetype = "dashed", linewidth = 1) + # Sombreado annotate("rect", xmin = 0.85, xmax = Inf, ymin = 0, ymax = -Inf, alpha = 0.3, fill = "#7a5d7e") + theme_classic(base_size = 25) + #scale_y_continuous(breaks = c(0,2,6,8,10)) + labs( x = "", y = expression(CCA~and~others~calcifiers~kg~CaCO[3]~m^{-2}~year^{-1}) ) + theme( axis.line = element_line(color = "black"), panel.border = element_rect(color = "black", fill = NA, size = 0.5) ) CCA_NCC ############################################################################ ######################### Micro ####################### mod.HGCMicro<- lmer(Micro~Status+(1|Reef), data=datH, REML = FALSE) Anova(mod.HGCMicro) lsmeans(mod.HG,pairwise ~ Status) Micro_NCC<-ggpredict(mod.HGCMicro, terms = c("Status" ), type = "fixed") %>% plot()+ xlab("") + labs(y=expression(CoralkgCaCO[3]~m{-2}~año{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media annotate("rect",xmin=0.85,xmax=Inf,ymin=0,ymax=-Inf, alpha=0.3, fill='#7a5d7e', )+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) Micro_NCC ############################################################################ ######################### Macro ####################### mod.HGCMacro<- lmer(Macro~Status+(1|Reef), data=datH, REML = FALSE) Anova(mod.HGCMacro) Macro_NCC<-ggpredict(mod.HGCMacro, terms = c("Status" ), type = "fixed") %>% plot()+ xlab("") + labs(y=expression(CoralkgCaCO[3]~m{-2}~año{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media annotate("rect",xmin=0.85,xmax=Inf,ymin=0,ymax=-Inf, alpha=0.3, fill='#7a5d7e', )+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) Macro_NCC ############################################################################ ######################### Parrotfish ####################### names(datH) View(datH) mod.HGCParrotfish<- lmer(Parrotfish~Status+(1|Reef), data=datH) summary(mod.HGCParrotfish2) Anova(mod.HGCParrotfish2) Parrotfish_NCC<-ggpredict(mod.HGCParrotfish, terms = c("Status" ), type = "fe") %>% plot(connect.lines = T, add.data = T, line.size=1, dot.size=3, dot.alpha = 0.5, colors = '#993404') + theme_classic(base_size = 25)+ xlab("") + labs(y=expression(CoralkgCaCO[3]~m{-2}~año{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media annotate("rect",xmin=0.85,xmax=Inf,ymin=0,ymax=-Inf, alpha=0.3, fill='#7a5d7e', )+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) Parrotfish_NCC ############################################################################ ######################### Urchin ####################### mod.HGCUrchin<- lmer(Urchin~Status+(1|Reef), data=datH) Anova(mod.HGCUrchin) Urchin_NCC<-ggpredict(mod.HGCUrchin, terms = c("Status" ), type = "fixed") %>% plot()+ xlab("") + labs(y=expression(CoralkgCaCO[3]~m{-2}~año{-1})) + guides()+ggtitle(" ") + geom_hline(yintercept =0, color= "red",linetype="dashed", lwd=1)+ ###media annotate("rect",xmin=0.85,xmax=Inf,ymin=0,ymax=-Inf, alpha=0.3, fill='#7a5d7e', )+ theme( axis.line = element_line(color = "black"), # Asegurar que los ejes se vean panel.border = element_rect(color = "black", fill = NA, size = 0.5) # Añadir borde alrededor )+ scale_colour_manual(values = c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF"),name = "Reef") + scale_fill_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) + scale_linetype_manual(values =c("#8B2500", "#CD6839", "#FF0000","#FF9932", "#FFCC65", "#FFEE99", "#54FF9F", "#CCFFFF", "#99EDFF","#79CDCD", "#528B8B","#0054FF")) Urchin_NCC #############################################################