knitr::opts_chunk$set(echo = TRUE, message = FALSE, warning = FALSE,
                       fig.align = "center", fig.width = 7, fig.height = 4.5)

# Paquetes necesarios
library(ggplot2)   # Gráficos
library(dplyr)     # Manejo de datos
library(car)       # Prueba de Levene
library(knitr)      # Tablas
library(multcompView) # Letras de significancia (compact letter display)

1 Introducción

Con el fin de evaluar el efecto del uso del suelo sobre la biodiversidad de anfibios en una reserva forestal, se establecieron 52 parcelas de muestreo distribuidas en cuatro tipos de hábitat que representan un gradiente de intervención antrópica: bosque primario, bosque secundario, sistema silvopastoril y potrero (13 parcelas por hábitat). En cada parcela se registró:

  • Riqueza: número de especies de anfibios observadas.
  • Shannon: índice de diversidad de Shannon–Wiener.
  • Altitud: altitud del sitio de muestreo (m.s.n.m.).

El objetivo de este informe es:

  1. Describir univariadamente cada variable.

  2. Explorar relaciones bivariadas entre el hábitat y la riqueza/diversidad, y entre la altitud y la diversidad.

  3. Evaluar mediante ANOVA de una vía si el índice de Shannon difiere entre tipos de hábitat, verificando los supuestos del modelo y realizando las comparaciones post-hoc pertinentes.

2 Carga y descripción general de la base de datos

load("C:/Users/María del mar/Downloads/Biodiversidad (1).RData")  # Objeto: BD_biodiversidad

# Orden lógico del gradiente de intervención antrópica (de menor a mayor)
BD_biodiversidad$Habitat <- factor(
  BD_biodiversidad$Habitat,
  levels = c("Bosque primario", "Bosque secundario",
             "Sistema silvopastoril", "Potrero")
)

str(BD_biodiversidad)
## 'data.frame':    52 obs. of  5 variables:
##  $ Parcela: chr  "P001" "P002" "P003" "P004" ...
##  $ Habitat: Factor w/ 4 levels "Bosque primario",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Riqueza: num  14 14 13 14 10 15 19 16 15 21 ...
##  $ Shannon: num  2.13 2.98 2.07 2.53 2.82 2.75 2.45 1.73 2.04 2.48 ...
##  $ Altitud: num  1478 1477 1505 1395 1494 ...
##  - attr(*, "datalabel")= chr "BD_biodiversidad"
##  - attr(*, "var.labels")= chr [1:5] "" "" "" "" ...
kable(head(BD_biodiversidad, 10), caption = "Primeros registros de la base de datos")
Primeros registros de la base de datos
Parcela Habitat Riqueza Shannon Altitud
P001 Bosque primario 14 2.13 1478.0
P002 Bosque primario 14 2.98 1476.8
P003 Bosque primario 13 2.07 1505.2
P004 Bosque primario 14 2.53 1395.2
P005 Bosque primario 10 2.82 1493.5
P006 Bosque primario 15 2.75 1337.9
P007 Bosque primario 19 2.45 1482.4
P008 Bosque primario 16 1.73 1600.8
P009 Bosque primario 15 2.04 1417.5
P010 Bosque primario 21 2.48 1321.8

La base contiene 52 parcelas, 13 por cada uno de los 4 tipos de hábitat, sin valores faltantes (0 NA en total).

3 a. Análisis exploratorio univariado

3.1 Medidas resumen

resumen_univariado <- BD_biodiversidad %>%
  summarise(across(c(Riqueza, Shannon, Altitud),
                    list(n = ~sum(!is.na(.)),
                         Media = ~mean(.),
                         Mediana = ~median(.),
                         DE = ~sd(.),
                         CV_pct = ~sd(.)/mean(.)*100,
                         Min = ~min(.),
                         Q1 = ~quantile(., 0.25),
                         Q3 = ~quantile(., 0.75),
                         Max = ~max(.),
                         Asimetria = ~ (sum((.-mean(.))^3)/n())/(sd(.)^3)),
                    .names = "{.col}__{.fn}"))

