PREGUNTO A CHATGPT: ¿¿¿Al calcular los density plots, se podía de algún modo ponderar por la abundancia relativa de la especie en lugar de solamente por la presencia???

Primeros pasos

Excelente pregunta 🌿 — y muy avanzada, porque justamente uno de los supuestos centrales del método de Broennimann et al. (2012) (el PCA-env que usa ecospat) es que trabaja con presencias, no con abundancias.
Pero sí hay formas de ponderar las densidades por abundancia o probabilidad de ocurrencia. Vamos a desglosarlo 👇

🔹 Contexto: cómo `ecospat` calcula las densidades

El paso clave es:

z.sp <- ecospat.grid.clim.dyn(
  glob = scores.bkg.env,
  glob1 = scores.bkg.env,
  sp = scores.occ.sp,
  R = 100
)

Internamente, esta función crea un mapa de densidad suavizada (kernel density) en el espacio PCA de clima, usando las coordenadas (scores.occ.sp) de presencia.

👉 Cada punto tiene peso 1, es decir, todas las ocurrencias se consideran igual de importantes.

🔹 Qué pasa si se tienen abundancias o probabilidades de presencia

Si tus datos incluyen una abundancia relativa (o algún tipo de peso de importancia, como biomasa, cobertura o esfuerzo de muestreo), se puede incorporar como un peso en el cálculo de la densidad fuera del flujo directo de ecospat.

El paquete ecospat no tiene un argumento explícito de “weights”, pero se puede adaptar el flujo así:

💡 Opción 1 — Usar un kernel density ponderado manualmente

Se puede calcular tu propio mapa de densidad en el espacio PCA usando, por ejemplo:

library(MASS)

# Supongamos que se tiene:
# scores.occ.sp  → coordenadas PCA (2 ejes)
# weights.sp     → abundancia relativa o peso de cada punto

dens.abund <- kde2d(
  x = scores.occ.sp[,1],
  y = scores.occ.sp[,2],
  w = weights.sp,   # <--- pondera por abundancia
  n = 100
)

Luego se puede escalar este mapa y usarlo para construir un objeto “tipo ecospat” (como z.sp) o simplemente comparar visualmente con otros nichos.

💡 Opción 2 — Expandir los puntos proporcionalmente a su abundancia

Si las abundancias son discretas (p. ej. conteos enteros), se puede expandir el dataset de ocurrencias repitiendo cada punto según su abundancia:

occ.expanded <- scores.occ.sp[rep(1:nrow(scores.occ.sp), times = round(weights.sp)), ]

Y luego usarlo normalmente en:

z.sp <- ecospat.grid.clim.dyn(glob = scores.bkg.env, sp = occ.expanded, R = 100)

Esto hace que los puntos más abundantes aporten más densidad a la superficie.

💡 Opción 3 — Si tus abundancias son probabilidades (0–1)

Puedes usar esas probabilidades como pesos al hacer el kernel density (como en Opción 1).
Esto sería conceptualmente equivalente a un nicho de probabilidad de presencia ponderado.

🔹 Limitación importante

⚠️ Cuidado con la interpretación:
Los tests de ecospat (equivalency y similarity) asumen que las densidades reflejan probabilidades de uso del espacio ambiental, no abundancia poblacional real. Si ponderás por abundancia, ya no estás modelando el “nicho realizado” estándar, sino algo más cercano a la intensidad de uso del espacio ambiental (ver Hirzel & Le Lay 2008, Broennimann et al. 2012 discuten esto).

✅ En resumen

Método	Se puede hacer	Cómo	Comentario
`ecospat.grid.clim.dyn` estándar	✅	Solo presencias (peso 1)	Modelo original
Ponderar por abundancia	⚙️ Manual	`kde2d` o expandir puntos	Cambia la interpretación
Ponderar por probabilidad	⚙️ Manual	Kernel density ponderado	Representa uso relativo del espacio ambiental

¿Querés que te muestre un ejemplo en R de cómo implementar el kernel density ponderado (Opción 1) con tus coordenadas PCA y pesos de abundancia?

Ejemplo de Kernel Density Ponderado

Perfecto 🙌 Te muestro un ejemplo completo y comentado en R para que veas cómo calcular un mapa de densidad ambiental ponderado por abundancia relativa, equivalente al que usa ecospat.grid.clim.dyn, pero incluyendo tus pesos.

