library(psych)
library(knitr)
library(kableExtra)
library(ggplot2)
library(dplyr)
library(car)
library(agricolae)
library(gridExtra)

cv <- function(x) (sd(x) / mean(x)) * 100

Este informe reúne tres estudios independientes (suelos y biomasa vegetal, fisiología de moluscos y biodiversidad de anfibios) con el fin de aplicar análisis exploratorio, ANOVA y pruebas post-ANOVA sobre datos biológicos y ambientales. En cada punto se sigue la misma lógica: describir las variables, evaluar diferencias entre grupos verificando los supuestos del modelo, y confirmar con pruebas post-hoc cuáles grupos difieren entre sí, interpretando cada resultado en su contexto.

Punto 1 – Salinidad

Para estudiar la relación entre ciertas características del suelo y la producción de biomasa (g) de una planta forrajera natural, se recolectaron 45 muestras en diferentes ambientes. En cada una se registró la biomasa producida (variable respuesta) y cuatro covariables del suelo: pH, salinidad, zinc y potasio.

45Muestras

4Covariables del suelo

1Variable respuesta (Biomasa)

a. Análisis Exploratorio Univariado

Panorama general de las variables

load("Salinidad.RData")

tabla_describe_sal <- describe(Salinidad)

kable(round(tabla_describe_sal, 2), 
      caption = "Estadísticos descriptivos",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE) %>%
  scroll_box(width = "100%")

Estadísticos descriptivos
	vars	n	mean	sd	median	trimmed	mad	min	max	range	skew	kurtosis	se
Biomasa	1	45	1082.17	546.29	991.83	1022.02	526.40	369.82	2337.33	1967.50	0.91	-0.02	81.44
pH	2	45	4.61	1.25	4.45	4.46	1.41	3.20	7.45	4.25	0.87	-0.03	0.19
Salinidad	3	45	30.27	3.72	30.00	30.14	4.45	24.00	38.00	14.00	0.31	-1.04	0.55
Zinc	4	45	17.83	8.27	19.24	18.40	6.62	0.21	31.29	31.08	-0.66	-0.06	1.23
Potasio	5	45	797.38	297.58	773.30	778.48	352.09	350.73	1441.67	1090.94	0.48	-1.00	44.36

No hay datos faltantes en ninguna de las 45 muestras. Las cinco variables están en escalas muy distintas entre sí (Biomasa y Potasio en cientos, pH en unidades de 3 a 7), lo que hace poco justo comparar directamente sus desviaciones estándar; por eso a continuación se calcula el coeficiente de variación de cada una.

¿Qué tan dispersas están las variables?

cv_tabla_sal <- data.frame(
  Variable = names(sapply(Salinidad, cv)),
  CV = round(sapply(Salinidad, cv), 2)
)