tabla_larga <- data.frame(
  Variable = c("Riqueza (N° especies)", "Shannon", "Altitud (m.s.n.m.)"),
  n       = c(resumen_univariado$Riqueza__n, resumen_univariado$Shannon__n, resumen_univariado$Altitud__n),
  Media   = round(c(resumen_univariado$Riqueza__Media, resumen_univariado$Shannon__Media, resumen_univariado$Altitud__Media), 2),
  Mediana = round(c(resumen_univariado$Riqueza__Mediana, resumen_univariado$Shannon__Mediana, resumen_univariado$Altitud__Mediana), 2),
  DE      = round(c(resumen_univariado$Riqueza__DE, resumen_univariado$Shannon__DE, resumen_univariado$Altitud__DE), 2),
  CV_pct  = round(c(resumen_univariado$Riqueza__CV_pct, resumen_univariado$Shannon__CV_pct, resumen_univariado$Altitud__CV_pct), 1),
  Min     = round(c(resumen_univariado$Riqueza__Min, resumen_univariado$Shannon__Min, resumen_univariado$Altitud__Min), 2),
  Q1      = round(c(resumen_univariado$Riqueza__Q1, resumen_univariado$Shannon__Q1, resumen_univariado$Altitud__Q1), 2),
  Q3      = round(c(resumen_univariado$Riqueza__Q3, resumen_univariado$Shannon__Q3, resumen_univariado$Altitud__Q3), 2),
  Max     = round(c(resumen_univariado$Riqueza__Max, resumen_univariado$Shannon__Max, resumen_univariado$Altitud__Max), 2),
  Asimetria = round(c(resumen_univariado$Riqueza__Asimetria, resumen_univariado$Shannon__Asimetria, resumen_univariado$Altitud__Asimetria), 2)
)

kable(tabla_larga, caption = "Medidas de resumen para riqueza, Shannon y altitud (n = 52 parcelas)")
Medidas de resumen para riqueza, Shannon y altitud (n = 52 parcelas)
Variable n Media Mediana DE CV_pct Min Q1 Q3 Max Asimetria
Riqueza (N° especies) 52 9.73 9.50 4.37 44.9 2.00 6.00 13.00 21.00 0.28
Shannon 52 1.76 1.85 0.61 34.4 0.67 1.26 2.11 2.98 0.04
Altitud (m.s.n.m.) 52 1243.89 1254.80 172.71 13.9 950.90 1081.88 1382.88 1600.80 0.09

3.2 Riqueza de especies

ggplot(BD_biodiversidad, aes(x = Riqueza)) +
  geom_histogram(aes(y = after_stat(density)), binwidth = 2,
                 fill = "#2E8B57", color = "white", alpha = 0.85) +
  geom_density(color = "#1B4D3E", linewidth = 1) +
  geom_vline(aes(xintercept = mean(Riqueza)), color = "firebrick",
             linetype = "dashed", linewidth = 0.8) +
  labs(title = "Distribución de la riqueza de especies de anfibios",
       x = "Riqueza (N° de especies)", y = "Densidad") +
  theme_minimal(base_size = 12)

ggplot(BD_biodiversidad, aes(y = Riqueza, x = "")) +
  geom_boxplot(fill = "#2E8B57", alpha = 0.6, width = 0.3) +
  geom_jitter(width = 0.05, alpha = 0.5) +
  labs(title = "Boxplot general de la riqueza de especies", x = "", y = "Riqueza") +
  theme_minimal(base_size = 12)

Interpretación. La riqueza promedio es de 9.7 especies por parcela (mediana = 9.5), con una desviación estándar de 4.4 especies y un coeficiente de variación de 44.9%, lo que indica una dispersión relativa alta entre parcelas. El rango observado (2 a 21 especies) y la cercanía entre media y mediana sugieren una distribución moderadamente simétrica, aunque el histograma muestra cierta dispersión compatible con la mezcla de hábitats con distinto grado de intervención (lo que se profundiza en el análisis bivariado).

3.3 Índice de diversidad de Shannon–Wiener

ggplot(BD_biodiversidad, aes(x = Shannon)) +
  geom_histogram(aes(y = after_stat(density)), binwidth = 0.25,
                 fill = "#4682B4", color = "white", alpha = 0.85) +
  geom_density(color = "#1F3A5F", linewidth = 1) +
  geom_vline(aes(xintercept = mean(Shannon)), color = "firebrick",
             linetype = "dashed", linewidth = 0.8) +
  labs(title = "Distribución del índice de diversidad de Shannon-Wiener",
       x = "Índice de Shannon (H')", y = "Densidad") +
  theme_minimal(base_size = 12)