Supongamos que ya tienes:

scores.occ.sp: coordenadas de tus ocurrencias en el espacio PCA (2 ejes).
weights.sp: un vector numérico con abundancia relativa o peso (de la misma longitud que el número de ocurrencias).

🧩 Paso 1 — Datos de ejemplo

# Simulemos datos tipo PCA y abundancias
set.seed(42)
scores.occ.sp <- data.frame(
  PC1 = rnorm(200, 0, 1),
  PC2 = rnorm(200, 0, 1)
)
weights.sp <- runif(200, 0.1, 5)  # abundancia relativa o pesos

🧮 Paso 2 — Calcular densidad ponderada con `kde2d.weighted`

MASS::kde2d() no tiene argumento de peso por defecto, así que usamos una función personalizada:

# Función para kernel density 2D con pesos
kde2d.weighted <- function(x, y, w = NULL, n = 100, lims = c(range(x), range(y))) {
  if (is.null(w)) w <- rep(1, length(x))
  w <- w / sum(w)  # normalizar los pesos
  gx <- seq(lims[1], lims[2], length.out = n)
  gy <- seq(lims[3], lims[4], length.out = n)
  h <- c(bandwidth.nrd(x), bandwidth.nrd(y))
  z <- outer(
    gx, gy,
    Vectorize(function(i, j) {
      sum(w * dnorm((i - x) / h[1]) * dnorm((j - y) / h[2]))
    })
  )
  z <- z / sum(z)  # normalizar la densidad
  list(x = gx, y = gy, z = z)
}

🌍 Paso 3 — Calcular la densidad ponderada

dens.abund <- kde2d.weighted(
  x = scores.occ.sp$PC1,
  y = scores.occ.sp$PC2,
  w = weights.sp,
  n = 100
)

🎨 Paso 4 — Visualizar la densidad

library(ggplot2)

dens_df <- expand.grid(PC1 = dens.abund$x, PC2 = dens.abund$y)
dens_df$z <- as.vector(dens.abund$z)

ggplot(dens_df, aes(PC1, PC2, fill = z)) +
  geom_raster() +
  scale_fill_viridis_c(option = "plasma") +
  coord_equal() +
  labs(
    title = "Densidad ponderada por abundancia en el espacio PCA",
    fill = "Densidad"
  ) +
  theme_minimal()

🔹 Interpretación:
Las zonas más “calientes” (valores altos de z) son las condiciones ambientales donde la especie tiene mayor abundancia relativa, no solo presencia.

💡 Opcional: comparar con densidad sin ponderar

dens.unweighted <- kde2d.weighted(scores.occ.sp$PC1, scores.occ.sp$PC2, n = 100)
mean(abs(dens.abund$z - dens.unweighted$z))

Esto te muestra cuánto cambia la forma del nicho al ponderar por abundancia.

¿Querés que te muestre cómo podrías usar este mapa de densidad ponderado dentro del flujo de ecospat (por ejemplo, reemplazando el z$z.occ de un objeto ecospat.grid.clim.dyn)?

Integrar mapa ponderado en flujo de ecospat

Perfecto 🌿 — te muestro cómo integrar tu densidad ponderada por abundancia dentro del flujo de trabajo de ecospat, de modo que puedas seguir usando las funciones de comparación de nichos (ecospat.niche.overlap, ecospat.niche.equivalency.test, etc.) sin romper su estructura.

🧭 Objetivo

Queremos reemplazar el componente de densidad estándar (z$z.occ) de un objeto ecospat.grid.clim.dyn() con nuestra versión ponderada por abundancia.

🔹 Paso 1 — Crear el nicho estándar en el espacio PCA

(Usando ecospat.grid.clim.dyn con presencias normales)

z.sp <- ecospat.grid.clim.dyn(
  glob = scores.bkg.env,       # Fondo ambiental (scores del PCA)
  glob1 = scores.bkg.env,      # Mismo fondo, salvo que uses subsets
  sp = scores.occ.sp,          # Ocurrencias (sin ponderar)
  R = 100                      # Resolución de la grilla
)

Esto te genera algo así como:

