Distribuciones geográficas

Juliana Herrera-Pérez & Fabricio Villalobos

Área de distribución geográfica

“…el espacio donde las condiciones ecológicas favorecen, real o potencialmente y en varios niveles, las interacciones no efímeras de los individuos de una especie” Mota-Vargas & Rojas-Soto 2012

Puntos, polígonos y mapas de distribución

Los registros (colectas georeferenciadas) de las especies son los datos primarios de biodiversidad, a partir de los cuáles podemos estimar las áreas de distribución de estas y describir/evaluar los patrones de diversidad que emergen de su agregación (traslape; e.g., gradiente geográfico de riqueza)

  • En este ejemplo/ejercicio veremos cómo obtener dichos registros directamente desde R
  • También, veremos cómo generar áreas de distribución (extenciones de presencia) a partir de estos registros, creando polígonos de diferentes tipos (mínimo, alpha y alpha dinámico)

Paquetes necesarios:

library(rgbif)
library(dplyr)
library(ggplot2)
library(sf)
library(rnaturalearth)
library(alphahull)
library(rangeBuilder)
library(janitor)

Obteniendo registros de presencia de GBIF

Para comenzar podemos escoger una especie y aplicar la función occ_data

El objeto sp_data es una lista con datos sobre los resultados obtenidos en GBIF (incluyendo algunos metadatos)

Para trabajar únicamente con la tabla de registros hay que seleccionar el objeto data dentro del mismo

sp_data <- occ_data(scientificName = "Musonycteris harrisoni", limit = 500)[[2]]

Checar el nombre de las columnas (para después buscar únicamente las de posición geográfica: lat/long)

names(sp_data)[1:30]
 [1] "key"                    "scientificName"         "decimalLatitude"       
 [4] "decimalLongitude"       "issues"                 "datasetKey"            
 [7] "publishingOrgKey"       "installationKey"        "hostingOrganizationKey"
[10] "publishingCountry"      "protocol"               "lastCrawled"           
[13] "lastParsed"             "crawlId"                "projectId"             
[16] "basisOfRecord"          "individualCount"        "occurrenceStatus"      
[19] "sex"                    "taxonKey"               "kingdomKey"            
[22] "phylumKey"              "classKey"               "orderKey"              
[25] "familyKey"              "genusKey"               "speciesKey"            
[28] "acceptedTaxonKey"       "acceptedScientificName" "kingdom"

1.Crear otro objeto a partir sp_data únicamente con long/lat

2.Quedarse únicamente con los puntos/registros individuales (i.e., excluir duplicados)

3.Transformarlo en un objeto espacial

sp_p1<-sp_data%>%
  select(decimalLongitude,decimalLatitude,species)%>%
  mutate(lat=decimalLatitude,lon=decimalLongitude)%>%
  distinct() %>%
  na.omit() %>% 
  sf::st_as_sf(coords = c('decimalLongitude','decimalLatitude'),crs="EPSG: 4326")

Graficar (poner en un mapa) los puntos de presencia de nuestra especie

ggplot()+ geom_sf(data=sp_p1, col="blue",pch=19)

Agregar el mapa del mundo para saber qué onda!

wrld <- ne_countries(scale = "small",returnclass = "sf")
ggplot()+
  geom_sf(data=wrld)+geom_sf(data=sp_p1,col="blue",pch=19,size=1)+coord_sf(expand = F) +
  labs(x="Longitud decimal ",
       y="Latitud decimal",
       title=expression(paste("Puntos reportados ", italic("Musonycteris harrisoni"))))+
  theme(plot.title = element_text(hjust = 0.5))

Hay algo claramente equivocado, ¿cierto? Los puntos/registros necesitan ser “curados” (limpiados)

Eliminar los puntos con mala georeferencia (en este caso, puntos obvios en el “viejo mundo”)

sp_p1<-sp_data%>%
  select(decimalLongitude,decimalLatitude,species)%>%
  mutate(lat=decimalLatitude,lon=decimalLongitude)%>%
  distinct() %>% na.omit() %>% 
  sf::st_as_sf(coords = c('decimalLongitude','decimalLatitude'),crs="EPSG: 4326")%>%
  filter(lat> 0.5) %>% filter(lat< 22)

Ahora sí, mapeamos de nuevo pero sólamente en la región de interés (México)

mex_map <- filter(wrld,name=="Mexico")
ggplot()+geom_sf(data=mex_map)+
  geom_sf(data=sp_p1,col="blue",pch=19,size=1)+coord_sf(expand = F)