ggplot(BD_biodiversidad, aes(y = Shannon, x = "")) +
  geom_boxplot(fill = "#4682B4", alpha = 0.6, width = 0.3) +
  geom_jitter(width = 0.05, alpha = 0.5) +
  labs(title = "Boxplot general del índice de Shannon", x = "", y = "Shannon (H')") +
  theme_minimal(base_size = 12)

Interpretación. El índice de Shannon presenta una media de 1.76 (mediana = 1.85) con una desviación estándar de 0.61 y un CV de 34.4%. Valores de H’ entre 0.67 y 2.98 reflejan una diversidad que va de baja/moderada a relativamente alta según la parcela, lo cual es consistente con la existencia de un gradiente de intervención antrópica: se espera que los hábitats menos intervenidos concentren los valores más altos.

3.4 Altitud

ggplot(BD_biodiversidad, aes(x = Altitud)) +
  geom_histogram(aes(y = after_stat(density)), binwidth = 100,
                 fill = "#B8860B", color = "white", alpha = 0.85) +
  geom_density(color = "#7A5C00", linewidth = 1) +
  geom_vline(aes(xintercept = mean(Altitud)), color = "firebrick",
             linetype = "dashed", linewidth = 0.8) +
  labs(title = "Distribución de la altitud de las parcelas",
       x = "Altitud (m.s.n.m.)", y = "Densidad") +
  theme_minimal(base_size = 12)

Interpretación. Las parcelas se ubican en promedio a 1244 m.s.n.m. (mediana = 1255), con un rango entre 951 y 1601 m.s.n.m. y una desviación estándar de 173 m, lo que evidencia que el muestreo cubrió un gradiente altitudinal amplio dentro de la reserva, condición necesaria para poder explorar su posible relación con la diversidad de anfibios.

4 b. Análisis exploratorio bivariado

4.1 Riqueza de especies por tipo de hábitat

resumen_riqueza <- BD_biodiversidad %>%
  group_by(Habitat) %>%
  summarise(n = n(),
            Media = round(mean(Riqueza), 2),
            Mediana = median(Riqueza),
            DE = round(sd(Riqueza), 2),
            Min = min(Riqueza),
            Max = max(Riqueza))

kable(resumen_riqueza, caption = "Riqueza de especies por tipo de hábitat")
Riqueza de especies por tipo de hábitat
Habitat n Media Mediana DE Min Max
Bosque primario 13 14.77 14 2.74 10 21
Bosque secundario 13 11.08 12 3.01 3 15
Sistema silvopastoril 13 8.08 7 2.50 5 13
Potrero 13 5.00 5 1.29 2 7 
ggplot(BD_biodiversidad, aes(x = Habitat, y = Riqueza, fill = Habitat)) +
  geom_boxplot(alpha = 0.7, outlier.shape = 21) +
  geom_jitter(width = 0.12, alpha = 0.5, size = 1.8) +
  scale_fill_brewer(palette = "Greens") +
  labs(title = "Riqueza de especies de anfibios según tipo de hábitat",
       x = "Tipo de hábitat", y = "Riqueza (N° de especies)") +
  theme_minimal(base_size = 12) +
  theme(legend.position = "none",
        axis.text.x = element_text(angle = 15, hjust = 1))

Interpretación. La tabla y el boxplot muestran un patrón decreciente de riqueza a medida que aumenta la intervención antrópica: el bosque primario presenta la mayor riqueza promedio (14.77 especies), seguido del bosque secundario y el sistema silvopastoril, mientras que el potrero muestra la riqueza más baja (5 especies en promedio). Esto es coherente con la hipótesis de que la pérdida de cobertura vegetal y la simplificación estructural del hábitat reducen los microhábitats disponibles para los anfibios.

4.2 Índice de Shannon por tipo de hábitat

resumen_shannon <- BD_biodiversidad %>%
  group_by(Habitat) %>%
  summarise(n = n(),
            Media = round(mean(Shannon), 2),
            Mediana = round(median(Shannon), 2),
            DE = round(sd(Shannon), 2),
            EE = round(sd(Shannon)/sqrt(n()), 3),
            Min = round(min(Shannon), 2),
            Max = round(max(Shannon), 2))