str(z.sp)
# List of 8
# $ z.occ : matriz 100x100 con la densidad suavizada (presencias)
# $ z.glob: matriz 100x100 con la densidad del fondo
# $ x, $y : coordenadas del grid PCA
# etc.

🔹 Paso 2 — Calcular tu densidad ponderada (como antes)

Usamos la función personalizada kde2d.weighted() que definimos antes 👇

dens.abund <- kde2d.weighted(
  x = scores.occ.sp$PC1,
  y = scores.occ.sp$PC2,
  w = weights.sp,   # tus abundancias o pesos relativos
  n = 100,
  lims = c(range(z.sp$x), range(z.sp$y))
)

🔹 Paso 3 — Reemplazar la densidad en el objeto `z.sp`

# Reemplazamos la densidad de ocurrencias por la ponderada
z.sp$z.occ <- dens.abund$z / sum(dens.abund$z)  # normalizamos para mantener escala

Puedes comprobarlo visualmente:

par(mfrow = c(1, 2))
image(z.sp$x, z.sp$y, z.sp$z.occ, main = "Densidad ponderada")
image(z.sp$x, z.sp$y, z.sp$z.glob, main = "Fondo ambiental")

🔹 Paso 4 — Usarlo en los tests de `ecospat`

Ahora z.sp tiene exactamente la misma estructura que un objeto creado con ecospat.grid.clim.dyn, pero con tu densidad ajustada.

Se puede usar sin problema en:

# Comparación de solapamiento
ecospat.niche.overlap(z.sp, z.other, cor = TRUE)

# Test de equivalencia
ecospat.niche.equivalency.test(z.sp, z.other, rep = 100)

# Test de similitud
ecospat.niche.similarity.test(z.sp, z.other, rep = 100)

⚠️ Consideraciones importantes

Aspecto	Implicación
Interpretación	Ahora el “nicho” refleja la intensidad de uso ambiental ponderada por abundancia, no solo presencia.
Comparabilidad	Si se compara con otra especie/grupo, conviene que ambos nichos estén construidos bajo el mismo criterio (también ponderados o ambos sin ponderar).
Escala	Si tus pesos no están normalizados, conviene hacerlo (`w <- w / sum(w)`) antes del kernel.
Reproducibilidad	Guardá la función `kde2d.weighted()` en tu script para asegurar consistencia.

✅ En resumen

Con este método:

Se conserva toda la estructura ecospat (por compatibilidad con sus funciones).
Se puede incorporar abundancia relativa o probabilidad de presencia.
La interpretación pasa de “nicho ocupado” → a “nicho ponderado por intensidad de uso”.

¿Querés que te deje una versión del código completa en bloque (lista para copiar y ejecutar con tus datos reales)?

Código completo para nicho ponderado

Perfecto 👌 — aquí tienes un bloque completo de código en R, listo para usar con tus datos.
Este flujo te permite generar nichos ecológicos ponderados por abundancia e integrarlos al marco de trabajo de ecospat para hacer comparaciones de nicho (solapamiento, equivalencia, similitud, etc.).

🧩 Código completo para nicho ponderado por abundancia

# -----------------------------------------------------------
# 1️⃣ Cargar librerías
# -----------------------------------------------------------
library(ecospat)
library(MASS)   # para kde2d
library(dplyr)

# -----------------------------------------------------------
# 2️⃣ Función para estimar densidad 2D ponderada
# -----------------------------------------------------------
kde2d.weighted <- function(x, y, w = NULL, n = 100, lims = c(range(x), range(y))) {
  if (is.null(w)) w <- rep(1, length(x))
  w <- w / sum(w)  # normalizamos pesos
  
  # Grilla de densidad
  gx <- seq(lims[1], lims[2], length.out = n)
  gy <- seq(lims[3], lims[4], length.out = n)
  
  # Matriz de densidad ponderada
  z <- matrix(0, n, n)
  for (i in seq_along(x)) {
    kx <- dnorm(gx, mean = x[i], sd = diff(range(x)) / n)
    ky <- dnorm(gy, mean = y[i], sd = diff(range(y)) / n)
    z <- z + w[i] * outer(kx, ky)
  }
  
  list(x = gx, y = gy, z = z / sum(z))
}

# -----------------------------------------------------------
# 3️⃣ Crear PCA ambiental (o usar el existente)
# -----------------------------------------------------------

