0.1 Introducción general

Este documento reúne los tres puntos del taller: (1) la relación entre la biomasa vegetal y las características del suelo (pH, salinidad, zinc y potasio), (2) el consumo de oxígeno en dos tipos de moluscos bajo distintas concentraciones de agua de mar, y (3) el efecto del uso del suelo sobre la biodiversidad de anfibios en una reserva forestal. Cada punto incluye un análisis exploratorio univariado y bivariado, seguido de la evaluación de un modelo ANOVA, la verificación de sus supuestos y las comparaciones post-hoc correspondientes.


1 Características del suelo y biomasa vegetal

Se analizarán 45 muestras para explorar la relación entre la biomasa de una planta forrajera y las características del suelo: pH, salinidad, zinc y potasio.

1.1 Carga de los datos

load("Salinidad.RData")

1.2 a. Análisis exploratorio univariado

Primero se examina la estructura de los datos y se calculan estadísticos descriptivos para cada variable.

str(Salinidad)
## 'data.frame':    45 obs. of  5 variables:
##  $ Biomasa  : num  765 954 828 755 896 ...
##  $ pH       : num  5 4.7 4.2 4.4 5.55 5.5 4.25 4.45 4.75 4.6 ...
##  $ Salinidad: int  33 35 32 30 33 33 36 30 38 30 ...
##  $ Zinc     : num  16.5 14 15.3 17.3 22.3 ...
##  $ Potasio  : num  1442 1299 1154 1045 522 ...

1.2.1 Estadísticos descriptivos

sal_descriptivos <- data.frame(
  Variable = names(Salinidad),
  n = sapply(Salinidad, function(x) sum(!is.na(x))),
  Faltantes = sapply(Salinidad, function(x) sum(is.na(x))),
  Media = sapply(Salinidad, mean, na.rm = TRUE),
  Mediana = sapply(Salinidad, median, na.rm = TRUE),
  DE = sapply(Salinidad, sd, na.rm = TRUE),
  Minimo = sapply(Salinidad, min, na.rm = TRUE),
  Maximo = sapply(Salinidad, max, na.rm = TRUE)
)

tabla_bonita(
  sal_descriptivos,
  digits = 2,
  row.names = FALSE,
  caption = "Tabla 1. Estadísticos descriptivos de las variables estudiadas."
)
Tabla 1. Estadísticos descriptivos de las variables estudiadas.
Variable n Faltantes Media Mediana DE Minimo Maximo
Biomasa 45 0 1082.17 991.83 546.29 369.82 2337.33
pH 45 0 4.61 4.45 1.25 3.20 7.45
Salinidad 45 0 30.27 30.00 3.72 24.00 38.00
Zinc 45 0 17.83 19.24 8.27 0.21 31.29
Potasio 45 0 797.38 773.30 297.58 350.73 1441.67

1.2.2 Distribución de las variables

sal_largos <- stack(Salinidad)
names(sal_largos) <- c("Valor", "Variable")

ggplot(sal_largos, aes(x = Valor)) +
  geom_histogram(
    bins = 8,
    fill = "steelblue",
    color = "white"
  ) +
  facet_wrap(~ Variable, scales = "free", ncol = 2) +
  labs(
    title = "Distribución de las características estudiadas",
    x = "Valor de la variable",
    y = "Número de muestras",
    caption = "Biomasa en gramos; pH sin unidades. Las demás unidades no se especifican en la guía."
  ) +
  tema_taller

1.2.3 Diagramas de caja

ggplot(sal_largos, aes(x = "", y = Valor)) +
  geom_boxplot(
    width = 0.4,
    fill = "lightblue",
    color = "gray30",
    outlier.color = "red",
    outlier.size = 2
  ) +
  facet_wrap(~ Variable, scales = "free_y", ncol = 2) +
  labs(
    title = "Dispersión y posibles valores atípicos",
    x = NULL,
    y = "Valor de la variable",
    caption = "Los puntos rojos indican valores fuera de los límites de 1,5 veces el rango intercuartílico."
  ) +
  tema_taller +
  theme(
    axis.text.x = element_blank(),
    axis.ticks.x = element_blank()
  )

1.2.4 Interpretación