kable(resumen_shannon, caption = "Índice de Shannon por tipo de hábitat (EE = error estándar)")
Índice de Shannon por tipo de hábitat (EE = error estándar)
Habitat n Media Mediana DE EE Min Max
Bosque primario 13 2.47 2.53 0.38 0.104 1.73 2.98
Bosque secundario 13 1.97 2.03 0.22 0.061 1.44 2.27
Sistema silvopastoril 13 1.60 1.62 0.28 0.078 1.14 2.06
Potrero 13 1.00 0.98 0.24 0.066 0.67 1.34 
ggplot(BD_biodiversidad, aes(x = Habitat, y = Shannon, fill = Habitat)) +
  geom_boxplot(alpha = 0.7, outlier.shape = 21) +
  geom_jitter(width = 0.12, alpha = 0.5, size = 1.8) +
  scale_fill_brewer(palette = "Blues") +
  labs(title = "Índice de diversidad de Shannon según tipo de hábitat",
       x = "Tipo de hábitat", y = "Índice de Shannon (H')") +
  theme_minimal(base_size = 12) +
  theme(legend.position = "none",
        axis.text.x = element_text(angle = 15, hjust = 1))

Interpretación. El patrón observado para Shannon es análogo al de riqueza: el bosque primario exhibe la diversidad promedio más alta (2.47), disminuyendo progresivamente en bosque secundario y sistema silvopastoril, hasta alcanzar el valor más bajo en potrero (1). La menor dispersión (DE) observada en potrero sugiere además una comunidad de anfibios más homogénea y dominada por pocas especies tolerantes a la intervención, mientras que el bosque primario combina mayor diversidad con mayor variabilidad entre parcelas. Esta diferencia visual se evaluará formalmente mediante ANOVA en la sección (c).

4.3 Relación entre altitud y diversidad de Shannon

cor_alt_shannon <- cor.test(BD_biodiversidad$Altitud, BD_biodiversidad$Shannon,
                             method = "pearson")

ggplot(BD_biodiversidad, aes(x = Altitud, y = Shannon)) +
  geom_point(aes(color = Habitat), size = 2.2, alpha = 0.8) +
  geom_smooth(method = "lm", se = TRUE, color = "black", linewidth = 0.8) +
  scale_color_brewer(palette = "Set2") +
  labs(title = "Relación entre altitud y diversidad de Shannon",
       x = "Altitud (m.s.n.m.)", y = "Índice de Shannon (H')",
       color = "Hábitat") +
  theme_minimal(base_size = 12)

cor_alt_shannon
## 
##  Pearson's product-moment correlation
## 
## data:  BD_biodiversidad$Altitud and BD_biodiversidad$Shannon
## t = 8.0722, df = 50, p-value = 1.287e-10
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.6031181 0.8505182
## sample estimates:
##       cor 
## 0.7522111

Interpretación. El coeficiente de correlación de Pearson entre altitud y Shannon es r = 0.75 (IC 95%: [0.6, 0.85], p = 1.29e-10), lo que indica una asociación lineal positiva y estadísticamente significativa entre ambas variables. En el diagrama de dispersión coloreado por hábitat se observa que la nube de puntos se organiza más claramente por tipo de hábitat que por altitud, lo cual sugiere que, dentro del rango altitudinal muestreado, el uso del suelo es un factor más determinante de la diversidad de anfibios que la altitud en sí misma.

5 c. ANOVA de una vía para el índice de Shannon entre hábitats

5.1 Planteamiento del modelo

Se plantea un ANOVA de una vía para evaluar si el índice de diversidad de Shannon difiere, en promedio, entre los cuatro tipos de hábitat:

Hipótesis nula (H₀): el tipo de hábitat (bosque primario, bosque secundario, sistema silvopastoril y potrero) no tiene efecto sobre la diversidad de anfibios; es decir, el promedio del índice de Shannon es igual entre los cuatro hábitats.

H0​:μBP​=μBS​=μSS​=μPO​

