datos$abund_total <- rowSums(especies)
ggplot(datos, aes(x=Parcela, y=abund_total)) +
geom_boxplot(fill="skyblue") +
theme_minimal() +
labs(title="Abundancia total por parcela",
x="Parcela",
y="Total de individuos")

# Fase 2 — Diversidad alfa por parcela
library(ggplot2)
# asegurar que las abundancias sean numéricas
especies <- data.frame(lapply(especies, as.numeric))
# riqueza de especies
riqueza <- rowSums(especies > 0)
# total de individuos por parcela
N <- rowSums(especies)
# proporciones
p <- especies / N
# índice de Shannon
shannon <- -rowSums(p * log(p), na.rm = TRUE)
# índice de Simpson
simpson <- 1 - rowSums(p^2, na.rm = TRUE)
# tabla de resultados
div_alpha <- data.frame(riqueza, shannon, simpson)
# gráfico
ggplot(div_alpha, aes(x = riqueza, y = shannon)) +
geom_point(color = "black", size = 3) +
theme_minimal() +
labs(title = "Diversidad alfa",
x = "Riqueza de especies",
y = "Índice de Shannon")

# Fase 3 — Diversidad beta (Bray-Curtis + NMDS con símbolos por etapa)
library(vegan)
# asegurar que especies sean numéricas
especies <- data.frame(lapply(especies, as.numeric))
# distancia Bray-Curtis
dist_bc <- vegdist(especies, method = "bray")
# NMDS
nmds <- metaMDS(dist_bc, k = 2, trymax = 100)
## Run 0 stress 0.2359893
## Run 1 stress 0.2416836
## Run 2 stress 0.24089
## Run 3 stress 0.2401874
## Run 4 stress 0.2456725
## Run 5 stress 0.2384584
## Run 6 stress 0.2420171
## Run 7 stress 0.2497067
## Run 8 stress 0.2423942
## Run 9 stress 0.244837
## Run 10 stress 0.2435351
## Run 11 stress 0.2396573
## Run 12 stress 0.2407004
## Run 13 stress 0.2370278
## Run 14 stress 0.2465409
## Run 15 stress 0.2451368
## Run 16 stress 0.2405292
## Run 17 stress 0.2456837
## Run 18 stress 0.2461727
## Run 19 stress 0.2396027
## Run 20 stress 0.2454913
## Run 21 stress 0.2472796
## Run 22 stress 0.2474893
## Run 23 stress 0.2402888
## Run 24 stress 0.2426131
## Run 25 stress 0.2532829
## Run 26 stress 0.2393044
## Run 27 stress 0.2443609
## Run 28 stress 0.2451614
## Run 29 stress 0.2442635
## Run 30 stress 0.2392652
## Run 31 stress 0.2459351
## Run 32 stress 0.240987
## Run 33 stress 0.2412861
## Run 34 stress 0.2442569
## Run 35 stress 0.2476911
## Run 36 stress 0.2449062
## Run 37 stress 0.2387351
## Run 38 stress 0.2419828
## Run 39 stress 0.2473441
## Run 40 stress 0.2417272
## Run 41 stress 0.2454807
## Run 42 stress 0.2407297
## Run 43 stress 0.2440366
## Run 44 stress 0.2418088
## Run 45 stress 0.2450816
## Run 46 stress 0.2411047
## Run 47 stress 0.2404376
## Run 48 stress 0.2417591
## Run 49 stress 0.2441451
## Run 50 stress 0.2439626
## Run 51 stress 0.2425817
## Run 52 stress 0.2410538
## Run 53 stress 0.2420953
## Run 54 stress 0.2431994
## Run 55 stress 0.243366
## Run 56 stress 0.2487255
## Run 57 stress 0.2406456
## Run 58 stress 0.236798
## Run 59 stress 0.2436913
## Run 60 stress 0.239722
## Run 61 stress 0.237142
## Run 62 stress 0.2505865
## Run 63 stress 0.2389298
## Run 64 stress 0.2485323
## Run 65 stress 0.2443844
## Run 66 stress 0.2446347
## Run 67 stress 0.2401312
## Run 68 stress 0.2442051
## Run 69 stress 0.2437379
## Run 70 stress 0.2463793
## Run 71 stress 0.2367978
## Run 72 stress 0.244891
## Run 73 stress 0.2457738
## Run 74 stress 0.2467106
## Run 75 stress 0.2464753
## Run 76 stress 0.2514881
## Run 77 stress 0.2390793
## Run 78 stress 0.2441824
## Run 79 stress 0.2418518
## Run 80 stress 0.2408587
## Run 81 stress 0.2427618
## Run 82 stress 0.243522
## Run 83 stress 0.2427404
## Run 84 stress 0.241556
## Run 85 stress 0.2446409
## Run 86 stress 0.2378338
## Run 87 stress 0.2426781
## Run 88 stress 0.2478905
## Run 89 stress 0.2452677
## Run 90 stress 0.2434193
## Run 91 stress 0.2477175
## Run 92 stress 0.2390524
## Run 93 stress 0.244438
## Run 94 stress 0.240365
## Run 95 stress 0.2460099
## Run 96 stress 0.2469245
## Run 97 stress 0.2399522
## Run 98 stress 0.2482688
## Run 99 stress 0.237515
## Run 100 stress 0.2405812
## *** Best solution was not repeated -- monoMDS stopping criteria:
## 21: no. of iterations >= maxit
## 79: stress ratio > sratmax
# coordenadas
coord <- as.data.frame(scores(nmds))
# traer etapa y normalizar texto
coord$Etapa <- tolower(trimws(as.character(datos$Etapa)))
# crear vector de símbolos (por defecto círculo)
pch_vec <- rep(16, nrow(coord)) # ● pionera
pch_vec[coord$Etapa == "intermedia"] <- 17 # ▲
pch_vec[coord$Etapa == "tardia"] <- 15 # ■
pch_vec[coord$Etapa == "tardía"] <- 15 # ■ (con acento)
pch_vec[coord$Etapa == "madura"] <- 3 # +
# gráfico
plot(coord$NMDS1, coord$NMDS2,
pch = pch_vec,
col = "darkblue",
xlab = "NMDS1",
ylab = "NMDS2",
main = "NMDS Bray-Curtis por etapa sucesional")
# leyenda
legend("topright",
legend = c("Pionera","Intermedia","Tardía","Madura"),
pch = c(16,17,15,3),
col = "darkblue")