kable(cv_tabla_sal, 
      caption = "Coeficiente de variación (%)",
      col.names = c("Variable", "CV (%)"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Coeficiente de variación (%)
Variable	CV (%)
Biomasa	50.48
pH	27.22
Salinidad	12.29
Zinc	46.40
Potasio	37.32

La biomasa (CV ≈ 50%) y el zinc (CV ≈ 46%) son, por mucho, las variables más dispersas del conjunto; la salinidad (CV ≈ 12%) es la más estable entre muestras.

Formas de distribución y presencia de atípicos

p1 <- ggplot(Salinidad, aes(x = Biomasa)) +
  geom_histogram(bins = 10, fill = "#622B51", color = "white") +
  theme_minimal() + labs(title = "Biomasa", x = "Biomasa (g)", y = "Frecuencia")

p2 <- ggplot(Salinidad, aes(x = pH)) +
  geom_histogram(bins = 10, fill = "#84477C", color = "white") +
  theme_minimal() + labs(title = "pH", x = "pH", y = "Frecuencia")

p3 <- ggplot(Salinidad, aes(x = Salinidad)) +
  geom_histogram(bins = 10, fill = "#534166", color = "white") +
  theme_minimal() + labs(title = "Salinidad", x = "Salinidad", y = "Frecuencia")

p4 <- ggplot(Salinidad, aes(x = Zinc)) +
  geom_histogram(bins = 10, fill = "#7C92A1", color = "white") +
  theme_minimal() + labs(title = "Zinc", x = "Zinc", y = "Frecuencia")

p5 <- ggplot(Salinidad, aes(x = Potasio)) +
  geom_histogram(bins = 10, fill = "#401C39", color = "white") +
  theme_minimal() + labs(title = "Potasio", x = "Potasio", y = "Frecuencia")

grid.arrange(p1, p2, p3, p4, p5, ncol = 3)

Biomasa y pH muestran una ligera cola hacia la derecha; Zinc, en cambio, tiene una cola hacia la izquierda, pista de que hay algunas muestras con valores muy bajos, algo que se confirma con los boxplots a continuación.

b1 <- ggplot(Salinidad, aes(y = Biomasa)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Biomasa") + theme(axis.text.x = element_blank())

b2 <- ggplot(Salinidad, aes(y = pH)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "pH") + theme(axis.text.x = element_blank())

b3 <- ggplot(Salinidad, aes(y = Salinidad)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Salinidad") + theme(axis.text.x = element_blank())

b4 <- ggplot(Salinidad, aes(y = Zinc)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Zinc") + theme(axis.text.x = element_blank())

b5 <- ggplot(Salinidad, aes(y = Potasio)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Potasio") + theme(axis.text.x = element_blank())

grid.arrange(b1, b2, b3, b4, b5, ncol = 3)

Solo el Zinc muestra puntos por fuera de los bigotes de la caja; el resto de variables no presenta valores fuera de rango. A continuación se identifican esos casos puntualmente.

outliers_zinc_sal <- boxplot.stats(Salinidad$Zinc)$out

outliers_tabla_sal <- data.frame(
  Muestra = seq_along(outliers_zinc_sal),
  Zinc = round(outliers_zinc_sal, 4)
)

kable(outliers_tabla_sal, 
      caption = "Valores atípicos detectados en Zinc",
      col.names = c("N°", "Valor de Zinc"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Valores atípicos detectados en Zinc
N°	Valor de Zinc
1	0.3729
2	0.2703
3	0.3205
4	0.2648
5	0.2105

Son 5 muestras con Zinc entre 0.21 y 0.37, muy por debajo del resto, siendo un grupo de valores particular dentro de la base.

En conjunto, la biomasa y el zinc son las variables con mayor variabilidad relativa (CV ≈ 50% y 46%), mientras que la salinidad es la más estable (CV ≈ 12%); el pH tiene una dispersión intermedia (CV ≈ 27%) y el potasio se ubica también en un nivel medio-alto (CV ≈ 37%). El zinc es la única variable con valores atípicos identificados.

b. Análisis Exploratorio Bivariado

Matriz de correlaciones

matriz_cor_sal <- cor(Salinidad[, c("Biomasa", "pH", "Salinidad", "Zinc", "Potasio")], 
                   method = "pearson")

kable(round(matriz_cor_sal, 3), 
      caption = "Matriz de correlación de Pearson",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Matriz de correlación de Pearson
	Biomasa	pH	Salinidad	Zinc	Potasio
Biomasa	1.000	0.928	-0.067	-0.781	-0.073
pH	0.928	1.000	-0.045	-0.720	0.032
Salinidad	-0.067	-0.045	1.000	-0.427	-0.020
Zinc	-0.781	-0.720	-0.427	1.000	0.079
Potasio	-0.073	0.032	-0.020	0.079	1.000

La covariable más relacionada con la biomasa

cor_biomasa_sal <- data.frame(
  Variable = c("pH", "Salinidad", "Zinc", "Potasio"),
  Correlacion = round(matriz_cor_sal["Biomasa", c("pH", "Salinidad", "Zinc", "Potasio")], 3)
)

kable(cor_biomasa_sal, 
      caption = "Correlación de Pearson entre Biomasa y cada covariable", 
      col.names = c("Variable", "r (Pearson)"), 
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Correlación de Pearson entre Biomasa y cada covariable
Variable	r (Pearson)
pH	0.928
Salinidad	-0.067
Zinc	-0.781
Potasio	-0.073

El pH (r = 0.928) y el Zinc (r = -0.781) destacan como las covariables de mayor peso: la primera con relación fuerte y positiva, la segunda fuerte pero negativa. Salinidad (r = -0.067) y Potasio (r = -0.073) prácticamente no se relacionan linealmente con la biomasa. Por tener la correlación más fuerte, el pH es la covariable elegida para el análisis de varianza del literal c).

Comprobación visual: diagramas de dispersión

s1 <- ggplot(Salinidad, aes(x = pH, y = Biomasa)) +
  geom_point(color = "#622B51", size = 2) + theme_minimal() +
  labs(title = "Biomasa vs pH", x = "pH", y = "Biomasa")

s2 <- ggplot(Salinidad, aes(x = Salinidad, y = Biomasa)) +
  geom_point(color = "#534166", size = 2) + theme_minimal() +
  labs(title = "Biomasa vs Salinidad", x = "Salinidad", y = "Biomasa")

s3 <- ggplot(Salinidad, aes(x = Zinc, y = Biomasa)) +
  geom_point(color = "#84477C", size = 2) + theme_minimal() +
  labs(title = "Biomasa vs Zinc", x = "Zinc", y = "Biomasa")

s4 <- ggplot(Salinidad, aes(x = Potasio, y = Biomasa)) +
  geom_point(color = "#7C92A1", size = 2) + theme_minimal() +
  labs(title = "Biomasa vs Potasio", x = "Potasio", y = "Biomasa")

grid.arrange(s1, s2, s3, s4, ncol = 2)

Los diagramas confirman visualmente lo que ya decía la matriz: la nube de puntos de Biomasa vs. pH sube de forma consistente, la de Biomasa vs. Zinc baja, y las otras dos (Salinidad y Potasio) se ven dispersas sin ninguna tendencia clara.

c. ANOVA y Postanova

Construcción de los niveles de pH

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

terciles_tabla_sal <- data.frame(
  Percentil = names(terciles_sal),
  pH = round(as.numeric(terciles_sal), 2)
)

kable(terciles_tabla_sal, 
      caption = "Puntos de corte (terciles) del pH",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Puntos de corte (terciles) del pH
Percentil	pH
0%	3.20
33.33333%	3.88
66.66667%	4.90
100%	7.45

Salinidad$pH_nivel <- cut(Salinidad$pH, 
                            breaks = terciles_sal, 
                            labels = c("Bajo", "Medio", "Alto"), 
                            include.lowest = TRUE)

conteo_tabla_sal <- as.data.frame(table(Salinidad$pH_nivel))

kable(conteo_tabla_sal, 
      caption = "Número de observaciones por nivel de pH",
      col.names = c("Nivel", "N"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Número de observaciones por nivel de pH
Nivel	N
Bajo	15
Medio	15
Alto	15

Los terciles dejan tres grupos parejos (15 muestras cada uno), condición conveniente para las comparaciones que siguen.

Ajuste del modelo ANOVA

anova_sal <- aov(Biomasa ~ pH_nivel, data = Salinidad)

anova_tabla_sal <- as.data.frame(summary(anova_sal)[[1]])
anova_tabla_sal <- round(anova_tabla_sal, 4)

kable(anova_tabla_sal, 
      caption = "Tabla ANOVA: Biomasa ~ pH_nivel",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Tabla ANOVA: Biomasa ~ pH_nivel
	Df	Sum Sq	Mean Sq	F value	Pr(>F)
pH_nivel	2	7712683	3856341.6	29.8928	0
Residuals	42	5418235	129005.6	NA	NA

El nivel de pH resulta con un efecto altamente significativo sobre la biomasa (F(2,42) = 29.89, p < 0.001).

¿Es válido el modelo? Verificación de supuestos

shapiro_res_sal <- shapiro.test(residuals(anova_sal))
levene_res_sal <- leveneTest(Biomasa ~ pH_nivel, data = Salinidad)

supuestos_tabla_sal <- data.frame(
  Prueba = c("Shapiro-Wilk (normalidad)", "Levene (homogeneidad de varianzas)"),
  Estadistico = c(round(shapiro_res_sal$statistic, 4), round(levene_res_sal$`F value`[1], 4)),
  p_valor = c(round(shapiro_res_sal$p.value, 4), round(levene_res_sal$`Pr(>F)`[1], 4))
)

kable(supuestos_tabla_sal, 
      caption = "Verificación de supuestos del ANOVA",
      col.names = c("Prueba", "Estadístico", "p-valor"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Verificación de supuestos del ANOVA
Prueba	Estadístico	p-valor
Shapiro-Wilk (normalidad)	0.9732	0.3749
Levene (homogeneidad de varianzas)	8.6753	0.0007

Con p = 0.375, Shapiro-Wilk no rechaza la normalidad de los residuales. Levene sí rechaza la homogeneidad de varianzas (p = 0.0007), una limitación del análisis que conviene tener presente al leer los resultados que siguen.

Diferencias entre niveles: prueba LSD

posanova_sal <- LSD.test(anova_sal, trt = "pH_nivel")

kable(as.data.frame(posanova_sal$statistics), 
      caption = "Estadísticos generales del modelo",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Estadísticos generales del modelo
MSerror	Df	Mean	CV	t.value	LSD
129005.6	42	1082.173	33.19004	2.018082	264.6747

kable(as.data.frame(posanova_sal$parameters), 
      caption = "Parámetros de la prueba LSD",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Parámetros de la prueba LSD
test	p.ajusted	name.t	ntr	alpha
Fisher-LSD	none	pH_nivel	3	0.05

grupos_tabla_sal <- posanova_sal$groups
grupos_tabla_sal$Nivel <- rownames(grupos_tabla_sal)
grupos_tabla_sal <- grupos_tabla_sal[, c("Nivel", "Biomasa", "groups")]

kable(grupos_tabla_sal, 
      caption = "Agrupamiento LSD — Nivel de pH",
      col.names = c("Nivel de pH", "Biomasa promedio", "Grupo"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = FALSE)

Agrupamiento LSD — Nivel de pH
Nivel de pH	Biomasa promedio	Grupo
Alto	1605.3853	a
Medio	1048.1093	b
Bajo	593.0233	c

Los tres niveles reciben letras distintas (a, b, c): Alto (1605.39 g), Medio (1048.11 g) y Bajo (593.02 g) difieren entre sí, en línea con la correlación positiva ya observada entre pH y biomasa.

Síntesis y hallazgos principales

Univariado: Biomasa (CV ≈ 50%) y Zinc (CV ≈ 46%) son las variables más dispersas; Salinidad (CV ≈ 12%) la más estable. Zinc es la única con valores atípicos (5 muestras cercanas a cero).
Bivariado: pH (r = 0.928) y Zinc (r = -0.781) son las covariables con relación fuerte a la biomasa; Salinidad y Potasio no muestran relación lineal relevante.
ANOVA: el nivel de pH genera diferencias altamente significativas en la biomasa (F(2,42) = 29.89, p < 0.001). La normalidad se cumple, pero no la homogeneidad de varianzas (Levene, p = 0.0007). LSD separa los tres niveles entre sí.
Lectura general: el pH del suelo es la variable con mayor peso explicativo sobre la producción de biomasa forrajera, seguido del Zinc; suelos menos ácidos y con menos zinc rinden más biomasa.

Punto 2 – Moluscos

Dos tipos de moluscos (A y B) fueron expuestos a tres concentraciones distintas de agua de mar (100%, 75% y 50%), registrando el consumo de oxígeno como la proporción de O2 por unidad de peso seco del organismo. El experimento sigue una estructura factorial 2×3 (tipo de molusco por concentración de agua), con 8 réplicas por combinación y 48 observaciones en total.

48Observaciones

2Tipos de molusco

3Concentraciones de agua

1Variable respuesta (Consumo O2)

load("moluscos.RData")

BD_moluscos <- BD_moluscos %>%
  mutate(
    c_agua  = factor(c_agua, levels = c(100, 75, 50)),
    molusco = factor(molusco)
  )

a. Análisis Exploratorio Univariado

Balance del diseño experimental

tabla_molusco_mol <- as.data.frame(table(BD_moluscos$molusco))
kable(tabla_molusco_mol,
      caption = "Frecuencia por tipo de molusco",
      col.names = c("Molusco", "N"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Frecuencia por tipo de molusco
Molusco	N
A	24
B	24

tabla_agua_mol <- as.data.frame(table(BD_moluscos$c_agua))
kable(tabla_agua_mol,
      caption = "Frecuencia por concentración de agua de mar",
      col.names = c("Concentración (%)", "N"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Frecuencia por concentración de agua de mar
Concentración (%)	N
100	16
75	16
50	16

El diseño queda perfectamente balanceado: 24 observaciones por tipo de molusco y 16 por cada concentración de agua de mar.

El consumo de oxígeno en cifras

resumen_o2_mol <- data.frame(
  Estadistico = c("Mínimo", "1er Cuartil", "Mediana", "Media", "3er Cuartil",
                   "Máximo", "Desviación estándar"),
  Valor = round(c(min(BD_moluscos$cons_o), quantile(BD_moluscos$cons_o, 0.25),
                    median(BD_moluscos$cons_o), mean(BD_moluscos$cons_o),
                    quantile(BD_moluscos$cons_o, 0.75), max(BD_moluscos$cons_o),
                    sd(BD_moluscos$cons_o)), 2)
)

kable(resumen_o2_mol,
      caption = "Estadísticos descriptivos del consumo de O2",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Estadísticos descriptivos del consumo de O2
Estadistico	Valor
Mínimo	1.80
1er Cuartil	6.31
Mediana	9.70
Media	9.30
3er Cuartil	11.23
Máximo	18.80
Desviación estándar	3.68

ggplot(BD_moluscos, aes(x = cons_o)) +
  geom_histogram(bins = 10, fill = "#622B51", color = "white") +
  labs(title = "Distribución del consumo de O2", x = "Consumo de O2", y = "Frecuencia") +
  theme_minimal()

ggplot(BD_moluscos, aes(y = cons_o)) +
  geom_boxplot(fill = "#D9CFE0") +
  labs(title = "Boxplot del consumo de O2", y = "Consumo de O2") +
  theme_minimal()

La media (9.30) y la mediana (9.70) del consumo de O2 son cercanas, con una dispersión considerable (sd = 3.68, CV ≈ 40%). El boxplot marca un valor atípico (18.8), el máximo de la muestra. Este resumen todavía mezcla las seis combinaciones posibles de molusco y concentración; el literal b) desagrega ese comportamiento.

b. Análisis Exploratorio Bivariado

Consumo de oxígeno por combinación de factores

ggplot(BD_moluscos, aes(x = c_agua, y = cons_o, fill = molusco)) +
  geom_boxplot() +
  scale_fill_manual(values = c("#622B51", "#7C92A1")) +
  labs(title = "Consumo de O2 según concentración de agua y tipo de molusco",
       x = "Concentración de agua de mar (%)",
       y = "Consumo de O2",
       fill = "Molusco") +
  theme_minimal()

medias_mol <- BD_moluscos %>%
  group_by(c_agua, molusco) %>%
  summarise(media_cons_o = round(mean(cons_o), 2), .groups = "drop")

kable(medias_mol,
      caption = "Media de consumo de O2 por combinación molusco x concentración",
      col.names = c("Concentración (%)", "Molusco", "Media Consumo O2"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Media de consumo de O2 por combinación molusco x concentración
Concentración (%)	Molusco	Media Consumo O2
100	A	9.94
100	B	7.41
75	A	7.89
75	B	6.10
50	A	12.18
50	B	12.33

En 100% y 75%, el molusco A consume en promedio más oxígeno que el B (9.94 frente a 7.41, y 7.89 frente a 6.10); en 50%, ambos quedan prácticamente igualados (12.18 frente a 12.33).

¿Se comportan igual los dos moluscos? Lectura del gráfico de interacción

ggplot(medias_mol, aes(x = c_agua, y = media_cons_o, color = molusco, group = molusco)) +
  geom_point(size = 3) +
  geom_line(linewidth = 1) +
  scale_color_manual(values = c("#401C39", "#84477C")) +
  labs(title = "Gráfico de interacción: concentración x tipo de molusco",
       x = "Concentración de agua de mar (%)",
       y = "Media de consumo de O2",
       color = "Molusco") +
  theme_minimal()

Las dos líneas dibujan la misma forma de “V” frente a la concentración; bajan de 100% a 75% y suben con fuerza en 50%, lo que ya sugiere visualmente que el patrón general no depende del tipo de molusco. Que la brecha entre A y B se cierre justo en la concentración más baja podría insinuar una interacción leve, algo que el ANOVA permite confirmar o descartar formalmente.

c. ANOVA de dos vías y Postanova

Modelo con interacción: Consumo ~ Molusco × Concentración

modelo_mol <- aov(cons_o ~ molusco * c_agua, data = BD_moluscos)

anova_tabla_mol <- as.data.frame(summary(modelo_mol)[[1]])
anova_tabla_mol <- round(anova_tabla_mol, 4)

kable(anova_tabla_mol,
      caption = "Tabla ANOVA: Consumo de O2 ~ Molusco * Concentración de agua",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Tabla ANOVA: Consumo de O2 ~ Molusco * Concentración de agua
	Df	Sum Sq	Mean Sq	F value	Pr(>F)
molusco	1	23.2269	23.2269	2.6508	0.1110
c_agua	2	230.8160	115.4080	13.1712	0.0000
molusco:c_agua	2	15.3563	7.6781	0.8763	0.4238
Residuals	42	368.0112	8.7622	NA	NA

Comprobación de supuestos del modelo

shapiro_res_mol <- shapiro.test(residuals(modelo_mol))
levene_res_mol  <- leveneTest(cons_o ~ molusco * c_agua, data = BD_moluscos)

supuestos_tabla_mol <- data.frame(
  Prueba = c("Shapiro-Wilk (normalidad)", "Levene (homogeneidad de varianzas)"),
  Estadistico = c(round(shapiro_res_mol$statistic, 4), round(levene_res_mol$`F value`[1], 4)),
  p_valor = c(round(shapiro_res_mol$p.value, 4), round(levene_res_mol$`Pr(>F)`[1], 4))
)

kable(supuestos_tabla_mol,
      caption = "Verificación de supuestos del ANOVA",
      col.names = c("Prueba", "Estadístico", "p-valor"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Verificación de supuestos del ANOVA
Prueba	Estadístico	p-valor
Shapiro-Wilk (normalidad)	0.9582	0.0857
Levene (homogeneidad de varianzas)	0.1723	0.9715

Shapiro-Wilk (p = 0.086) y Levene (p = 0.972) están ambos por encima de 0.05: ninguno de los dos supuestos se rechaza, así que el modelo es válido tal como fue planteado.

Lectura de efectos principales y de la interacción

El tipo de molusco no resulta significativo (F(1,42) = 2.65, p = 0.111), y tampoco la interacción molusco × concentración (F(2,42) = 0.88, p = 0.424): el efecto de la concentración es, entonces, el mismo para los dos moluscos, pese a la diferencia visual del literal b). La concentración de agua sí tiene un efecto muy marcado (F(2,42) = 13.17, p < 0.001), y es el único factor sobre el que tiene sentido hacer comparaciones post-hoc.

Comparación entre concentraciones (LSD)

posanova_mol <- LSD.test(modelo_mol, "c_agua", console = FALSE)

kable(as.data.frame(posanova_mol$statistics),
      caption = "Estadísticos generales del modelo",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Estadísticos generales del modelo
MSerror	Df	Mean	CV	t.value	LSD
8.76217	42	9.304792	31.8126	2.018082	2.112028

grupos_tabla_mol <- posanova_mol$groups
grupos_tabla_mol$Concentracion <- rownames(grupos_tabla_mol)
grupos_tabla_mol <- grupos_tabla_mol[, c("Concentracion", "cons_o", "groups")]

kable(grupos_tabla_mol,
      caption = "Agrupamiento LSD — Concentración de agua de mar",
      col.names = c("Concentración (%)", "Consumo O2 promedio", "Grupo"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Agrupamiento LSD — Concentración de agua de mar
Concentración (%)	Consumo O2 promedio	Grupo
50	12.25062	a
100	8.67125	b
75	6.99250	b

50% queda en un grupo aparte (“a”), con el consumo más alto (12.25); 100% (8.67) y 75% (6.99) comparten grupo (“b”), sin diferencia entre sí. La caída de 100% a 75% no llega a ser significativa, pero el salto al llegar a 50% sí lo es. Es un indicio de que una dilución fuerte del agua de mar exige a los moluscos, sin importar el tipo, un consumo de oxígeno notoriamente mayor.

Cierre del análisis

Univariado: diseño balanceado (24 obs. por molusco, 16 por concentración). Consumo de O2 con media 9.30, dispersión considerable (CV ≈ 40%) y un valor atípico (18.8).
Bivariado: patrón en “V” del consumo frente a la concentración, presente en ambos moluscos; la brecha entre A y B es clara en 100% y 75%, pero desaparece en 50%.
ANOVA: solo la concentración de agua resulta significativa (F(2,42) = 13.17, p < 0.001); molusco (p = 0.111) e interacción (p = 0.424) no lo son. Ambos supuestos del modelo se cumplen. LSD separa 50% de las otras dos concentraciones.
Lectura general: la concentración de agua de mar es el factor que gobierna el consumo de oxígeno, con un efecto igual para los dos tipos de molusco; diluir el agua hasta 50% dispara el consumo, coherente con un mayor esfuerzo osmorregulador ante el estrés hipoosmótico.

Punto 3 – Biodiversidad

Con el fin de evaluar el efecto del uso del suelo sobre la biodiversidad de anfibios en una reserva forestal, se establecieron parcelas de muestreo 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 una de las 52 parcelas (13 por hábitat) se registró la riqueza de especies, el índice de diversidad de Shannon-Wiener y la altitud (m.s.n.m.) del sitio de muestreo.

52Parcelas

4Tipos de hábitat

13Parcelas por hábitat

3Variables registradas

a. Análisis Exploratorio Univariado

Descripción general y dispersión de las variables

load("Biodiversidad.RData")

tabla_describe_bio <- describe(BD_biodiversidad[, c("Riqueza", "Shannon", "Altitud")])

kable(round(tabla_describe_bio, 2),
      caption = "Estadísticos descriptivos",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE) %>%
  scroll_box(width = "100%")

Estadísticos descriptivos
	vars	n	mean	sd	median	trimmed	mad	min	max	range	skew	kurtosis	se
Riqueza	1	52	9.73	4.37	9.50	9.60	5.19	2.00	21.00	19.00	0.28	-0.72	0.61
Shannon	2	52	1.76	0.61	1.85	1.76	0.64	0.67	2.98	2.31	0.04	-0.84	0.08
Altitud	3	52	1243.89	172.71	1254.80	1240.59	215.64	950.90	1600.80	649.90	0.09	-1.26	23.95

cv_tabla_bio <- data.frame(
  Variable = c("Riqueza", "Shannon", "Altitud"),
  CV = round(sapply(BD_biodiversidad[, c("Riqueza", "Shannon", "Altitud")], cv), 2)
)

kable(cv_tabla_bio,
      caption = "Coeficiente de variación (%)",
      col.names = c("Variable", "CV (%)"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Coeficiente de variación (%)
Variable	CV (%)
Riqueza	44.92
Shannon	34.44
Altitud	13.88

Las 52 parcelas están completas, sin datos faltantes, y repartidas en partes iguales entre los 4 hábitats (13 cada uno). En cuanto a dispersión relativa, Riqueza (CV ≈ 45%) es la más variable, seguida de Shannon (CV ≈ 34%); la Altitud (CV ≈ 14%) es comparativamente estable.

Distribución y presencia de valores atípicos

b1b <- ggplot(BD_biodiversidad, aes(y = Riqueza)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Riqueza") + theme(axis.text.x = element_blank())

b2b <- ggplot(BD_biodiversidad, aes(y = Shannon)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Shannon") + theme(axis.text.x = element_blank())

b3b <- ggplot(BD_biodiversidad, aes(y = Altitud)) +
  geom_boxplot(fill = "#D9CFE0") + theme_minimal() +
  labs(title = "Altitud") + theme(axis.text.x = element_blank())

grid.arrange(b1b, b2b, b3b, ncol = 3)

Ninguna de las tres cajas muestra puntos por fuera de los bigotes: no hay valores atípicos evidentes en ninguna variable.

Forma de las distribuciones

p1b <- ggplot(BD_biodiversidad, aes(x = Riqueza)) +
  geom_histogram(bins = 10, fill = "#622B51", color = "white") +
  theme_minimal() + labs(title = "Riqueza de especies", x = "Riqueza", y = "Frecuencia")

p2b <- ggplot(BD_biodiversidad, aes(x = Shannon)) +
  geom_histogram(bins = 10, fill = "#84477C", color = "white") +
  theme_minimal() + labs(title = "Índice de Shannon", x = "Shannon", y = "Frecuencia")

p3b <- ggplot(BD_biodiversidad, aes(x = Altitud)) +
  geom_histogram(bins = 10, fill = "#534166", color = "white") +
  theme_minimal() + labs(title = "Altitud", x = "Altitud (m.s.n.m.)", y = "Frecuencia")

grid.arrange(p1b, p2b, p3b, ncol = 3)

Riqueza y Altitud tienen una ligera cola hacia la derecha; Shannon es prácticamente simétrica. Al estar mezclados aquí los 4 hábitats, parte de esa forma responde a la combinación de grupos distintos, lo que se desagrega en el literal b).

b. Análisis Exploratorio Bivariado

Relación entre altitud y diversidad

ggplot(BD_biodiversidad, aes(x = Altitud, y = Shannon)) +
  geom_point(aes(color = Habitat), size = 2) +
  scale_color_manual(values = c("#401C39", "#622B51", "#84477C", "#7C92A1")) +
  geom_smooth(method = "lm", se = FALSE, color = "black", linetype = "dashed") +
  theme_minimal() +
  labs(title = "Shannon vs Altitud", x = "Altitud (m.s.n.m.)", y = "Índice de Shannon")

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

cor_tabla_bio <- data.frame(
  Relacion = "Altitud vs Shannon",
  r = round(cor_alt_bio$estimate, 3),
  p_valor = signif(cor_alt_bio$p.value, 3)
)

kable(cor_tabla_bio,
      caption = "Correlación de Pearson entre Altitud y Shannon",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Correlación de Pearson entre Altitud y Shannon
Relacion	r	p_valor
Altitud vs Shannon	0.752	0

La correlación es positiva y fuerte (r = 0.752, p < 0.001): a mayor altitud, mayor diversidad. Pero los puntos se agrupan por color siguiendo el mismo gradiente que la altitud, lo que sugiere que buena parte de esa relación refleja que los hábitats menos intervenidos están también en las zonas más altas de la reserva, más que un efecto altitudinal puro.

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

BD_biodiversidad$Habitat <- factor(BD_biodiversidad$Habitat,
                                    levels = c("Bosque primario", "Bosque secundario",
                                               "Sistema silvopastoril", "Potrero"))

g1b <- ggplot(BD_biodiversidad, aes(x = Habitat, y = Riqueza, fill = Habitat)) +
  geom_boxplot() + theme_minimal() +
  scale_fill_manual(values = c("#401C39", "#622B51", "#84477C", "#D9CFE0")) +
  labs(title = "Riqueza de especies por hábitat", x = "", y = "Riqueza") +
  theme(axis.text.x = element_text(angle = 20, hjust = 1), legend.position = "none")

g2b <- ggplot(BD_biodiversidad, aes(x = Habitat, y = Shannon, fill = Habitat)) +
  geom_boxplot() + theme_minimal() +
  scale_fill_manual(values = c("#401C39", "#622B51", "#84477C", "#D9CFE0")) +
  labs(title = "Índice de Shannon por hábitat", x = "", y = "Shannon") +
  theme(axis.text.x = element_text(angle = 20, hjust = 1), legend.position = "none")

grid.arrange(g1b, g2b, ncol = 2)

tabla_grupo_bio <- BD_biodiversidad %>%
  group_by(Habitat) %>%
  summarise(
    n = n(),
    Riqueza_media = round(mean(Riqueza), 2),
    Riqueza_sd = round(sd(Riqueza), 2),
    Shannon_media = round(mean(Shannon), 2),
    Shannon_sd = round(sd(Shannon), 2)
  )

kable(tabla_grupo_bio,
      caption = "Riqueza y Shannon promedio por tipo de hábitat",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Riqueza y Shannon promedio por tipo de hábitat
Habitat	n	Riqueza_media	Riqueza_sd	Shannon_media	Shannon_sd
Bosque primario	13	14.77	2.74	2.47	0.38
Bosque secundario	13	11.08	3.01	1.97	0.22
Sistema silvopastoril	13	8.08	2.50	1.60	0.28
Potrero	13	5.00	1.29	1.00	0.24

Riqueza y Shannon bajan de forma ordenada a lo largo del gradiente de intervención: máximos en bosque primario, mínimos en potrero. A mayor intervención humana, menor biodiversidad de anfibios registrada.

c. ANOVA y Postanova

Modelo: Shannon ~ Hábitat

anova_bio <- aov(Shannon ~ Habitat, data = BD_biodiversidad)

anova_tabla_bio <- as.data.frame(summary(anova_bio)[[1]])
anova_tabla_bio <- round(anova_tabla_bio, 4)

kable(anova_tabla_bio,
      caption = "Tabla ANOVA: Shannon ~ Habitat",
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Tabla ANOVA: Shannon ~ Habitat
	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

El hábitat tiene un efecto altamente significativo sobre el Shannon (F(3,48) = 60.61, p < 0.001), al menos uno de los cuatro hábitats difiere de los demás.

Supuestos del modelo

shapiro_res_bio <- shapiro.test(residuals(anova_bio))
levene_res_bio <- leveneTest(Shannon ~ Habitat, data = BD_biodiversidad)

supuestos_tabla_bio <- data.frame(
  Prueba = c("Shapiro-Wilk (normalidad)", "Levene (homogeneidad de varianzas)"),
  Estadistico = c(round(shapiro_res_bio$statistic, 4), round(levene_res_bio$`F value`[1], 4)),
  p_valor = c(round(shapiro_res_bio$p.value, 4), round(levene_res_bio$`Pr(>F)`[1], 4))
)

kable(supuestos_tabla_bio,
      caption = "Verificación de supuestos del ANOVA",
      col.names = c("Prueba", "Estadístico", "p-valor"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Verificación de supuestos del ANOVA
Prueba	Estadístico	p-valor
Shapiro-Wilk (normalidad)	0.9785	0.4625
Levene (homogeneidad de varianzas)	1.2768	0.2930

Shapiro-Wilk (p = 0.463) y Levene (p = 0.293) están ambos por encima de 0.05: ninguno de los dos supuestos se rechaza, así que el modelo y su post-hoc son confiables.

Comparaciones entre hábitats (LSD)

posanova_bio <- LSD.test(anova_bio, trt = "Habitat")

kable(as.data.frame(posanova_bio$statistics),
      caption = "Estadísticos generales del modelo",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Estadísticos generales del modelo
MSerror	Df	Mean	CV	t.value	LSD
0.0817362	48	1.762308	16.22279	2.010635	0.2254674

kable(as.data.frame(posanova_bio$parameters),
      caption = "Parámetros de la prueba LSD",
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Parámetros de la prueba LSD
test	p.ajusted	name.t	ntr	alpha
Fisher-LSD	none	Habitat	4	0.05

grupos_tabla_bio <- posanova_bio$groups
grupos_tabla_bio$Habitat <- rownames(grupos_tabla_bio)
grupos_tabla_bio <- grupos_tabla_bio[, c("Habitat", "Shannon", "groups")]

kable(grupos_tabla_bio,
      caption = "Agrupamiento LSD — Tipo de hábitat",
      col.names = c("Hábitat", "Shannon promedio", "Grupo"),
      row.names = FALSE,
      table.attr = "class='tabla-verde'") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                full_width = FALSE)

Agrupamiento LSD — Tipo de hábitat
Hábitat	Shannon promedio	Grupo
Bosque primario	2.467692	a
Bosque secundario	1.974615	b
Sistema silvopastoril	1.603077	c
Potrero	1.003846	d

Los cuatro hábitats reciben letras distintas (a, b, c, d): Bosque primario (2.47), Bosque secundario (1.97), Sistema silvopastoril (1.60) y Potrero (1.00) difieren todos entre sí, confirmando el gradiente decreciente ya sugerido en el análisis bivariado.

Conclusiones generales

Univariado: distribuciones aproximadamente simétricas y sin valores atípicos en las tres variables; Riqueza y Shannon son bastante más dispersas que la Altitud.
Bivariado: Riqueza y Shannon bajan de forma ordenada según el gradiente de intervención; la Altitud correlaciona fuerte con Shannon (r = 0.752), aunque probablemente confundida con el propio hábitat.
ANOVA: el hábitat genera diferencias altamente significativas en el Shannon (F(3,48) = 60.61, p < 0.001), con ambos supuestos cumplidos. LSD separa los cuatro hábitats entre sí, con diversidad decreciente: Bosque primario > Bosque secundario > Sistema silvopastoril > Potrero.
Lectura general: la intervención antrópica del suelo tiene un efecto negativo, claro y estadísticamente robusto sobre la biodiversidad de anfibios, lo que respalda priorizar la conservación de bosque primario y limitar la conversión hacia sistemas silvopastoriles o potreros.

TALLER 1

Exploratorio, ANOVA y pruebas post-ANOVA

Elizabeth Rebolledo (202416732) y Maria Paula Giraldo (202418989)

25 de Septiembre del 2026

Punto 1 – Salinidad

a. Análisis Exploratorio Univariado

Panorama general de las variables

¿Qué tan dispersas están las variables?

Formas de distribución y presencia de atípicos

b. Análisis Exploratorio Bivariado

Matriz de correlaciones

La covariable más relacionada con la biomasa

Comprobación visual: diagramas de dispersión

c. ANOVA y Postanova

Construcción de los niveles de pH

Ajuste del modelo ANOVA

¿Es válido el modelo? Verificación de supuestos

Diferencias entre niveles: prueba LSD

Síntesis y hallazgos principales

Punto 2 – Moluscos

a. Análisis Exploratorio Univariado

Balance del diseño experimental

El consumo de oxígeno en cifras

b. Análisis Exploratorio Bivariado

Consumo de oxígeno por combinación de factores

¿Se comportan igual los dos moluscos? Lectura del gráfico de interacción

c. ANOVA de dos vías y Postanova

Modelo con interacción: Consumo ~ Molusco × Concentración

Comprobación de supuestos del modelo

Lectura de efectos principales y de la interacción

Comparación entre concentraciones (LSD)

Cierre del análisis

Punto 3 – Biodiversidad

a. Análisis Exploratorio Univariado

Descripción general y dispersión de las variables

Distribución y presencia de valores atípicos

Forma de las distribuciones

b. Análisis Exploratorio Bivariado

Relación entre altitud y diversidad

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

c. ANOVA y Postanova

Modelo: Shannon ~ Hábitat

Supuestos del modelo

Comparaciones entre hábitats (LSD)

Conclusiones generales