# Suponiendo que ya tenés tu PCA global:
# pca.env.global <- dudi.pca(environmental_data, scannf = FALSE, nf = 2)
# Y tus puntuaciones globales (fondo ambiental):
scores.bkg.env <- pca.env.global$li

# -----------------------------------------------------------
# 4️⃣ Ocurrencias y abundancias de tu especie
# -----------------------------------------------------------
# Ejemplo: tus ocurrencias en coordenadas del PCA
# (las mismas variables que usaste para la PCA)
scores.occ.sp <- scores.bkg.env[match(occ.sp$ID, rownames(scores.bkg.env)), ]

# Vector de abundancias relativas (normalizado)
weights.sp <- occ.sp$abundance / sum(occ.sp$abundance)

# -----------------------------------------------------------
# 5️⃣ Crear el objeto de nicho estándar (solo presencias)
# -----------------------------------------------------------
z.sp <- ecospat.grid.clim.dyn(
  glob = scores.bkg.env,
  glob1 = scores.bkg.env,
  sp = scores.occ.sp,
  R = 100
)

# -----------------------------------------------------------
# 6️⃣ Calcular densidad ponderada por abundancia
# -----------------------------------------------------------
dens.abund <- kde2d.weighted(
  x = scores.occ.sp$Axis1,  # o PC1
  y = scores.occ.sp$Axis2,  # o PC2
  w = weights.sp,
  n = 100,
  lims = c(range(z.sp$x), range(z.sp$y))
)

# -----------------------------------------------------------
# 7️⃣ Reemplazar densidad en el objeto ecospat
# -----------------------------------------------------------
z.sp$z.occ <- dens.abund$z / sum(dens.abund$z)

# -----------------------------------------------------------
# 8️⃣ Visualizar para comprobar
# -----------------------------------------------------------
par(mfrow = c(1, 2))
image(z.sp$x, z.sp$y, z.sp$z.glob, main = "Fondo ambiental", col = topo.colors(100))
image(z.sp$x, z.sp$y, z.sp$z.occ, main = "Densidad ponderada por abundancia", col = heat.colors(100))

# -----------------------------------------------------------
# 9️⃣ Usar el nicho ponderado en los análisis
# -----------------------------------------------------------
# Ejemplo: comparar con otra especie u origen
# z.other <- ecospat.grid.clim.dyn(...)

ecospat.niche.overlap(z.sp, z.other, cor = TRUE)
ecospat.niche.equivalency.test(z.sp, z.other, rep = 100)
ecospat.niche.similarity.test(z.sp, z.other, rep = 100)

📘 Notas importantes

weights.sp puede ser abundancia, frecuencia relativa, o probabilidad de presencia (si proviene de un modelo).
Si se tienen muchos ceros o abundancias muy sesgadas, se puede transformar antes (por ejemplo, log1p(weights.sp) o raíz cuadrada).
Si se quiere mantener la estructura comparable, hay que aplicar el mismo proceso para cada grupo o especie.

¿Querés que te deje una versión del código adaptada a tus objetos específicos (por ejemplo, z.NZ, z.Exotic, scores.bkg.env, etc.) para que puedas copiar y ejecutar directamente sobre tu análisis actual?

Código adaptado a mi caso

Perfecto 🌱 — te dejo el código completamente adaptado a tu contexto, usando tu estructura de trabajo (por lo que mencionaste antes: z.NZ, z.Exotic, scores.bkg.env, PCA-ocu, etc.).
Esto te va a permitir reemplazar los nichos construidos con presencia/ausencia por nichos ponderados por abundancia relativa, manteniendo compatibilidad total con las funciones de ecospat.

🧩 Versión adaptada a tu análisis: PCA-Ocu (Broennimann et al.)

# -----------------------------------------------------------
# 1️⃣ Cargar librerías
# -----------------------------------------------------------
library(ecospat)
library(MASS)

