rm(list=ls())
# Paquetes
library(ggplot2)
library(readxl)
library(gridExtra)

# Leer datos
ruta_archivo <- "data-cz-values-interactions-parana.xlsx"
datos <- read_excel(ruta_archivo)

# Tipos
datos$STATUS <- as.factor(datos$STATUS)
datos$sform  <- as.factor(datos$sform)
datos$C <- as.numeric(datos$C)
datos$Z <- as.numeric(datos$Z)

# Tests
kruskal.test(C ~ STATUS, data = datos)
kruskal.test(Z ~ STATUS, data = datos)
kruskal.test(C ~ sform, data = datos)
kruskal.test(Z ~ sform, data = datos)

# Calcular límites Y comunes para estandarizar escalas entre paneles
y_rng <- range(datos$Z, na.rm = TRUE)
pad <- diff(y_rng) * 0.05
y_limits <- c(y_rng[1] - pad, y_rng[2] + pad)

# Paletas
pal_status <- c("R" = "#B0BEC5", "M" = "#90A4AE")
pal_sform <- c("ANNUAL HERB" = "#a6bddb", 
               "CLIMBING PLANT" = "#d0d1e6", 
               "NN" = "#b8cde6", 
               "PARASITE BUSH" = "#9ebcda",
               "SHRUB" = "#8c96c6",
               "SHRUB/TREE" = "#6e91c6",
               "TREE" = "#4f79b6")
# tamaño común para etiquetas de ejes
axis_label_size <- 14

# Plot 1: STATUS — violin + narrow box + points
p1 <- ggplot(datos, aes(x = STATUS, y = Z, fill = STATUS)) +
  geom_violin(trim = FALSE, width = 0.9, scale = "width") +
  geom_boxplot(width = 0.12, outlier.shape = NA, alpha = 0.9) +
  geom_jitter(width = 0.12, height = 0, size = 1.6, alpha = 0.9) +
  scale_fill_manual(values = pal_status) +
  coord_cartesian(ylim = y_limits) +
  labs(
    title = "Within-Module Degree (Z) — Status",
    x = "Species Status (R = Resident, M = Migratory)",
    y = "Within-Module Degree (Z)"
  ) +
  theme_classic(base_size = 15, base_family = "Times New Roman") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 17),
    axis.title.x = element_text(size = axis_label_size, margin = margin(t = 10)),
    axis.title.y = element_text(size = axis_label_size, margin = margin(r = 10)),
    axis.text = element_text(color = "black"),
    legend.position = "none"
  )

# Plot 2: sform — violin + narrow box + points
p2 <- ggplot(datos, aes(x = sform, y = Z, fill = sform)) +
  geom_violin(trim = FALSE, width = 0.9, scale = "width") +
  geom_boxplot(width = 0.08, outlier.shape = NA, alpha = 0.9) +
  geom_jitter(width = 0.12, height = 0, size = 1.3, alpha = 0.9) +
  scale_fill_manual(values = pal_sform) +
  coord_cartesian(ylim = y_limits) +
  labs(
    title = "Within-Module Degree (Z) — Life form",
    x = "Life form",
    y = "Within-Module Degree (Z)"
  ) +
  theme_classic(base_size = 13, base_family = "Times New Roman") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 16),
    axis.title.x = element_text(size = axis_label_size, margin = margin(t = 10)),
    axis.title.y = element_text(size = axis_label_size, margin = margin(r = 10)),
    axis.text.x = element_text(angle = 45, hjust = 1, color = "black"),
    axis.text.y = element_text(color = "black"),
    legend.position = "none"
  )


# Combinar con panel widths estandarizados
combined <- arrangeGrob(p1, p2, ncol = 2, widths = c(1, 1))

# Guardar
ggsave("figure3_violin_boxpoints.png", combined, width = 16, height = 6, dpi = 600, units = "in")


###otro
# tamaño común para etiquetas de ejes
axis_label_size <- 16

# Mediana como línea horizontal (crossbar) y SIN etiquetas de n
p1 <- ggplot(datos, aes(x = STATUS, y = Z, fill = STATUS)) +
  geom_violin(trim = FALSE, width = 0.9, scale = "width", alpha = 0.6) +
  geom_jitter(width = 0.12, height = 0, size = 1.6, alpha = 0.8) +
  stat_summary(fun = median, geom = "crossbar", width = 0.25, fatten = 0, color = "black") +
  scale_fill_manual(values = pal_status) +
  coord_cartesian(ylim = y_limits) +
  labs(title = "Within-Module Degree (Z) — Status",
       x = "Species Status (R = Resident, M = Migratory)",
       y = "Within-Module Degree (Z)") +
  theme_classic(base_size = 15, base_family = "Times New Roman") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 17),
    axis.title.x = element_text(size = axis_label_size, margin = margin(t = 10)),
    axis.title.y = element_text(size = axis_label_size, margin = margin(r = 10)),
    axis.text = element_text(color = "black"),
    legend.position = "none"
  )

p2 <- ggplot(datos, aes(x = sform, y = Z, fill = sform)) +
  geom_violin(trim = FALSE, width = 0.9, scale = "width", alpha = 0.6) +
  geom_jitter(width = 0.12, height = 0, size = 1.3, alpha = 0.8) +
  stat_summary(fun = median, geom = "crossbar", width = 0.25, fatten = 0, color = "black") +
  scale_fill_manual(values = pal_sform) +
  coord_cartesian(ylim = y_limits) +
  labs(title = "Within-Module Degree (Z) — Life form",
       x = "Life form",
       y = "Within-Module Degree (Z)") +
  theme_classic(base_size = 13, base_family = "Times New Roman") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 16),
    axis.title.x = element_text(size = axis_label_size, margin = margin(t = 10)),
    axis.title.y = element_text(size = axis_label_size, margin = margin(r = 10)),
    axis.text.x = element_text(angle = 45, hjust = 1, color = "black"),
    axis.text.y = element_text(color = "black"),
    legend.position = "none"
  )
# Combinar con panel widths estandarizados
combined <- arrangeGrob(p1, p2, ncol = 2, widths = c(1, 1))

# Guardar
ggsave("figure3_violin_boxpointssincajaotrob.png", combined, width = 16, height = 6, dpi = 600, units = "in")