# Fase 4 — Curvas rango-abundancia por etapa sucesional
# asegurar que especies sean numéricas
especies <- data.frame(lapply(especies, as.numeric))
# obtener etapas únicas
etapas <- unique(datos$Etapa)
# hacer gráfico para cada etapa
for(e in etapas){
# filtrar datos de la etapa
sub_especies <- especies[datos$Etapa == e, ]
# abundancia total por especie
abund <- colSums(sub_especies)
# ordenar abundancias
rank_abund <- sort(abund, decreasing = TRUE)
# gráfico
plot(rank_abund,
type = "b",
pch = 19,
col = "black",
xlab = "Rango de especies",
ylab = "Abundancia",
main = paste("Curva rango-abundancia -", e))
}




# Fase 5 — Red de interacciones ecológicas con todas las especies
library(igraph)
# asegurar que especies sean numéricas
especies <- data.frame(lapply(especies, as.numeric))
# calcular abundancia total por especie
abund_especies <- colSums(especies)
# crear matriz de interacción entre todas las especies
matriz_red <- outer(abund_especies, abund_especies, FUN = function(x,y) (x+y)/2)
# convertir matriz en red
red <- graph_from_adjacency_matrix(matriz_red,
mode = "undirected",
weighted = TRUE,
diag = FALSE)
# graficar red
plot(red,
vertex.size = 15,
vertex.color = "orange",
vertex.label = names(abund_especies),
vertex.label.color = "black",
main = "Red ecológica de especies")

# Fase 6 — Pirámides de energía por etapa sucesional
par(mfrow=c(2,2)) # 4 gráficos en una sola figura
niveles <- c("Productores","Herbivoros","Carnivoros","Depredadores")
# Pionera
energia <- c(12000,1200,120,12)
bp <- barplot(energia,
names.arg=niveles,
col=c("darkgreen","yellowgreen","orange","red"),
main="Pionera",
ylab="Energia",
ylim=c(0,13000))
text(bp, energia+200, energia)
# Intermedia
energia <- c(10000,1000,100,10)
bp <- barplot(energia,
names.arg=niveles,
col=c("darkgreen","yellowgreen","orange","red"),
main="Intermedia",
ylab="Energia",
ylim=c(0,11000))
text(bp, energia+200, energia)
# Tardia
energia <- c(8000,800,80,8)
bp <- barplot(energia,
names.arg=niveles,
col=c("darkgreen","yellowgreen","orange","red"),
main="Tardia",
ylab="Energia",
ylim=c(0,9000))
text(bp, energia+200, energia)
# Madura
energia <- c(6000,600,60,6)
bp <- barplot(energia,
names.arg=niveles,
col=c("darkgreen","yellowgreen","orange","red"),
main="Madura",
ylab="Energia",
ylim=c(0,7000))
text(bp, energia+200, energia)

# Fase 7 — Relación N:P por etapa (igual estilo que Fase 1)
# etapas
etapas <- c("Pionera","Intermedia","Tardia","Madura")
# valores de ejemplo de relación N:P
relacion_np <- c(12, 10, 8, 6)
# gráfico tipo Fase 1
barplot(relacion_np,
names.arg = etapas,
col = "forestgreen",
main = "Relación N:P por etapa sucesional",
ylab = "Relación N:P")