Y ¿Cómo eliminamos los registros que están en el mar?

sp_p1<-sp_data%>%
  select(decimalLongitude,decimalLatitude,species)%>%
  mutate(lat=decimalLatitude,lon=decimalLongitude)%>%
  distinct() %>%
  na.omit() %>% 
  sf::st_as_sf(coords = c('decimalLongitude','decimalLatitude'),crs="EPSG: 4326")%>%
  filter(lat> 0.5) %>%
  filter(lat< 22)%>%
  filter(lon> -105.56611)

Polígono convexo mínimo

Una vez tenemos los datos curados, podemos crear nuestro mcp

sp1_mcp <- st_convex_hull(st_union(sp_p1) )

¿Cómo se ve?

sp1_mcp2 <- st_as_sf(sp1_mcp)

ggplot()+
  geom_sf(data=mex_map)+
  geom_sf(data=sp1_mcp,
          fill="blue")

Polígono alfa (alpha hull)

Usamos el paquete alphahull

NOTA: Esta función solo acepta tablas como entrada

sp_p2<-as.data.frame(st_coordinates(sp_p1))
sp1_alphahull <- ahull(sp_p2, alpha = 6)
Error in eval(expr, envir, enclos): shull: duplicate points found

Error: shull: duplicate points found

Falla porque encuentra puntos duplicados o, como en este caso, puntos en una línea recta (i.e, mismo X y/o mismo Y).

¿Cómo podemos identificar y solucionar este error?

sp_p3<-sp_p2 %>% select(X, Y)%>% 
mutate(X = ifelse(duplicated(sp_p2$X), X + rnorm(1, mean = 0, sd = 0.0001), X))%>% 
mutate(Y = ifelse(duplicated(sp_p2$Y), Y + rnorm(1, mean = 0, sd = 0.0001), Y))

Ahora si, podemos crear el Alpha Hull con un valor de alpha escogido (por la razón que crean relevante)

sp1_alphahull<- ahull(sp_p3, alpha = 1)

Para observar el alpha hull, necesitamos que el objeto sea de tipo espacial del paquete sf. Para eso usaremos una función independiente, disponible en su carpeta de trabajo

source(file = here::here("data","ah2sf.R"))
sp1_alphahull.poly <- ah2sf(sp1_alphahull)

¿Cómo se ve?

ggplot()+
  geom_sf(data=mex_map)+geom_sf(data=sp1_alphahull.poly,fill="blue")

Polígono alfa dinámico

Usamos el paquete rangeBuilder, el cual crea un polígono alpha hull con un valor de alpha “óptimo” basado en la distribución espacial de los puntos

sp_range <- getDynamicAlphaHull(
  sp_p3, #Tabla de puntos/registros de la especie
  coordHeaders = c("decimalLongitude", "decimalLatitude"),# x y y
  fraction = 0.95,   # la fracción mínima de registros que debe incluir el polígono
  partCount = 2,  # el máximo de polígonos disyuntos permitidos
  initialAlpha = 1, # Alpha inicial
  alphaIncrement = 0.5,
  alphaCap = 1000,
  clipToCoast = "terrestrial"  # solo la parte terrestre del polígono se mantendrá (se cortan las partes no-terrestres/acuáticas con base en un mapa descargado de naturalearth).
)
alpha <- sp_range[[2]]
alpha
[1] "alpha1.5"

Convertir el polígono alpha a un objeto sf

sp1_dynalpha <- st_make_valid(st_as_sf(sp_range[[1]]))

¿Cómo podemos visualizarlo?

ggplot()+ geom_sf(data=mex_map)+ 
  geom_sf(data=sp1_dynalpha, fill="blue")

¿Y ….Cómo se ven todos los polígonos?

ggplot()+
  geom_sf(data=mex_map)+ geom_sf(data=sp1_mcp,fill="red",alpha=0.1) +
  geom_sf(data=sp_range[[1]], fill="blue",alpha=0.5)+ 
  geom_sf(data=sp1_alphahull.poly,fill="yellow",alpha=0.5)

Finalmente, podemos salvar esos polígonos como shapefiles, para usarlos en otros software (e.g. ArcGIS) y eventualmente juntar los de varias especies para otros análisis (ejercicio siguiente)

st_write(sp1_mcp2, "sp1_min_convex.shp")
st_write(sp1_alphahull.poly, "sp1_alphahull.shp")
st_write(sp1_dynalpha, "sp1_dyn_alphahull.shp")