# -----------------------------------------------------------
# 2️⃣ Función para densidad 2D ponderada por abundancia
# -----------------------------------------------------------
kde2d.weighted <- function(x, y, w = NULL, n = 100, lims = c(range(x), range(y))) {
  if (is.null(w)) w <- rep(1, length(x))
  w <- w / sum(w)  # normalizamos pesos
  
  gx <- seq(lims[1], lims[2], length.out = n)
  gy <- seq(lims[3], lims[4], length.out = n)
  
  z <- matrix(0, n, n)
  for (i in seq_along(x)) {
    kx <- dnorm(gx, mean = x[i], sd = diff(range(x)) / n)
    ky <- dnorm(gy, mean = y[i], sd = diff(range(y)) / n)
    z <- z + w[i] * outer(kx, ky)
  }
  list(x = gx, y = gy, z = z / sum(z))
}

# -----------------------------------------------------------
# 3️⃣ Inputs: PCA y datos
# -----------------------------------------------------------
# Tu PCA ya creada con ocurrencias (PCA-Ocu):
# pca.env.global <- dudi.pca(environmental_data, scannf = FALSE, nf = 2)

# Puntuaciones globales (de fondo ambiental)
scores.bkg.env <- pca.env.global$li

# Puntuaciones de ocurrencias por grupo:
# (ya proyectadas sobre el PCA)
scores.occ.NZ     <- scores.bkg.env[match(occ.NZ$ID, rownames(scores.bkg.env)), ]
scores.occ.Exotic <- scores.bkg.env[match(occ.Exotic$ID, rownames(scores.bkg.env)), ]

# Pesos o abundancias relativas
weights.NZ     <- occ.NZ$abund / sum(occ.NZ$abund)
weights.Exotic <- occ.Exotic$abund / sum(occ.Exotic$abund)

# -----------------------------------------------------------
# 4️⃣ Crear nichos base (no ponderados)
# -----------------------------------------------------------
z.NZ <- ecospat.grid.clim.dyn(
  glob  = scores.bkg.env,
  glob1 = scores.bkg.env,
  sp    = scores.occ.NZ,
  R     = 100
)

z.Exotic <- ecospat.grid.clim.dyn(
  glob  = scores.bkg.env,
  glob1 = scores.bkg.env,
  sp    = scores.occ.Exotic,
  R     = 100
)

# -----------------------------------------------------------
# 5️⃣ Calcular densidades ponderadas
# -----------------------------------------------------------
dens.NZ <- kde2d.weighted(
  x = scores.occ.NZ$Axis1,   # o PC1 según tu PCA
  y = scores.occ.NZ$Axis2,   # o PC2
  w = weights.NZ,
  n = 100,
  lims = c(range(z.NZ$x), range(z.NZ$y))
)

dens.Exotic <- kde2d.weighted(
  x = scores.occ.Exotic$Axis1,
  y = scores.occ.Exotic$Axis2,
  w = weights.Exotic,
  n = 100,
  lims = c(range(z.Exotic$x), range(z.Exotic$y))
)

# -----------------------------------------------------------
# 6️⃣ Reemplazar densidades originales por ponderadas
# -----------------------------------------------------------
z.NZ$z.occ     <- dens.NZ$z / sum(dens.NZ$z)
z.Exotic$z.occ <- dens.Exotic$z / sum(dens.Exotic$z)

# -----------------------------------------------------------
# 7️⃣ Visualizar los nichos
# -----------------------------------------------------------
par(mfrow = c(1, 2))
image(z.NZ$x, z.NZ$y, z.NZ$z.occ, main = "Nicho ponderado - Nativo (NZ)", col = heat.colors(100))
image(z.Exotic$x, z.Exotic$y, z.Exotic$z.occ, main = "Nicho ponderado - Exótico", col = heat.colors(100))

# -----------------------------------------------------------
# 8️⃣ Comparaciones de nicho ponderado
# -----------------------------------------------------------
# Superposición
ecospat.niche.overlap(z.NZ, z.Exotic, cor = TRUE)

# Test de equivalencia
equiv.abund <- ecospat.niche.equivalency.test(z.NZ, z.Exotic, rep = 100)
equiv.abund

# Test de similitud (orden importa)
sim.abund <- ecospat.niche.similarity.test(z.NZ, z.Exotic, rep = 100, rand.type = 1)
sim.abund

# -----------------------------------------------------------
# 9️⃣ Interpretación visual rápida
# -----------------------------------------------------------
hist(sim.abund$sim$D, col = "lightblue",
     main = "Distribución permutada de D (Similitud ponderada)",
     xlab = "Índice de solapamiento D")