Las cinco variables presentaron 45 observaciones y ningún dato faltante. La biomasa tuvo una media de 1082,17 g, una mediana de 991,83 g y una desviación estándar de 546,29 g, con valores entre 369,82 y 2337,33 g. Su distribución mostró asimetría hacia valores altos, evidenciando variación considerable entre las muestras.

El pH presentó una media de 4,61 y valores entre 3,20 y 7,45. Predominaron las condiciones ácidas, aunque algunas muestras presentaron valores cercanos a la neutralidad o ligeramente alcalinos. La salinidad tuvo una media de 30,27 y una desviación estándar de 3,72, con valores entre 24 y 38.

El zinc presentó una media de 17,83, una mediana de 19,24 y una desviación estándar de 8,27. Se observaron valores muy bajos y una concentración principal de observaciones en valores intermedios. El potasio tuvo una media de 797,38 y una desviación estándar de 297,58, con valores entre 350,73 y 1441,67; su histograma mostró una distribución irregular.

En conjunto, las muestras abarcan distintas condiciones del suelo y una amplia variación de biomasa. Este análisis describe cada variable por separado y no establece relaciones causales.

1.3 b. Análisis exploratorio bivariado

Se examina la relación de la biomasa con el pH, la salinidad y el zinc mediante gráficos de dispersión.

sal_relaciones <- data.frame(
  Biomasa = rep(Salinidad$Biomasa, times = 3),
  Valor = c(
    Salinidad$pH,
    Salinidad$Salinidad,
    Salinidad$Zinc
  ),
  Variable = rep(
    c("pH", "Salinidad", "Zinc"),
    each = nrow(Salinidad)
  )
)

ggplot(sal_relaciones, aes(x = Valor, y = Biomasa)) +
  geom_point(
    color = "steelblue",
    size = 2,
    alpha = 0.8
  ) +
  geom_smooth(
    method = "lm",
    formula = y ~ x,
    se = FALSE,
    color = "black"
  ) +
  facet_wrap(~ Variable, scales = "free_x", ncol = 1) +
  labs(
    title = "Relación de la biomasa con las características del suelo",
    x = "Valor de la característica del suelo",
    y = "Biomasa (g)"
  ) +
  tema_taller

1.3.1 Correlación entre la biomasa y las covariables

Se utiliza la correlación de Pearson para describir la dirección y la intensidad de la asociación lineal. Como análisis complementario, se calcula la correlación de Spearman, basada en los rangos de los datos.

sal_covariables <- c("pH", "Salinidad", "Zinc")

sal_correlaciones <- data.frame(
  Variable = sal_covariables,
  Pearson = sapply(sal_covariables, function(v) {
    cor(Salinidad$Biomasa, Salinidad[[v]], method = "pearson")
  }),
  Spearman = sapply(sal_covariables, function(v) {
    cor(Salinidad$Biomasa, Salinidad[[v]], method = "spearman")
  })
)

sal_correlaciones$Pearson_absoluto <- abs(sal_correlaciones$Pearson)

sal_correlaciones <- sal_correlaciones[
  order(-sal_correlaciones$Pearson_absoluto),
]

tabla_bonita(
  sal_correlaciones,
  digits = 3,
  row.names = FALSE,
  caption = "Tabla 2. Correlaciones de la biomasa con las características del suelo."
)
Tabla 2. Correlaciones de la biomasa con las características del suelo.
Variable Pearson Spearman Pearson_absoluto
pH 0.928 0.878 0.928
Zinc -0.781 -0.625 0.781
Salinidad -0.067 -0.157 0.067

1.3.2 Interpretación

El pH presentó la mayor asociación con la biomasa, tanto mediante Pearson (r = 0,928) como mediante Spearman (rho = 0,878). La relación fue positiva: las muestras con mayor pH tendieron a presentar mayor biomasa.

El zinc mostró una asociación negativa con la biomasa (r = -0,781; rho = -0,625). La salinidad presentó una asociación lineal muy débil (r = -0,067), acompañada de una correlación de Spearman también débil (rho = -0,157). Estos resultados coinciden con las tendencias observadas en los gráficos de dispersión.

Se seleccionó el pH para la comparación entre niveles porque presentó la mayor asociación entre las covariables evaluadas. Estas correlaciones describen asociaciones individuales y no separan la contribución del pH de la de otras características del suelo.

1.4 c. Comparación de la biomasa entre niveles de pH