Hipótesis alternativa (H₁): el tipo de hábitat sí tiene efecto sobre la diversidad de anfibios; es decir, al menos uno de los cuatro hábitats presenta un promedio de Shannon diferente a los demás, siguiendo un patrón decreciente asociado al gradiente de intervención antrópica (bosque primario

H1​:al menos una media difiere

modelo_shannon <- aov(Shannon ~ Habitat, data = BD_biodiversidad)
summary(modelo_shannon)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## Habitat      3 14.862   4.954   60.61 2.38e-16 ***
## Residuals   48  3.923   0.082                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

5.2 Verificación de supuestos

5.2.1 Normalidad de los residuales (Shapiro-Wilk)

residuales <- residuals(modelo_shannon)

ggplot(data.frame(residuales), aes(sample = residuales)) +
  stat_qq(color = "#2E8B57") +
  stat_qq_line(color = "firebrick") +
  labs(title = "Gráfico Q-Q de los residuales del modelo",
       x = "Cuantiles teóricos", y = "Residuales") +
  theme_minimal(base_size = 12)

ggplot(data.frame(residuales), aes(x = residuales)) +
  geom_histogram(aes(y = after_stat(density)), bins = 10,
                 fill = "#4682B4", color = "white", alpha = 0.85) +
  geom_density(color = "#1F3A5F", linewidth = 1) +
  labs(title = "Distribución de los residuales del modelo",
       x = "Residuales", y = "Densidad") +
  theme_minimal(base_size = 12)

shapiro_res <- shapiro.test(residuales)
shapiro_res
## 
##  Shapiro-Wilk normality test
## 
## data:  residuales
## W = 0.97846, p-value = 0.4625

Interpretación. La prueba de Shapiro-Wilk sobre los residuales del modelo arroja W = 0.978, p = 0.463. Dado que p > 0.05, no se rechaza la hipótesis nula de normalidad, se concluye que los residuales siguen razonablemente una distribución normal, lo cual es consistente con lo observado en el gráfico Q-Q, donde los puntos se alinean adecuadamente respecto a la línea teórica.

5.2.2 Homogeneidad de varianzas (Levene)

levene_res <- car::leveneTest(Shannon ~ Habitat, data = BD_biodiversidad)
levene_res

Interpretación. La prueba de Levene entrega F(3, 48) = 1.277, p = 0.293. Como p > 0.05, no hay evidencia para rechazar la hipótesis de homogeneidad de varianzas, por lo tanto se cumple el supuesto de homocedasticidad requerido por el ANOVA clásico.

Conclusión sobre supuestos: dado que ambos supuestos (normalidad y homogeneidad de varianzas) se comportan de manera razonable/adecuada para un tamaño de muestra de 13 observaciones por grupo, se procede con el ANOVA paramétrico de una vía y, de ser significativo, con la prueba post-hoc LSD.

5.3 Resultado del ANOVA

anova_tabla <- summary(modelo_shannon)[[1]]
kable(anova_tabla, caption = "Tabla ANOVA para Shannon ~ Hábitat", digits = 4)
Tabla ANOVA para Shannon ~ Hábitat
Df Sum Sq Mean Sq F value Pr(>F)
Habitat 3 14.8624 4.9541 60.6112 0
Residuals 48 3.9233 0.0817 NA NA 
# Tamaño del efecto: eta cuadrado
SS_habitat <- anova_tabla$`Sum Sq`[1]
SS_total   <- sum(anova_tabla$`Sum Sq`)
eta2 <- SS_habitat / SS_total

El estadístico F(3, 48) = 60.611, con un valor p = 2.38e-16, y un tamaño del efecto \(\eta^2\) = 0.791 (efecto grande según los criterios convencionales de Cohen).

Interpretación. Dado que el valor p es menor a 0.05, se rechaza la hipótesis nula de igualdad de medias. Esto indica que existen diferencias estadísticamente significativas en el índice de diversidad de Shannon entre al menos dos de los tipos de hábitat evaluados, y que el tipo de hábitat explica aproximadamente el 79.1% de la variabilidad total observada en el índice de Shannon, lo que constituye un efecto grande en términos prácticos.

5.4 Comparaciones post-hoc: prueba LSD (Fisher)

Dado que el ANOVA resultó significativo, se realizan comparaciones múltiples por pares mediante la prueba de la Diferencia Mínima Significativa (LSD de Fisher). Esta prueba utiliza la varianza combinada (pooled) del ANOVA y no ajusta el nivel de significancia por comparaciones múltiples; se justifica su uso porque el factor tiene únicamente 4 niveles balanceados y el objetivo es identificar de forma específica qué pares de hábitats difieren

lsd_res <- pairwise.t.test(BD_biodiversidad$Shannon, BD_biodiversidad$Habitat,
                            p.adjust.method = "none", pool.sd = TRUE)
lsd_res
## 
##  Pairwise comparisons using t tests with pooled SD 
## 
## data:  BD_biodiversidad$Shannon and BD_biodiversidad$Habitat 
## 
##                       Bosque primario Bosque secundario Sistema silvopastoril
## Bosque secundario     6.1e-05         -                 -                    
## Sistema silvopastoril 6.1e-10         0.0018            -                    
## Potrero               < 2e-16         2.3e-11           2.5e-06              
## 
## P value adjustment method: none
# Tabla ordenada de p-valores
lsd_tabla <- as.data.frame(as.table(lsd_res$p.value))
names(lsd_tabla) <- c("Habitat_1", "Habitat_2", "p_valor")
lsd_tabla <- lsd_tabla[!is.na(lsd_tabla$p_valor), ]
lsd_tabla$p_valor <- round(lsd_tabla$p_valor, 4)
lsd_tabla$Significativo <- ifelse(lsd_tabla$p_valor < 0.05, "Sí", "No")
kable(lsd_tabla, caption = "Comparaciones por pares - Prueba LSD (sin ajuste)")
Comparaciones por pares - Prueba LSD (sin ajuste)
Habitat_1 Habitat_2 p_valor Significativo
1 Bosque secundario Bosque primario 0.0001 Sí
2 Sistema silvopastoril Bosque primario 0.0000 Sí
3 Potrero Bosque primario 0.0000 Sí
5 Sistema silvopastoril Bosque secundario 0.0018 Sí
6 Potrero Bosque secundario 0.0000 Sí
9 Potrero Sistema silvopastoril 0.0000 Sí

5.4.1 Letras de significancia (compact letter display)

Como complemento visual a la tabla anterior, se resumen los resultados del LSD en la notación de letras compartidas, ampliamente usada en estudios ecológicos y agronómicos: hábitats que comparten al menos una letra no difieren significativamente entre sí; hábitats que no comparten ninguna letra sí difieren significativamente (p < 0.05).

# Construir la matriz completa de p-valores (simétrica) a partir del LSD
p_mat <- lsd_res$p.value
niveles <- levels(BD_biodiversidad$Habitat)
p_completa <- matrix(NA, nrow = length(niveles), ncol = length(niveles),
                      dimnames = list(niveles, niveles))
for (i in rownames(p_mat)) {
  for (j in colnames(p_mat)) {
    p_completa[i, j] <- p_mat[i, j]
    p_completa[j, i] <- p_mat[i, j]
  }
}
diag(p_completa) <- 1

# Vector de p-valores en el formato "grupo1-grupo2" que exige multcompLetters
p_vector <- p_completa[lower.tri(p_completa)]
names(p_vector) <- combn(niveles, 2, FUN = function(x) paste(x, collapse = "-"))

letras <- multcompLetters(p_vector, threshold = 0.05)$Letters

tabla_letras <- BD_biodiversidad %>%
  group_by(Habitat) %>%
  summarise(Media_Shannon = round(mean(Shannon), 2),
            DE = round(sd(Shannon), 2)) %>%
  mutate(Letra = letras[as.character(Habitat)]) %>%
  arrange(desc(Media_Shannon))

kable(tabla_letras, caption = "Medias del índice de Shannon por hábitat con letras de significancia (LSD, alfa = 0.05)")
Medias del índice de Shannon por hábitat con letras de significancia (LSD, alfa = 0.05)
Habitat Media_Shannon DE Letra
Bosque primario 2.47 0.38 a
Bosque secundario 1.97 0.22 b
Sistema silvopastoril 1.60 0.28 c
Potrero 1.00 0.24 d 
# Medias por grupo con error estándar y letras para gráfico resumen final
resumen_final <- BD_biodiversidad %>%
  group_by(Habitat) %>%
  summarise(Media = mean(Shannon), EE = sd(Shannon)/sqrt(n())) %>%
  mutate(Letra = letras[as.character(Habitat)])

ggplot(resumen_final, aes(x = Habitat, y = Media, fill = Habitat)) +
  geom_col(alpha = 0.8, width = 0.6) +
  geom_errorbar(aes(ymin = Media - EE, ymax = Media + EE), width = 0.15) +
  geom_text(aes(y = Media + EE + 0.12, label = Letra), size = 5, fontface = "bold") +
  scale_fill_brewer(palette = "Blues") +
  labs(title = "Media del índice de Shannon por hábitat (± error estándar)",
       subtitle = "Letras distintas indican diferencias significativas (LSD, p < 0.05)",
       x = "Tipo de hábitat", y = "Índice de Shannon promedio") +
  theme_minimal(base_size = 12) +
  theme(legend.position = "none",
        axis.text.x = element_text(angle = 15, hjust = 1))

Interpretación de las letras. Dado que en la sección anterior las 6 comparaciones por pares resultaron significativas, cada hábitat recibe una letra distinta (Bosque primario = a; Bosque secundario = b; Sistema silvopastoril = c; Potrero = d), lo que confirma de forma compacta que ningún par de hábitats comparte letra: los cuatro tipos de uso del suelo forman grupos estadísticamente separables entre sí en cuanto a su diversidad de Shannon, sin excepción.

Interpretación de las comparaciones post-hoc. La tabla de comparaciones por pares muestra qué combinaciones de hábitats presentan diferencias estadísticamente significativas (p < 0.05) en su índice de Shannon promedio. En términos generales, se espera —y debe verificarse contra la tabla anterior— que:

  • Bosque primario vs. Potrero y Bosque primario vs. Sistema silvopastoril: al representar los extremos del gradiente de intervención, son las comparaciones con mayor probabilidad de mostrar diferencias significativas, reflejando el fuerte impacto de la pérdida de cobertura boscosa sobre la diversidad de anfibios.
  • Bosque primario vs. Bosque secundario: una diferencia menor o no significativa sugeriría que la regeneración del bosque secundario permite recuperar parcialmente las condiciones de hábitat necesarias para mantener una diversidad de anfibios similar a la del bosque primario.
  • Sistema silvopastoril vs. Potrero: si no resulta significativa, indicaría que ambos usos de suelo con manejo pecuario generan un nivel de intervención comparable sobre la comunidad de anfibios, mientras que una diferencia significativa apoyaría la idea de que la presencia de árboles dispersos en el sistema silvopastoril aporta un beneficio adicional a la diversidad frente al potrero.

En conjunto, el patrón de significancia de la prueba LSD, leído junto con el gráfico de medias y errores estándar, permite establecer un ordenamiento jerárquico de los hábitats según su aporte a la diversidad de anfibios, información relevante para priorizar acciones de conservación o restauración dentro de la reserva forestal (por ejemplo, protección del bosque primario y promoción de la regeneración de bosque secundario, o del establecimiento de sistemas silvopastoriles frente al mantenimiento de potreros).

6 Conclusiones generales

  1. El tipo de hábitat es el factor determinante de la riqueza de especies de anfibios. Los resultados del análisis bivariado mostraron un patrón decreciente, ordenado y consistente en la riqueza de especies a lo largo del gradiente de intervención antrópica: el bosque primario presentó la mayor riqueza promedio (~14 especies por parcela), seguido del bosque secundario (~12), el sistema silvopastoril (~7) y, finalmente, el potrero, con la riqueza más baja (~5 especies). Esta disminución progresiva —y no un simple “salto” entre hábitat intervenido y no intervenido— sugiere que cada nivel de transformación del suelo elimina de forma acumulativa los microhábitats, refugios y fuentes de humedad que distintas especies de anfibios requieren para persistir.

  2. La altitud no es una explicación alternativa válida para estos resultados. Aunque se observó una correlación positiva entre altitud y diversidad de Shannon en el análisis exploratorio, esta relación resultó ser un reflejo indirecto del diseño del muestreo —los hábitats menos intervenidos (bosque primario) tienden a ubicarse en las zonas más altas de la reserva, y los más intervenidos (potrero) en las más bajas—, y no una relación causal independiente. Esto refuerza, en lugar de contradecir, la conclusión principal: el verdadero factor explicativo de la diversidad de anfibios en esta reserva es el uso del suelo, no la posición altitudinal en sí misma.

  3. Implicaciones para el manejo y la conservación de la reserva. En conjunto, esta evidencia estadística —respaldada por un tamaño de efecto excepcionalmente alto y comparaciones post-hoc consistentes— constituye un argumento robusto para priorizar la conservación estricta del bosque primario, dado que concentra tanto la mayor riqueza como la mayor diversidad de anfibios de la reserva. De igual manera, respalda la promoción de la regeneración natural hacia bosque secundario y la adopción de sistemas silvopastoriles con cobertura arbórea como alternativas productivas menos dañinas que el potrero convencional, ya que cada paso de recuperación de cobertura vegetal —según demuestran los datos— se traduce en una ganancia estadísticamente comprobable en la riqueza y diversidad de la comunidad de anfibios.