abline(v = sim.abund$obs$D, col = "red", lwd = 2)

💬 Cómo interpretarlo

Componente	Qué mide	Interpretación si ponderás por abundancia
`D`, `I`	Superposición del nicho ponderado	Mide cuánto se solapan los ambientes intensamente utilizados por cada grupo.
`expansion`, `unfilling`, `stability`	Dinámica entre nichos	Reflejan proporciones del nicho ambiental ocupadas con distinta intensidad.
`p.D`, `p.I`	Significancia (tests)	Indican si el solapamiento observado es mayor o menor que el esperado por azar.

⚠️ Tips finales

Si tus abundancias tienen muchos ceros o valores extremos, se puede suavizar:
```
weights.NZ <- sqrt(occ.NZ$abund) / sum(sqrt(occ.NZ$abund))
```
Para mantener comparabilidad, hay que usar el mismo método de ponderación para todas las regiones o especies.
Se puede ajustar la resolución (R) o el suavizado (n) si los nichos quedan muy rugosos o dispersos.

¿Si quieres también te puedo explicar cómo calcular la media ponderada del nicho (por ejemplo, el “centro del nicho” en el PCA ponderado por abundancia)?

Eso sirve para visualizar desplazamientos de nicho en el espacio PCA.

Calcular media ponderada del nicho

Excelente 🌿 — eso es un paso muy interesante: calcular el centro del nicho ponderado (también llamado centroid shift o weighted niche centroid) te permite comparar directamente dónde se ubica el nicho principal de cada grupo (nativo, exótico, etc.) dentro del espacio ambiental del PCA.

Vamos con una explicación breve y luego el código adaptado a tu flujo de trabajo (z.NZ, z.Exotic, PCA-ocu, etc.) 👇

🧭 ¿Qué es el centro del nicho ponderado?

En un PCA ambiental (por ejemplo el PCA-ocu de Broennimann), cada punto representa una celda ambiental del espacio del fondo, y cada especie (o población) tiene una densidad de ocupación sobre ese espacio (z.occ).

El centro del nicho ponderado se calcula como el promedio de las coordenadas de cada celda del PCA, ponderadas por la densidad de ocupación (en este caso, por la abundancia relativa si se usa la versión ponderada):

![[Pasted image 20251030172629.webp]]

donde (w_i) son las probabilidades de ocupación (densidades normalizadas).

🧩 Código adaptado a tu análisis

# -----------------------------------------------------------
# 1️⃣ Función para calcular el centro del nicho ponderado
# -----------------------------------------------------------
niche_centroid <- function(z) {
  # Extrae coordenadas del grid y densidad (z.occ)
  xx <- rep(z$x, each = length(z$y))
  yy <- rep(z$y, times = length(z$x))
  dens <- as.vector(z$z.occ)
  
  # Normalizar (por si no lo está)
  dens <- dens / sum(dens, na.rm = TRUE)
  
  # Calcular centroides ponderados
  cx <- sum(xx * dens, na.rm = TRUE)
  cy <- sum(yy * dens, na.rm = TRUE)
  
  return(c(Cx = cx, Cy = cy))
}

# -----------------------------------------------------------
# 2️⃣ Aplicar a tus nichos ponderados
# -----------------------------------------------------------
cent.NZ     <- niche_centroid(z.NZ)
cent.Exotic <- niche_centroid(z.Exotic)

cent.NZ
cent.Exotic

# -----------------------------------------------------------
# 3️⃣ Calcular distancia entre centros de nicho
# -----------------------------------------------------------
dist.centroids <- sqrt((cent.NZ[1] - cent.Exotic[1])^2 + (cent.NZ[2] - cent.Exotic[2])^2)
dist.centroids

🎨 Visualización en el espacio PCA

Se pueden graficar ambos nichos y sus centroides para ver el desplazamiento:

# Plot de los nichos en el PCA
par(mfrow = c(1,1))
image(z.NZ$x, z.NZ$y, z.NZ$z.occ, col = topo.colors(100), main = "Centros de nicho ponderado", xlab = "PCA 1", ylab = "PCA 2")
contour(z.Exotic$x, z.Exotic$y, z.Exotic$z.occ, add = TRUE, drawlabels = FALSE, col = "red")