Se seleccionó el pH porque presentó la mayor correlación absoluta de Pearson con la biomasa (r = 0,928), frente al zinc (r = -0,781) y la salinidad (r = -0,067). La correlación de Spearman mantuvo este orden de asociación.

1.4.1 Clasificación del pH en tres niveles

Se utilizaron los terciles de la distribución del pH para establecer los niveles bajo, medio y alto. Este criterio busca obtener grupos de tamaños aproximadamente similares, aunque los valores repetidos pueden producir diferencias en el número de muestras por grupo. Los niveles son relativos a esta muestra y no representan categorías universales de acidez del suelo.

sal_anova <- Salinidad

sal_cortes <- quantile(
  sal_anova$pH,
  probs = c(0, 1/3, 2/3, 1)
)

sal_cortes
##        0% 33.33333% 66.66667%      100% 
##  3.200000  3.883333  4.900000  7.450000
sal_anova$Nivel_pH <- cut(
  sal_anova$pH,
  breaks = sal_cortes,
  labels = c("Bajo", "Medio", "Alto"),
  include.lowest = TRUE,
  right = TRUE
)

table(sal_anova$Nivel_pH)
## 
##  Bajo Medio  Alto 
##    15    15    15

1.4.2 Modelo y verificación de supuestos

Se ajusta un ANOVA de una vía para comparar la biomasa media entre los tres niveles de pH. Se examinan la normalidad de los residuales y la homogeneidad de varianzas.

sal_modelo <- aov(Biomasa ~ Nivel_pH, data = sal_anova)

shapiro.test(residuals(sal_modelo))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(sal_modelo)
## W = 0.97316, p-value = 0.3749
car::leveneTest(
  Biomasa ~ Nivel_pH,
  data = sal_anova,
  center = median
)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value   Pr(>F)    
## group  2  8.6753 0.000702 ***
##       42                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

1.4.3 Gráficos de diagnóstico

par(mfrow = c(1, 2))

plot(sal_modelo, which = 1)
plot(sal_modelo, which = 2)

par(mfrow = c(1, 1))

1.4.4 Elección del análisis

La prueba de Shapiro-Wilk no mostró evidencia suficiente para rechazar la normalidad de los residuales (W = 0,97316; p = 0,3749). Sin embargo, la prueba de Levene detectó heterogeneidad de varianzas entre los grupos (F = 8,6753; p = 0,000702), consistente con la dispersión observada en el gráfico de residuales.

Por esta razón, se utiliza un ANOVA de Welch para comparar la biomasa media entre los tres niveles de pH, con un nivel de significancia de 0,05.

1.4.5 ANOVA de Welch

La hipótesis nula establece que la biomasa media es igual en los tres niveles de pH. La hipótesis alternativa establece que al menos una de estas medias es diferente.

sal_welch <- oneway.test(
  Biomasa ~ Nivel_pH,
  data = sal_anova,
  var.equal = FALSE
)

sal_welch
## 
##  One-way analysis of means (not assuming equal variances)
## 
## data:  Biomasa and Nivel_pH
## F = 34.727, num df = 2.000, denom df = 24.736, p-value = 6.58e-08

1.4.6 Comparaciones entre niveles de pH

El ANOVA de Welch detectó diferencias significativas en la biomasa media entre los niveles de pH (F = 34,727; gl = 2 y 24,736; p = 6,58 × 10^-8). Por tanto, se rechaza la hipótesis de igualdad de las tres medias.

Para identificar cuáles niveles difieren, se utiliza la prueba de Games-Howell. Se selecciona esta alternativa a la LSD convencional debido a la heterogeneidad de varianzas detectada por Levene.

sal_posthoc <- rstatix::games_howell_test(
  data = sal_anova,
  formula = Biomasa ~ Nivel_pH,
  conf.level = 0.95
)

tabla_bonita(
  sal_posthoc,
  digits = 6,
  caption = "Tabla 3. Comparaciones Games-Howell entre niveles de pH."
)
Tabla 3. Comparaciones Games-Howell entre niveles de pH.
.y. group1 group2 estimate conf.low conf.high p.adj p.adj.signif
Biomasa Bajo Medio 455.0859 264.9979 645.1739 0.000010 ****
Biomasa Bajo Alto 1012.3620 632.5370 1392.1870 0.000009 ****
Biomasa Medio Alto 557.2761 164.4401 950.1121 0.005071 **

