Ajuste inicial de los datos

Los modelos de distribución de especies son por tanto representaciones cartográficas de la idoneidad de un espacio para la presencia de una especie en función de las variables empleadas para generar dicha representación (Mateo, Felicísimo & Jesús Muñoz, 2011). Un aspecto importante a tener en cuenta para aplicar este modelo correctamente es tener un directorio de trabajo ordenado, donde se almacenarán los datos de ocurrencia, los predictores ambientales y los modelos candidatos. Cada uno de estos elementos se irán explicando a medida que se avance en la rutina.
Los datos de ocurrencia de especies pueden adquirirse de dos maneras: descargándolos de la página https://www.gbif.org/es/occurrence/search o mediante el paquete rgbif. Con ayuda de esta última fuente es más fácil acceder a una mayor cantidad de datos, los cuales ya tienen el formato para ser trabajados en R y pueden filtrarse con mayor facilidad las variables de interés. 
library(rgbif)
## Warning: package 'rgbif' was built under R version 3.6.3
#library(tidyverse)

#Buscamos los registros para la especie Panthera onca y filtramos para trabajar solo con aquellos registros con coordenadas geográficas y sin errores Geoespaciales.
ponca <- occ_data(scientificName = 'Panthera onca', 
                limit = 200000, 
                hasCoordinate = TRUE, 
                hasGeospatialIssue = FALSE)

# Los datos descargados por este paquete tienen 168 variables. Para correr este modelo solamente necesitamos el nombre científico, la longitud y la latitud, en ese orden.
ponca_occ <- as.data.frame(na.omit(ponca$data[ , c("scientificName", "decimalLongitude", "decimalLatitude")]))
View(ponca_occ)
ponca_occ$scientificName <- "Panthera onca"
names(ponca_occ) <- c("Species", "Longitude", "Latitude")

write.csv(ponca_occ, 'occ_Panthera_onca.csv', row.names = F)
Ahora que se tienen los datos de la especie, es necesario dividirlos en datos de entrenamiento, los cuales serán empleados para el desarrollo de cada modelo, y en datos de prueba, con los que se podrá determinar la exactitud de cada uno de estos. El primer paso para este ajuste de los datos es la instalación y llamada del paquete kuenm. 
if(!require(devtools)){
    install.packages("devtools")
}
## Loading required package: devtools
## Loading required package: usethis
if(!require(kuenm)){
    devtools::install_github("marlonecobos/kuenm")
}
## Loading required package: kuenm
library(kuenm)
Ahora se realiza la división de los datos en los grupos de entrenamiento y prueba. Para esto se utiliza la función kuenm_occsplit. 
set.seed(2020)
split <- kuenm_occsplit(occ = ponca_occ, train.proportion = 0.7, method = "random", save = TRUE, "ponca")
Ahora seleccionamos las variables ambientales que queremos en el modelo. Para este caso se tomaron las variables bioclimáticas bio1 y bio2 que fueron descargadas y ubicadas en el directorio de trabajo. 
library(dismo)
## Loading required package: raster
## Loading required package: sp
bioclim.data <- getData(name = "worldclim",
                        var= "bio",
                        res= 2.5)
Con las líneas de código anterior se descargaron las 19 variables climaticas Bioclim. En la carpeta del directorio de trabajo se eliminaron todas las variables exceptuando las bio1 y bio2, y con ayuda de ArcGIS se cambió el formato de .bil (Raster) a .asc (Ascii). 
vs <- kuenm_varcomb(var.dir = "wc2-5", out.dir = "M_variables", min.number = 2, in.format = "ascii", out.format = "ascii")
## Warning in dir.create(out.dir): 'M_variables' already exists
## 
## A total of 1 sets of variables resulted from combinations of 2 variables will be written.
## Warning in dir.create(sub_paths[k]): 'M_variables\Set_1' already exists

Proceso de calibración de los modelos