# Añadir centroides
points(cent.NZ[1], cent.NZ[2], pch = 19, col = "blue", cex = 1.5)
points(cent.Exotic[1], cent.Exotic[2], pch = 19, col = "red", cex = 1.5)
arrows(cent.NZ[1], cent.NZ[2], cent.Exotic[1], cent.Exotic[2],
       length = 0.1, lwd = 2, col = "darkgray")

legend("topright", legend = c("Nativo", "Exótico"),
       pch = 19, col = c("blue", "red"))

📊 Interpretación

Resultado	Qué indica
`cent.NZ`, `cent.Exotic`	Coordenadas del centro del nicho en el espacio PCA (ponderadas por abundancia).
`dist.centroids`	Distancia euclídea entre los centros → una medida del desplazamiento del nicho.
Gráfico	Visualiza la magnitud y dirección del desplazamiento.

🧠 Opcional: representar el desplazamiento como vector ambiental

Si quieres traducir ese desplazamiento a las variables ambientales originales del PCA, se puede “retroproyectar” los centroides al espacio original:

# Retroproyección al espacio de variables originales
centroids_PCA <- rbind(cent.NZ, cent.Exotic)
colnames(centroids_PCA) <- c("Axis1", "Axis2")

centroids_env <- centroids_PCA %*% t(pca.env.global$c1)  # c1 = loadings de variables en el PCA
centroids_env

Esto te da una idea de qué variables ambientales explican más el cambio en el centro del nicho (por ejemplo, temperatura, precipitación, etc.).

¿Querés que te agregue un bloque final para calcular y graficar el vector de desplazamiento del nicho en el espacio ambiental original (por variable, no solo en el PCA)?
Eso te permitiría ver, por ejemplo, si el desplazamiento del exótico fue hacia temperaturas más altas, menor precipitación, etc.

Graficar vector de desplazamiento de nicho

Perfecto 🌎 vamos a extender el análisis — vas a poder ver hacia dónde se desplazó el nicho exótico en términos de las variables ambientales originales (no solo en ejes de PCA).

Supongamos que ya tienes:

pca.env.global → tu PCA global (el que usaste en el PCA-ocu).
z.NZ, z.Exotic → tus objetos de nicho (densidades) construidos con ecospat.grid.clim.dyn.

🧭 Paso 1. Calcular centroides ponderados en el espacio PCA

Esto es igual que antes:

niche_centroid <- function(z) {
  xx <- rep(z$x, each = length(z$y))
  yy <- rep(z$y, times = length(z$x))
  dens <- as.vector(z$z.occ)
  dens <- dens / sum(dens, na.rm = TRUE)
  cx <- sum(xx * dens, na.rm = TRUE)
  cy <- sum(yy * dens, na.rm = TRUE)
  return(c(Cx = cx, Cy = cy))
}

cent.NZ     <- niche_centroid(z.NZ)
cent.Exotic <- niche_centroid(z.Exotic)

🔁 Paso 2. Retroproyectar los centroides al espacio de variables originales

El objeto pca.env.global (resultado de dudi.pca) tiene las cargas de las variables en pca.env.global$c1, que permiten reconstruir valores aproximados en el espacio original.

# Crear matriz con los centroides en el espacio PCA
centroids_PCA <- rbind(cent.NZ, cent.Exotic)
rownames(centroids_PCA) <- c("Nativo", "Exótico")

# Retroproyectar a las variables originales
# (cada fila: un grupo; cada columna: una variable ambiental)
centroids_env <- centroids_PCA %*% t(pca.env.global$c1)

centroids_env

Esto te da, para cada grupo, el valor medio ponderado estimado de cada variable ambiental (en las mismas unidades estandarizadas del PCA).

📊 Paso 3. Calcular y visualizar el vector de desplazamiento ambiental

# Calcular desplazamiento por variable
shift_env <- centroids_env["Exótico", ] - centroids_env["Nativo", ]
shift_env

# Convertir a data frame para graficar
df_shift <- data.frame(
  Variable = names(shift_env),
  Desplazamiento = as.numeric(shift_env)
)

# Barplot del desplazamiento
library(ggplot2)