1.4.7 Biomasa por nivel de pH

sal_biomasa_grupos <- split(
  sal_anova$Biomasa,
  sal_anova$Nivel_pH
)

sal_resumen_grupos <- data.frame(
  Nivel_pH = names(sal_biomasa_grupos),
  n = sapply(sal_biomasa_grupos, length),
  Media_g = sapply(sal_biomasa_grupos, mean),
  DE_g = sapply(sal_biomasa_grupos, sd),
  Mediana_g = sapply(sal_biomasa_grupos, median)
)

tabla_bonita(
  sal_resumen_grupos,
  digits = 2,
  row.names = FALSE,
  caption = "Tabla 4. Biomasa por nivel de pH."
)
Tabla 4. Biomasa por nivel de pH.
Nivel_pH n Media_g DE_g Mediana_g
Bajo 15 593.02 164.83 545.54
Medio 15 1048.11 244.91 1039.64
Alto 15 1605.39 547.60 1422.84
ggplot(sal_anova, aes(x = Nivel_pH, y = Biomasa)) +
  geom_boxplot(
    fill = "lightblue",
    width = 0.5,
    outlier.shape = NA
  ) +
  geom_point(
    position = position_jitter(width = 0.10, height = 0, seed = 123),
    alpha = 0.7,
    size = 2,
    color = "steelblue"
  ) +
  stat_summary(
    fun = mean,
    geom = "point",
    shape = 18,
    size = 4,
    color = "red"
  ) +
  labs(
    title = "Biomasa vegetal según el nivel de pH",
    x = "Nivel de pH",
    y = "Biomasa (g)",
    caption = "Los puntos azules representan muestras y los rombos rojos indican las medias."
  ) +
  tema_taller

1.4.8 Interpretación y conclusión

El ANOVA de Welch detectó diferencias en la biomasa media entre los niveles de pH (F = 34,727; gl = 2 y 24,736; p < 0,001). Las comparaciones Games-Howell indicaron diferencias significativas entre los tres pares de niveles.

La biomasa media del nivel medio superó a la del nivel bajo en 455,09 g, con un intervalo de confianza del 95 % de 265,00 a 645,17 g. El nivel alto superó al bajo en 1012,36 g (IC del 95 %: 632,54 a 1392,19 g) y al medio en 557,28 g (IC del 95 %: 164,44 a 950,11 g). Todos los valores de p ajustados fueron inferiores a 0,05.

El gráfico mostró un incremento de la biomasa media desde el nivel bajo hasta el alto, acompañado de mayor dispersión en este último grupo. Por tanto, en las muestras estudiadas, los niveles más altos de pH se asociaron con mayor biomasa vegetal.

La interpretación supone independencia entre las muestras, condición que depende del diseño de muestreo y no se verifica mediante Shapiro-Wilk o Levene. Además, los niveles de pH son categorías relativas a esta muestra y su construcción reduce la información de la variable continua.

Dado que el pH se seleccionó utilizando estos mismos datos, los resultados deben considerarse exploratorios. El análisis no permite atribuir causalidad ni descartar la influencia de otras características del suelo.


2 Consumo de oxígeno en moluscos según la concentración de agua de mar

Se evaluó el consumo de oxígeno (proporción de O2 por unidad de peso seco) en dos tipos de moluscos (A y B), sometidos a tres concentraciones distintas de agua de mar (100%, 75% y 50%). El objetivo es determinar, mediante un ANOVA de dos vías, si el tipo de molusco, la concentración de agua de mar y la interacción entre ambos factores afectan significativamente el consumo de oxígeno.

2.1 Carga de los datos

load("moluscos.RData")

BD_moluscos$molusco <- factor(BD_moluscos$molusco)
BD_moluscos$c_agua  <- factor(BD_moluscos$c_agua, levels = c(100, 75, 50),
                               labels = c("100%", "75%", "50%"))

str(BD_moluscos)
## tibble [48 × 3] (S3: tbl_df/tbl/data.frame)
##  $ c_agua : Factor w/ 3 levels "100%","75%","50%": 1 1 1 1 1 1 1 1 1 1 ...
##  $ molusco: Factor w/ 2 levels "A","B": 1 1 1 1 1 1 1 1 2 2 ...
##  $ cons_o : num [1:48] 7.16 8.26 6.78 14 13.6 11.1 8.93 9.66 6.14 6.14 ...