El proceso de calibración es un proceso en el cual el objetivo es determinar cuál es la combinación de parámetros que mejor representa el fenómeno de interés al encontrar el mejor ajuste de los datos. El número de modelos candidatos dependerá de tres elementos: 1) los diferentes sets de variables, que en este caso solo es 1 ya que solo se trabajan con dos variables, 2) los diferentes parámetros multiplicadores de regularización de MaxEnt y 3) combinaciones de clases de entidades en el que encontramos cinco diferentes clases (lineal, cuadrático, producto, umbral y bisagra). Es necesario recalcar que para cada combinación se construyen dos modelos, uno con el set total de datos y otro solo con los datos de entrenamidento. 
occ_joint <- "ponca_joint.csv"
occ_tra <- "ponca_train.csv"
M_var_dir <- "M_variables"
batch_cal <- "Candidate_models"
out_dir <- "Candidate_Models_P"
reg_mult <- c(0.1, 0.25, 0.5, 1, 2, 4)
f_clas <- c("lq", "lqp", "q", "h", "t", "lpth")
args <- NULL
maxent_path <- "D:/ENM_HUMBOLDT/maxent"
wait <- FALSE
run <- TRUE
Ya con los parámetros definidos el paso siguiente es aplicarlos para construir los diferentes modelos candidatos. Para esto se utilizará la función kuenm_cal. 
kuenm_cal(occ.joint = occ_joint, occ.tra = occ_tra, M.var.dir = M_var_dir, batch = batch_cal,
          out.dir = out_dir, reg.mult = reg_mult, f.clas = f_clas, args = args,
          maxent.path = maxent_path, wait = wait, run = run)
Existe un total de 36 modelos de distribución que combinan todos los parámetros definidos para Panthera onca, pero como cada modelo se realiza primero con el set completo de datos y luego con los datos de entrenamiento, hay un total de 72 modelos. Ahora se deben elegir los modelos más adecuados y que nos brinden información de calidad sobre el nicho ecológico de la especie.

Evaluación y selección de los mejores modelos

Para medir el rendimiento de los modelos se utilizarán tasas de omisión, que indican qué tan bien creados son los modelos al predecir con los datos de entrenamiento, los datos de prueba. El umbral que se utilizará para esta evaluación será del 5%. También se evaluará la complejidad del modelo, teniendo en cuenta que un modelo es menos favorable entre mayor sea su complejidad. Esta característica se medirá con el Criterio de Información de Akaike corregido (AICc), el delta AICc y los pesos AICc. 
occ_test <- "ponca_test.csv"
out_eval <- "Calibration_results"
threshold <- 5
rand_percent <- 50
iterations <- 100
kept <- TRUE
selection <- "OR_AICc"
paral_proc <- FALSE
Al tener los argumentos de la función definidos, se procede a evaluar la totalidad de modelos candidatos y luego seleccionar los mejores modelos de acuerdo con los tres criterios mencionados aneriormente. Para esto se utilizará la función kuenm_ceval. Desafortunadamente, dado que la base de datos contiene un elevado número de registros no me fue posible construir los modelos. Aplicando las líneas de código siguientes, que aparecen comentadas, generaría los modelos finales y cada uno de ellos contaría con los parámetros de selección del mejor modelo. 
#cal_eval <- kuenm_ceval(path = out_dir, occ.joint = occ_joint, occ.tra = occ_tra, occ.test = occ_test, batch = batch_cal, out.eval = out_eval, threshold = threshold, rand.percent = rand_percent, iterations = iterations, kept = kept, selection = selection, parallel.proc = paral_proc)
Para completar la actividad, se hará el mismo ejercicio que se realizó hasta este punto, pero con datos de las especies de Amblyomma americanum y Peltophryne empusa disponibles en la dirección de gitHub https://github.com/marlonecobos/kuenm#installing-the-package. Los resultados se presentan en nuevo Rmarkdown.