ggplot(df_shift, aes(x = reorder(Variable, Desplazamiento), y = Desplazamiento)) +
  geom_bar(stat = "identity", fill = "steelblue") +
  coord_flip() +
  labs(title = "Desplazamiento del centro del nicho (Exótico vs Nativo)",
       x = "Variable ambiental",
       y = "Cambio en el centro (Exótico - Nativo)") +
  theme_minimal()

🧠 Interpretación

Valores positivos → el nicho exótico se desplazó hacia valores mayores de esa variable ambiental.
Valores negativos → el nicho exótico está asociado a valores menores.
La magnitud indica cuánto contribuye cada variable al desplazamiento del nicho.

Por ejemplo:

Si bio1 (temperatura media anual) es positiva → el nicho exótico se movió hacia ambientes más cálidos.
Si bio12 (precipitación anual) es negativa → el nicho exótico ocupa ambientes más secos.

🌈 (Opcional) Mostrar ambos centroides en el espacio de variables

centroids_df <- as.data.frame(centroids_env)
centroids_df$Grupo <- rownames(centroids_df)

# Gráfico tipo radar
library(fmsb)

radar_data <- rbind(
  apply(centroids_df[, -ncol(centroids_df)], 2, max),
  apply(centroids_df[, -ncol(centroids_df)], 2, min),
  centroids_df[centroids_df$Grupo == "Nativo", -ncol(centroids_df)],
  centroids_df[centroids_df$Grupo == "Exótico", -ncol(centroids_df)]
)

rownames(radar_data) <- c("Max", "Min", "Nativo", "Exótico")

radarchart(radar_data, axistype = 1, 
           pcol = c("blue", "red"), plwd = 2, plty = 1,
           title = "Centros de nicho por variable ambiental")
legend("bottomright", legend = c("Nativo", "Exótico"),
       col = c("blue", "red"), lty = 1, lwd = 2, bty = "n")

← Volver al jardín

Primeros pasos

🔹 Contexto: cómo ecospat calcula las densidades

🔹 Qué pasa si se tienen abundancias o probabilidades de presencia

💡 Opción 1 — Usar un kernel density ponderado manualmente

💡 Opción 2 — Expandir los puntos proporcionalmente a su abundancia

💡 Opción 3 — Si tus abundancias son probabilidades (0–1)

🔹 Limitación importante

✅ En resumen

Ejemplo de Kernel Density Ponderado

🧩 Paso 1 — Datos de ejemplo

🧮 Paso 2 — Calcular densidad ponderada con kde2d.weighted

🌍 Paso 3 — Calcular la densidad ponderada

🎨 Paso 4 — Visualizar la densidad

💡 Opcional: comparar con densidad sin ponderar

Integrar mapa ponderado en flujo de ecospat

🧭 Objetivo

🔹 Paso 1 — Crear el nicho estándar en el espacio PCA

🔹 Paso 2 — Calcular tu densidad ponderada (como antes)

🔹 Paso 3 — Reemplazar la densidad en el objeto z.sp

🔹 Paso 4 — Usarlo en los tests de ecospat

⚠️ Consideraciones importantes

✅ En resumen

Código completo para nicho ponderado

🧩 Código completo para nicho ponderado por abundancia

📘 Notas importantes

Código adaptado a mi caso

🧩 Versión adaptada a tu análisis: PCA-Ocu (Broennimann et al.)

💬 Cómo interpretarlo

⚠️ Tips finales

Calcular media ponderada del nicho

🧭 ¿Qué es el centro del nicho ponderado?

🧩 Código adaptado a tu análisis

🎨 Visualización en el espacio PCA

📊 Interpretación

🧠 Opcional: representar el desplazamiento como vector ambiental

Graficar vector de desplazamiento de nicho

🧭 Paso 1. Calcular centroides ponderados en el espacio PCA

🔁 Paso 2. Retroproyectar los centroides al espacio de variables originales

📊 Paso 3. Calcular y visualizar el vector de desplazamiento ambiental

🧠 Interpretación

🌈 (Opcional) Mostrar ambos centroides en el espacio de variables

🔹 Contexto: cómo `ecospat` calcula las densidades

🧮 Paso 2 — Calcular densidad ponderada con `kde2d.weighted`

🔹 Paso 3 — Reemplazar la densidad en el objeto `z.sp`

🔹 Paso 4 — Usarlo en los tests de `ecospat`