2.2 a. Análisis exploratorio univariado

2.2.1 Variables categóricas

table(BD_moluscos$molusco)
## 
##  A  B 
## 24 24
table(BD_moluscos$c_agua)
## 
## 100%  75%  50% 
##   16   16   16

El diseño experimental está balanceado: hay el mismo número de observaciones para cada tipo de molusco y para cada nivel de concentración, lo cual es deseable para que el ANOVA sea más robusto y fácil de interpretar.

2.2.2 Estadísticos descriptivos: consumo de O2

mol_descriptivos <- data.frame(
  Variable = "Consumo de O2",
  n = length(BD_moluscos$cons_o),
  Media = mean(BD_moluscos$cons_o),
  Mediana = median(BD_moluscos$cons_o),
  DE = sd(BD_moluscos$cons_o),
  CV_pct = 100 * sd(BD_moluscos$cons_o) / mean(BD_moluscos$cons_o),
  Minimo = min(BD_moluscos$cons_o),
  Maximo = max(BD_moluscos$cons_o)
)

tabla_bonita(
  mol_descriptivos,
  digits = 2,
  row.names = FALSE,
  caption = "Tabla 5. Estadísticos descriptivos del consumo de O2."
)
Tabla 5. Estadísticos descriptivos del consumo de O2.
Variable n Media Mediana DE CV_pct Minimo Maximo
Consumo de O2 48 9.3 9.7 3.68 39.58 1.8 18.8

2.2.3 Distribución del consumo de O2

ggplot(BD_moluscos, aes(x = cons_o)) +
  geom_histogram(bins = 8, fill = "steelblue", color = "white") +
  labs(
    title = "Distribución del consumo de O2",
    x = "Consumo de O2 / peso seco", y = "Número de parcelas"
  ) +
  tema_taller

ggplot(BD_moluscos, aes(x = "", y = cons_o)) +
  geom_boxplot(fill = "lightblue", color = "gray30",
               outlier.color = "red", outlier.size = 2, width = 0.4) +
  labs(title = "Boxplot general - Consumo de O2", x = NULL, y = "Consumo de O2") +
  tema_taller +
  theme(axis.text.x = element_blank(), axis.ticks.x = element_blank())

2.2.4 Interpretación

El consumo de oxígeno presenta una media cercana a 9,3 y una desviación estándar de 3,7, con un coeficiente de variación de aproximadamente 40%, lo que indica una variabilidad considerable entre individuos. El histograma y el boxplot muestran una leve asimetría hacia la derecha (algunos valores más altos alejan la cola superior de la distribución).

2.3 b. Análisis exploratorio bivariado

2.3.1 Consumo de O2 según la concentración de agua de mar

tabla_bonita(
  aggregate(cons_o ~ c_agua, data = BD_moluscos, FUN = mean),
  digits = 2, row.names = FALSE,
  caption = "Tabla 6. Consumo medio de O2 según concentración."
)
Tabla 6. Consumo medio de O2 según concentración.
c_agua cons_o
100% 8.67
75% 6.99
50% 12.25
ggplot(BD_moluscos, aes(x = c_agua, y = cons_o)) +
  geom_boxplot(fill = "lightblue", color = "gray30", outlier.color = "red") +
  labs(
    title = "Consumo de O2 según concentración de agua de mar",
    x = "Concentración", y = "Consumo de O2"
  ) +
  tema_taller

El consumo promedio no disminuye de forma lineal con la concentración: baja levemente de 100% a 75%, y luego aumenta marcadamente al 50%, sugiriendo una posible respuesta de estrés osmótico ante la dilución del agua de mar.

2.3.2 ¿Es el mismo patrón para cada tipo de molusco?

mol_medias_cruzadas <- aggregate(cons_o ~ c_agua + molusco, data = BD_moluscos, FUN = mean)

tabla_bonita(
  mol_medias_cruzadas,
  digits = 2, row.names = FALSE,
  caption = "Tabla 7. Consumo medio de O2 por concentración y tipo de molusco."
)
Tabla 7. Consumo medio de O2 por concentración y tipo de molusco.
c_agua molusco cons_o
100% A 9.94
75% A 7.89
50% A 12.18
100% B 7.41
75% B 6.10
50% B 12.33
ggplot(BD_moluscos, aes(x = c_agua, y = cons_o, fill = molusco)) +
  geom_boxplot(color = "gray30") +
  scale_fill_manual(values = c("A" = "steelblue", "B" = "lightblue")) +
  labs(
    title = "Consumo de O2 por concentración y tipo de molusco",
    x = "Concentración", y = "Consumo de O2", fill = "Molusco"
  ) +
  tema_taller

ggplot(mol_medias_cruzadas, aes(x = c_agua, y = cons_o, color = molusco, group = molusco)) +
  geom_line(linewidth = 1) +
  geom_point(size = 3) +
  scale_color_manual(values = c("A" = "steelblue", "B" = "tomato")) +
  labs(
    title = "Gráfico de interacción: concentración x tipo de molusco",
    x = "Concentración", y = "Consumo medio de O2", color = "Molusco"
  ) +
  tema_taller

Las líneas del gráfico de interacción son aproximadamente paralelas entre el molusco A y el B, lo que sugiere de forma exploratoria que el efecto de la concentración sobre el consumo de oxígeno es similar en ambos tipos de molusco (no se observa un cruce marcado que indique una interacción fuerte). Esto se confirmará formalmente con el ANOVA.

2.4 c. ANOVA de dos vías

mol_modelo <- aov(cons_o ~ molusco * c_agua, data = BD_moluscos)
summary(mol_modelo)
##                Df Sum Sq Mean Sq F value   Pr(>F)    
## molusco         1   23.2   23.23   2.651    0.111    
## c_agua          2  230.8  115.41  13.171 3.63e-05 ***
## molusco:c_agua  2   15.4    7.68   0.876    0.424    
## Residuals      42  368.0    8.76                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Los resultados muestran que:

  • El factor molusco no tiene un efecto significativo sobre el consumo de O2 (p > 0.05).
  • El factor concentración sí tiene un efecto significativo sobre el consumo de O2 (p < 0.001).
  • La interacción molusco:concentración no es significativa (p > 0.05), confirmando que el efecto de la concentración es equivalente para ambos tipos de molusco.

2.4.1 Verificación de supuestos

Normalidad de los residuales

shapiro.test(residuals(mol_modelo))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(mol_modelo)
## W = 0.95824, p-value = 0.08571
ggplot(data.frame(residuos = residuals(mol_modelo)), aes(sample = residuos)) +
  stat_qq(color = "steelblue") +
  stat_qq_line(color = "red") +
  labs(title = "QQ-plot de los residuales", x = "Cuantiles teóricos", y = "Residuales") +
  tema_taller

Homogeneidad de varianzas

leveneTest(cons_o ~ molusco * c_agua, data = BD_moluscos)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  5  0.1723 0.9715
##       42

Ambos supuestos se cumplen (p > 0,05 en las dos pruebas), por lo que los resultados del ANOVA son confiables.

2.4.2 Comparaciones post-hoc (prueba de LSD)

Dado que la concentración fue el único efecto significativo con más de dos niveles, las comparaciones post-hoc se realizan sobre este factor:

mol_lsd <- LSD.test(mol_modelo, "c_agua")

tabla_bonita(
  mol_lsd$groups,
  digits = 3,
  caption = "Tabla 8. Comparaciones LSD del consumo de O2 según concentración."
)
Tabla 8. Comparaciones LSD del consumo de O2 según concentración.
cons_o groups
50% 12.251 a
100% 8.671 b
75% 6.992 b

La prueba de LSD muestra que el consumo de O2 al 50% de concentración es significativamente distinto (más alto) que al 100% y al 75%, mientras que estos dos últimos no difieren significativamente entre sí (grupos con la misma letra en la tabla no difieren).

2.5 Conclusiones

El consumo de oxígeno en los moluscos estudiados depende significativamente de la concentración de agua de mar, pero no del tipo de molusco considerado de forma aislada, y el patrón de respuesta a la concentración es equivalente entre el molusco A y el B (interacción no significativa). Específicamente, la dilución del agua de mar al 50% genera un aumento significativo en el consumo de oxígeno respecto a las concentraciones de 100% y 75%, lo cual es consistente con una respuesta fisiológica de mayor actividad metabólica ante el estrés osmótico generado por la menor salinidad del medio.


3 Efecto del uso del suelo sobre la biodiversidad de anfibios

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 (13 por hábitat) en cuatro tipos de hábitat que representan un gradiente de intervención antrópica: bosque primario, bosque secundario, sistema silvopastoril y potrero. En cada parcela se registró la riqueza de especies, el índice de diversidad de Shannon-Wiener y la altitud del sitio.

3.1 Carga de los datos

load("Biodiversidad.RData")

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] "" "" "" "" ...

3.2 a. Análisis exploratorio univariado

3.2.1 Estadísticos descriptivos

bio_vars <- c("Riqueza", "Shannon", "Altitud")

bio_descriptivos <- data.frame(
  Variable = bio_vars,
  n = sapply(BD_biodiversidad[bio_vars], length),
  Media = sapply(BD_biodiversidad[bio_vars], mean),
  Mediana = sapply(BD_biodiversidad[bio_vars], median),
  DE = sapply(BD_biodiversidad[bio_vars], sd),
  CV_pct = sapply(BD_biodiversidad[bio_vars], function(x) 100 * sd(x) / mean(x)),
  Minimo = sapply(BD_biodiversidad[bio_vars], min),
  Maximo = sapply(BD_biodiversidad[bio_vars], max)
)

tabla_bonita(
  bio_descriptivos,
  digits = 2, row.names = FALSE,
  caption = "Tabla 9. Estadísticos descriptivos de riqueza, Shannon y altitud."
)
Tabla 9. Estadísticos descriptivos de riqueza, Shannon y altitud.
Variable n Media Mediana DE CV_pct Minimo Maximo
Riqueza 52 9.73 9.50 4.37 44.92 2.00 21.00
Shannon 52 1.76 1.85 0.61 34.44 0.67 2.98
Altitud 52 1243.89 1254.80 172.71 13.88 950.90 1600.80

3.2.2 Distribución de las variables

bio_largos <- data.frame(
  Valor = c(BD_biodiversidad$Riqueza, BD_biodiversidad$Shannon, BD_biodiversidad$Altitud),
  Variable = rep(c("Riqueza", "Shannon", "Altitud"), each = nrow(BD_biodiversidad))
)
ggplot(bio_largos, aes(x = Valor)) +
  geom_histogram(bins = 8, fill = "steelblue", color = "white") +
  facet_wrap(~ Variable, scales = "free", ncol = 3) +
  labs(
    title = "Distribución de la riqueza, la diversidad de Shannon y la altitud",
    x = "Valor de la variable", y = "Número de parcelas"
  ) +
  tema_taller

ggplot(bio_largos, aes(x = "", y = Valor)) +
  geom_boxplot(fill = "lightblue", color = "gray30", outlier.color = "red", width = 0.4) +
  facet_wrap(~ Variable, scales = "free_y", ncol = 3) +
  labs(title = "Dispersión y posibles valores atípicos", x = NULL, y = "Valor de la variable") +
  tema_taller +
  theme(axis.text.x = element_blank(), axis.ticks.x = element_blank())

3.2.3 Interpretación

La riqueza presenta la mayor variabilidad relativa (CV ≈ 45%), seguida de la diversidad de Shannon (CV ≈ 34%). La altitud es la variable más homogénea (CV ≈ 14%). Las tres variables muestran distribuciones razonablemente simétricas, sin asimetrías marcadas.

3.3 b. Análisis exploratorio bivariado

3.3.1 Riqueza y diversidad de Shannon según el tipo de hábitat

ggplot(BD_biodiversidad, aes(x = Habitat, y = Riqueza)) +
  geom_boxplot(fill = "lightblue", color = "gray30", outlier.color = "red") +
  labs(title = "Riqueza de especies por tipo de hábitat", x = NULL, y = "Riqueza") +
  tema_taller +
  theme(axis.text.x = element_text(angle = 20, hjust = 1))

ggplot(BD_biodiversidad, aes(x = Habitat, y = Shannon)) +
  geom_boxplot(fill = "lightblue", color = "gray30", outlier.color = "red") +
  labs(title = "Diversidad de Shannon por tipo de hábitat", x = NULL, y = "Índice de Shannon") +
  tema_taller +
  theme(axis.text.x = element_text(angle = 20, hjust = 1))

Tanto la riqueza como la diversidad de Shannon disminuyen de forma escalonada siguiendo el gradiente de intervención antrópica: bosque primario > bosque secundario > sistema silvopastoril > potrero.

3.3.2 Relación entre altitud y diversidad

ggplot(BD_biodiversidad, aes(x = Altitud, y = Shannon)) +
  geom_point(color = "steelblue", size = 2, alpha = 0.8) +
  geom_smooth(method = "lm", formula = y ~ x, se = FALSE, color = "black") +
  labs(
    title = "Relación entre altitud y diversidad de Shannon",
    x = "Altitud (m.s.n.m.)", y = "Índice de Shannon"
  ) +
  tema_taller

cor.test(BD_biodiversidad$Altitud, BD_biodiversidad$Shannon, method = "pearson")
## 
##  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

Existe una correlación positiva fuerte entre la altitud y la diversidad de Shannon (r ≈ 0,75, p < 0,001): las parcelas ubicadas a mayor altitud tienden a presentar mayor diversidad. Esto probablemente refleja que, en esta reserva, los hábitats mejor conservados (bosque primario) se concentran en las zonas más altas.

3.4 c. ANOVA de una vía: Shannon según hábitat

bio_modelo <- aov(Shannon ~ Habitat, data = BD_biodiversidad)
summary(bio_modelo)
##             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

El hábitat tiene un efecto altamente significativo sobre la diversidad de Shannon (p < 0,001).

3.4.1 Verificación de supuestos

Normalidad de los residuales

shapiro.test(residuals(bio_modelo))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(bio_modelo)
## W = 0.97846, p-value = 0.4625
ggplot(data.frame(residuos = residuals(bio_modelo)), aes(sample = residuos)) +
  stat_qq(color = "steelblue") +
  stat_qq_line(color = "red") +
  labs(title = "QQ-plot de los residuales", x = "Cuantiles teóricos", y = "Residuales") +
  tema_taller

Homogeneidad de varianzas

leveneTest(Shannon ~ Habitat, data = BD_biodiversidad)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  3  1.2768  0.293
##       48

Ambos supuestos se cumplen (p > 0,05 en las dos pruebas), por lo que el resultado del ANOVA es confiable.

3.4.2 Comparaciones post-hoc (prueba de LSD)

bio_lsd <- LSD.test(bio_modelo, "Habitat")

tabla_bonita(
  bio_lsd$groups,
  digits = 3,
  caption = "Tabla 10. Comparaciones LSD de la diversidad de Shannon entre hábitats."
)
Tabla 10. Comparaciones LSD de la diversidad de Shannon entre hábitats.
Shannon groups
Bosque primario 2.468 a
Bosque secundario 1.975 b
Sistema silvopastoril 1.603 c
Potrero 1.004 d
ggplot(BD_biodiversidad, aes(x = Habitat, y = Shannon)) +
  geom_boxplot(fill = "lightblue", width = 0.5, outlier.shape = NA) +
  geom_point(
    position = position_jitter(width = 0.10, height = 0, seed = 123),
    alpha = 0.7, size = 2, color = "steelblue"
  ) +
  stat_summary(fun = mean, geom = "point", shape = 18, size = 4, color = "red") +
  labs(
    title = "Diversidad de Shannon según el tipo de hábitat",
    x = NULL, y = "Índice de Shannon",
    caption = "Los puntos azules representan parcelas y los rombos rojos indican las medias."
  ) +
  tema_taller +
  theme(axis.text.x = element_text(angle = 20, hjust = 1))

La prueba de LSD muestra que los cuatro tipos de hábitat difieren significativamente entre sí en su diversidad de Shannon (ninguna comparación por pares resulta no significativa; grupos con distinta letra en la tabla difieren).

3.5 Conclusiones

El uso del suelo afecta de manera significativa la diversidad de anfibios en la reserva. La diversidad de Shannon disminuye de forma gradual y estadísticamente distinguible a lo largo de todo el gradiente de intervención antrópica: el bosque primario presenta la mayor diversidad, seguido del bosque secundario, el sistema silvopastoril y, por último, el potrero con la menor diversidad. Adicionalmente, se observa una asociación positiva entre la altitud y la diversidad, consistente con la distribución espacial de los hábitats mejor conservados en las zonas más altas de la reserva.