areas_vida_1.qmd

---
format:
  html:
    toc: true
    toc-title: "Contenidos"
    toc-expand: 1
    code-overflow: wrap
editor_options:
  markdown:
    wrap: sentence
    canonical: true
theme: cosmo
lang: es
---
```{r}
#| include: false
source("_common.R")

data <- read_csv2("docs/datos/zebras.csv")
zebras <- data[, c("event-id", "individual-local-identifier",
                "location-long", "location-lat", "study-local-timestamp")]
dos_zebras <- zebras %>%
              rename(id = "event-id", identifier = "individual-local-identifier",
              long = "location-long",
              lat = "location-lat",
              timestamp = "study-local-timestamp") |>
          filter(identifier == "Z3864" | identifier == "Z6405") |>
          dplyr::select(identifier, long, lat, timestamp) |> na.omit()
zebras.proj <- SpatialPointsDataFrame(coords = as.data.frame(cbind(dos_zebras$long, dos_zebras$lat)),
                data = dos_zebras,
                proj4string = CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))
zebras.transf <- spTransform(zebras.proj,
                  CRS("+proj=lcc +lat_1=20 +lat_2=-23 +lat_0=0 +lon_0=25 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
```

# Áreas de vida {#sec-areas-vida}

La estimación de áreas de vida es una parte fundamental, pero no la única, en el estudio de la ecología y etología de los animales, y los artículos científicos de @kays_terrestrial_2015, @abrahms_emerging_2021, y @silva_autocorrelation-informed_2022 son una excelente fuente de información si estás involucrándote en el mundo de la ecología de movimiento. Si bien existen diversos métodos por los cuales se pueden estimar dichas áreas, en este módulo aprenderás sobre métodos tradicionales ampliamente utilizados y que serán comparados más adelante con nuevos métodos en el paquete _ctmm_.

## _adehabitatHR_

El paquete [_adeHabitatHR_](https://cran.r-project.org/web/packages/adehabitatHR/vignettes/adehabitatHR.pdf) es uno de cuatro paquetes que hace unos años formaban parte del paquete ya deprecado _adehabitat_, y se especializa en el análisis de áreas de vida. Este paquete, sin embargo, se enfrenta a un riesgo a futuro ya que otros paquetes utilizados dentro de la disciplina del análisis espacial, dependían en gran medida de los paquetes _raster_, _rgdal_, _rgeos_, y _sp_, que fueron reemplazados en el 2023 por _sf_ y _terra_. _sp_ es el único que se adaptó a esta transición y adoptó funciones del paquete _sf_ Sin embargo, _adehabitatHR_ no se adaptó a esta transformación y depende de la instalación de estos paquetes que ya no reciben mantenimiento.

A pesar de esto, aprender sobre _adeHabitatHR_ es aún necesario ya que no existe, a la fecha, un paquete que incorpore todas las utilidades que éste lo hace.

Ahora estimarás el área de vida de las zebras que fueron extraídas en el objeto `zebras.transf` en la @sec-transf-proy.

```{r}
#| code-fold: true
#| code-summary: "Obten el objeto `zebras.transf`"

zebras.transf <- spTransform(zebras.proj,
                  CRS("+proj=lcc +lat_1=20 +lat_2=-23 +lat_0=0 +lon_0=25 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
```

## Minimum Convex Polygon (MCP)

### Zebras en África
Este método es de los más sencillo para dibujar los límites de distribución de un animal. Aunque su uso original Carl O. Mohr (1947) fue destinado a la identificación de animales recapturados en una malla de captura.

El uso del MCP es limitado ya que este describe el alcance de la distribución de ubicaciones de un individuo, mas no la verdadera área de vida del animal bajo estudio. Para una discusión más detallada sobre esto puedes leer @signer_fresh_2021, y un debate profundo sobre rangos y ocurrencia de distribución en @fleming_estimating_2016.

A pesar de lo anteriormente mencionado, obtener el área del polígono y el poder ubicarlo sobre un mapa, nos permite observar a _grosso modo_ el espacio y hábitat que ocupa cierto ejemplar.

```{r}
library(adehabitatHR)
library(sf)
zebras.mcp <- mcp(zebras.transf[, "identifier"], percent= 100,
              unin="m", unout= "ha")
```

:::::: {.callout-note}
## Ejercicio
Ahora intenta con distintos porcentajes y observa las diferencias en tamaños de área.
:::



Si bien puedes utilizar funciones base de R para graficar estos polígonos, aprovecha la interacción de estos paquetes con _sf_ y _ggplot2_ para obtener un gráfico más estético.

```{r}
library(ggplot2)

ggplot() +
      geom_sf(data = st_as_sf(zebras.mcp), aes(fill = id),
              alpha = 0.2) +
      geom_sf(data = st_as_sf(zebras.transf), aes(colour = identifier)) +
      coord_sf() +
      theme_bw()
```

También puedes exportar estos polígonos con extensión `kml` para que puedas importarlos en Google Earth o cualquier software de tu preferencia.
```{r eval = FALSE}
st_write(st_as_sf(zebras.mcp), "zebras.kml", delete_layer = TRUE)
```

### Animales hipotéticos

Como caso de estudio, e intentando que mejores en el flujo de trabajo de estos análisis, creé un set de datos de tres mamíferos hipotéticos en la zona Norte de Loreto - Perú, y lo puedes [descargar aquí](https://mega.nz/folder/VwVQEbJK#lxYXI88_90s6hUihhwNXgg).

```{r include = FALSE}
mamiferos <- read.csv("docs/datos/mamiferos.csv")
```

```{r eval = FALSE}
mamiferos <- read.csv("mamiferos.csv")
```

Como siempre, haz una exploración rápida de este objeto y grafica los puntos satelitales rápidamente.

```{r}
plot(mamiferos[, c("longitude", "latitude")], pch=20)
```

:::::: {.callout-note}
## Ejercicio
Realiza una limpieza de este set de datos, elimina outliers y datos faltantes, y crea un objeto que contenga únicamente a "juancho". Finalmente, verifica que el proceso haya sido exitoso.
:::

```{r}
# Verificación visual
library(plotly)

mamiferos_outliers <- mamiferos %>%
  ggplot(aes(x= longitude, y=latitude, colour = name))+
  geom_point(alpha = 0.5)

ggplotly(mamiferos_outliers)
```

Con los gráficos interactivos de _plotly_ puedes usar el cursor para observar la posición exacta de los _outliers_. Si bien esto es de utilidad para un set de datos pequeño como este, puedes automatizar tu código para deshacerte de todos aquellos puntos menores a la longitud -76.

```{r}
#Limpiar el data frame
juancho <- mamiferos[-which(mamiferos$longitude < -76),] |>
                na.omit() |>
                filter(name == "juancho")
```

Ahora que has filtrado el set de datos, ya puedes realizar la proyeccion, transformación y cálculo del área de vida y núcleo de juancho.

```{r}
juancho.proj <- SpatialPointsDataFrame(coords = as.data.frame(cbind
                        (juancho$longitude, juancho$latitude)),
                        data = juancho, proj4string = CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))
```

A diferencia de las zebras, sabemos que nuestro lugar de muestreo es pequeño, y en una región tropical, por ende necesitamos un método de proyección adecuado. Entre tantos métodos existentes, ¿cuál elegirías?
En estos casos, lo aconsejable es utilizar [Universal Transverse Mercator](https://gisgeography.com/utm-universal-transverse-mercator-projection/) (UTM). No obstante, para hacer uso de este método de proyección, aún nos falta conocer un dato importante, la zona y el hemisferio.

Podemos copiar y pegar cualquiera de los puntos de nuestro set de datos en [esta página](http://rcn.montana.edu/resources/Converter.aspx) para poder conocer la información faltante.

```{r}
juancho.transf <- spTransform(juancho.proj,
                           CRS("+proj=utm +south +zone=18 +ellps=WGS84"))
```

Si deseas conocer las coordenadas UTM de todos los puntos, puedes crear un nuevo objeto y ejecutar las siguientes líneas.
```{r}
juancho.utm <- coordinates(juancho.transf)
head(juancho.utm)
```

:::::: {.callout-note}
## Ejercicio
Estima el área de vida (90%) y núcleo (50%) de juancho mediante el método MCP y realiza un gráfico mostrando ambos polígonos y puntos satelitales.
:::


```{r}
juancho.mcp.90 <- mcp(juancho.transf[, "name"], percent = 90,
                          unin = "m", unout = "ha")
juancho.mcp.50 <- mcp(juancho.transf[, "name"], percent = 50,
                      unin = "m", unout = "ha")
```

Tal como lo hiciste para el set de datos de zebras, puedes transformar los objetos _SpatialPolygonsDataFrame_ a un objeto _sf_ para graficarlo fácilmente.
```{r}
ggplot() +
  geom_sf(data = st_as_sf(juancho.mcp.90),
          color = "NA",
          fill = "black",
          alpha = 0.2) +
  geom_sf(data = st_as_sf(juancho.mcp.50),
          color = "NA",
          fill = "purple",
          alpha = 0.5) +
  geom_sf(data = st_as_sf(juancho.transf), alpha = 0.5) +
  coord_sf() +
  labs(title = "Áreas de vida y núcleo de juancho estimadas mediante MCP") +
  theme_bw()
```

## Kernel Density Estimate (KDE)

En KDE, una colina tri-dimensional o kernel es formada junto a cada punto, la forma y altura de esta colina depende del ancho de banda (_bandwidth_).
En este curso calcularás el KDE mediante el “_fixed kernel_” o href, “_least squares cross-validation_” o lcsv, e intentarás ajustar el bandwidth manualmente.
La teoría sobre estos métodos, sus pros y contras pueden ser encontrados en estos links: [1], [2].

### _Fixed kernel_ (href)

```{r}
juancho.khref <- kernelUD(juancho.transf[, "name"], h = "href")

# Observa este raster
image(juancho.khref)
```

Si quieres exportarlo para su uso en otros softwares como ArcGis puedes seguir estos pasos.

```{r warnings = FALSE, output = FALSE}
juanchopix <- estUDm2spixdf(juancho.khref)
st_write(st_as_sf(juanchopix), "juancho.shp", delete_layer = TRUE)
# Mira otros formatos disponibles con st_drivers()
```

También puedes plotearlo en R.
```{r eval = FALSE}
# Tambien lo puedes graficar en R
plot(juanchopix)
plot(juancho.transf, add = T, cex = 0.1)
```

Utilizando estos objetos, ya puedes estimar el área de vida de juancho mediante el método KDE.

```{r}
# Estimo el area usando el objeto juancho.khref segun una secuencia de porcentajes
kernel.area(juancho.khref, percent = seq (20, 95, 5),
            unin="m", unout="ha")
```

Ahora puedes utilizar la función `getverticeshr` para obtener los polígonos que te permitirán graficar el área de vida y núcleo de juancho.

```{r}
# Polígonos para el plot final
juancho.KDE90 <- getverticeshr(juancho.khref, percent = 90)
juancho.KDE50 <- getverticeshr(juancho.khref, percent = 50)

# Usa ggplot para graficar

ggplot() +
  geom_sf(data = st_as_sf(juancho.KDE90),
          color = "NA",
          fill = "black",
          alpha = 0.2) +
  geom_sf(data = st_as_sf(juancho.KDE50),
          color = "NA",
          fill = "purple",
          alpha = 0.5) +
  geom_sf(data = st_as_sf(juancho.transf), alpha = 0.5) +
  coord_sf() +
  labs(title = "Áreas de vida y núcleo de juancho estimada mediante KDE") +
  theme_bw()
```

:::::: {.callout-note}
## Ejericicio
La función `getvolumeUD` también te permite calcular el área de vida, pero ¿cuál crees que sea la diferencia entre ambos métodos?
:::


## Movebank

Con el advenimiento de nuevas tecnologías para el rastreo de animales, algunos problemas se hicieron evidentes debido al cumulo de información generada por los nuevos dispositivos. Esto conllevó al desarrollo del sistema [_Movebank_](https://www.movebank.org/cms/movebank-main), para poder mejorar el manejo del gran volumen y diversidad de información (@kays_movebank_2022).

[Movebank](https://www.movebank.org/cms/movebank-main) es una plataforma web que permite almacenar datos de telemetría (GPS, Argos, VHF, etc)  de manera privada y/o pública. A la fecha, cuentan con un total de 6 billones de locaciones generadas por más de 8 mil estudios que compreden a aproximadamente 1500 taxones. Además, almacena una gran diversidad de información asociada a cada proyecto (persona de contacto, grants, cita, licencia, etc), animal (Nombre, especie, peso, edad, etc) y a los sensores del collar (temperatura externa, altitud, frecuencia cardíaca, DOP, etc.

Estas variables han sido estandarizadas mediante lenguaje persistente, mediante una lista de vocabulario extensa que puede ser leído por máquinas y que está disponible en el Natural Environment Research Council's Vocabulary Server (mira los detalles [aquí](https://www.movebank.org/cms/movebank-content/movebank-attribute-dictionary)).

El trabajo que realizarás a continuación utilizará una porción de los datos recolectados por @castellanos_pilot_2022 de un zorro andino en Ecuador que puedes obtenerla [aquí](https://mega.nz/folder/F01TibSK#DHHrerhI5tiCFSbYISLLjA). Sin embargo, si quieres trabajar con las más de 6 mil posiciones satelitales, puedes acceder al [repositorio de Movebank](doi:10.5441/001/1.1q2c075p), o ingresar a [www.movebank.org](www.movebank.org) y buscar el proyecto "_Home range and movement patterns of the Andean Fox in Cotopaxi National Park Ecuador_".

## Autocorrelated Kernel Density Estimate (AKDE)

Lo descrito en esta sección es una parte muy pequeña de todo lo que puedes lograr con este paquete. Puedes encontrar mucho más en la [página oficial de ctmm](https://ctmm-initiative.github.io/ctmm/index.html). Además, el [grupo de google](https://groups.google.com/g/ctmm-user) de ctmm es muy activo y encontrarás respuestas a posibles problemas que puedas tener con este paquete. La teoría sobre cada función ha sido publicada en varias revistas, y estas funciones se expandirán aún más en próximos años.

### Error de posición

El primer paso en este flujo de trabajo es estimar el error de posición GPS. La explicación matemática de este parámetro se detalla en @fleming_comprehensive_2020, acompañada por varias medidas de error de una gran variedad de collares. En caso de que no poseas información de calibración del collar o si no fuiste capaz de obtenerla, esta [viñeta de R](https://ctmm-initiative.github.io/ctmm/articles/error.html) te explica un poco más sobre este proceso.

Aquuí usarás el error de calibración oportunista que se obtuvo tras la muerte del zorro andino Mashca.

```{r include = FALSE}
calibracion <- read.csv("docs/datos/mashca_calibration.csv")
calibracion <- calibracion[, c(3:5,8,9,15, 18,19)] |> na.omit()
```

```{r eval = FALSE}
calibracion <- read.csv("mashca_calibration.csv")
```

::::::: {.callout-note}
## Ejercicio
Remueve los NA's que se encuentran en las columnas de longitud y latitud del objeto `calibracion`.
:::

Ahora plotea este objeto utilizando cualquiera de las columnas de posición.
```{r}
plot(calibracion$utm.easting, calibracion$utm.northing)
```

Convierte este objeto a una nueva clase de objeto con la función `as.telemetry`. Asegúrate de saber cual es la zona horaria de tu lugar de estudio y del [código EPSG](https://spatialreference.org/ref/epsg/32717/), o puedes utilizar los argumentos `projstring` que utilizamos en @sec-transf-proy para la transformaciónn a UTM.

```{r}
library(ctmm)
calib.ctmm <- as.telemetry(calibracion, timezone = "America/Bogota", projection =
                                CRS("+init=epsg:32717"))
summary(calib.ctmm)
```

Ahora plotea este objeto, y compáralo con los gráficos anteriores. ¿Qué diferencias notas?
```{r}
plot(calib.ctmm, error = 2)
```

Finalmente, debemos aplicar el error de locación al objeto `calib.ctmm`.

```{r}
calib.model <- uere.fit(calib.ctmm)
summary(calib.model)
```

De esta manera, se puede estimar un error de posición de GPS de aproximadamente 36 metros. En este tutorial no removerás outliers, pero siempre es necesario aplicarlo a tu trabajo y encontrarás una explicación muy sencilla y directa sobre como hacerlo en la [página oficial de ctmm](https://ctmm-initiative.github.io/ctmm/articles/error.html#outlier-detection).

Para mostrar los beneficios del uso del método AKDE, he separado el set de datos de mashca aproximadamente en dos mitades.

```{r include = FALSE}
mashca_first <- read.csv("docs/datos/mashca_first_half.csv")
mashca_second <- read.csv("docs/datos/mashca_second_half.csv")
```

```{r eval = FALSE}
mashca_first <- read.csv("mashca_first_half.csv")
mashca_second <- read.csv("mashca_second_half.csv")
```

Como es habitual, limpia estos objetos e importalos a _ctmm_ tal como lo hiciste para la información de calibración.

```{r}
mashca_first <- mashca_first[, c(3:5,8,9,15,18,19, 22)] %>%
    na.omit() %>%
  filter(gps.dop < 3, gps.fix.type.raw == "val. GPS-3D")
mashca_second <- mashca_second[, c(3:5,8,9,15,18,19, 22)] %>%
  na.omit() %>%
  filter(gps.dop < 3, gps.fix.type.raw == "val. GPS-3D")

# Importa a ctmm
mashca_first_ctmm <- as.telemetry(mashca_first, timezone = "America/Bogota",
                            projection = CRS("+init=epsg:32717"))

mashca_second_ctmm <- as.telemetry(mashca_second, timezone = "America/Bogota",
                                  projection = CRS("+init=epsg:32717"))
```
Puedes utilizar _ctmm_ con listas, lo que mejorara tu flujo de trabajo cuando analices varios individuos.

```{r}
# Une ambos objetos en una lista
mashca_merged <- list(mashca_first_ctmm, mashca_second_ctmm)
names(mashca_merged) <- c("first","second")

# Aplica el modelo de calibracion a la lista
uere(mashca_merged) <- calib.model

# Observa los datos
plot(mashca_merged[[2]], error=2)
```

### Variogramas

El método implementado en _ctmm_ requiere que los animales presenten un comportamiento de rango-residencia indicado por la función de semivarianza (SVF) de un proceso de movimiento estocástico, el cual puede ser visualizado mediante un variograma. Los detalles sobre este método se encuentran en @fleming_fine-scale_2014.

Debido a que los variogramas incorporan el efecto de tasa de muestreo para evaluar distintos comportamientos, es recomendable evaluar el intervalo en el cual los puntos de locación fueron tomados.
```{r}
# Obtener grafico de intervalos
dt.plot(mashca_merged[[2]])
ablines = abline(h = c(1, 2, 4) %#% "hours", col = "red")
# Aplica este intervalo en un nuevo objeto
dt = c(c(2, 4) %#% "hours")

# Evalua el variograma
vg <- c()
for (i in 1:length(mashca_merged)){
vg[[i]] <- variogram(mashca_merged[[i]], dt)
}

# Puedes mirar este variograma como un gráfico estático o con zoom()
plot(vg[[2]], fraction = 0.9, level=c(0.5,0.95))
```

### _Model fitting_

```{r include = FALSE}
load("docs/datos/mashca_fits.RData")
```
Tras incorporar el error de posición a los datos y seleccionar un variograma, debes realizar una aproximación empírica a lo que el modelo "pareceria ser". Esto puede lograrse con la función `guess`.

```{r eval = FALSE}
# guess
guess <- c()
for (i in 1:length(mashca_merged)){
guess[[i]] <- ctmm.guess(mashca_merged[[i]],
                   variogram = vg[[i]],
                   interactive = FALSE)
guess[[i]]$error <- TRUE
}
```
También puedes hacerlo manualmente si es que tienes el suficiente conocimiento sobre los parámetros que estás eligiendo.
```{r eval = FALSE}
ctmm.guess(mashca_merged[[1]],
           variogram = vg[[1]])
```

```{r eval = FALSE}
# Ajusta un modelo
fits_first <- ctmm.select(mashca_merged[[1]], CTMM = guess[[1]], trace = 3,
                    cores = 0)
fits_second <- ctmm.select(mashca_merged[[2]], CTMM = guess[[2]], trace = 3,
                    cores = 0)
```

Puedes inspeccinar ambos objetos con `summary` y verás algo como esto:

```{r echo = FALSE}
# FITS
summary(fits_second)
```

En breve, ambos objetos muestran que no hay suficiente señal en la información para poder estimar velocidades no lineales ya que los grados de libertad de la velocidad son iguales a 0 y el modelo elegido fue OU anisotropic (lee @fleming_fine-scale_2014). Sin embargo, hay suficiente información para estimar el área de vida de Mashca ya que los grados de libertad del área son 419.78.

### AKDE

Para estimar la utilización de distribución de la primera mitad de datos de mashca:

```{r eval = FALSE}
mashca_first_akde <- akde(mashca_merged[[1]], fits_first,
                    trace = 2)
summary(mashca_first_akde, level.UD = 0.95)
```
```{r include = FALSE}
mashca_first_akde <- akde(mashca_merged[[1]], fits_first,
                    trace = 2)
mashca_second_akde <- akde(mashca_merged[[2]], fits_second,
                    trace = 2)
```
:::::: {.callout-note}
## Ejercicio
Calcula el AKDE para la segunda mitad de los datos de mashca. ¿Hay alguna diferencia?
:::

Finalmente, puedes realizar un plot de estas áreas.

```{r}
par(mfrow=c(2,2))
plot(mashca_merged[[1]], UD = mashca_first_akde,
      col.grid=NA, level.UD= 0.5, error=2, main = "Primera mitad AKDE core mashca")
plot(mashca_merged[[2]], UD = mashca_second_akde,
      col.grid=NA, level.UD= 0.5, error=2, main = "Segunda mitad AKDE core mashca")

plot(mashca_merged[[1]], UD = mashca_first_akde,
      col.grid=NA, level.UD= 0.95, error=2, main = "Primera mitad AKDE core mashca")

plot(mashca_merged[[2]], UD = mashca_second_akde,
      col.grid=NA, level.UD= 0.95, error=2, main = "Segunda mitad AKDE core mashca")
```

### AKDE vs KDE

Realiza una comparación de ambos métodos y analiza sus semejanzas y diferencias. ¿Es uno de estos métodos mejor que el otro?

```{r include = FALSE}
mashca.first.proj <- SpatialPointsDataFrame(coords = as.data.frame(cbind
                        (mashca_first$location.long, mashca_first$location.lat)),
                        data = mashca_first, proj4string = CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

mashca.first.transf <- spTransform(mashca.first.proj,
                           CRS("+proj=utm +south +zone=17 +ellps=WGS84"))

mashca.first.khref <- kernelUD(mashca.first.transf, h = "href")
kernel.area(mashca.first.khref, percent = seq (50, 95, 15),
            unin="m", unout="km")

# Obten los poligonos
mashca.first.KDE95 <- getverticeshr(mashca.first.khref, percent = 95)
mashca.first.KDE50 <- getverticeshr(mashca.first.khref, percent = 50)
```

```{r include = FALSE, eval = FALSE}
# ctmm package meses up with the base R plot function for KDE
# I tried to unload it but is not working either
par(mfrow=c(2,2))
plot(mashca_merged[[1]], UD = mashca_first_akde,
      col.grid=NA, level.UD= 0.5, error=2, main = "Primera mitad AKDE 50% mashca")
plot.default(mashca.first.KDE50, col = "grey", axes = TRUE, main = "Primera mitad KDE 50% mashca")
plot(mashca.first.transf, col = "black", add = TRUE, pch = 21, cex = 0.2)

plot(mashca_merged[[1]], UD = mashca_first_akde,
      col.grid=NA, level.UD= 0.95, error=2, main = "Primera mitad AKDE 95% mashca")
plot(mashca.first.KDE95, col = "grey", axes = TRUE, main = "Primera mitad AKDE 95% mashca")
plot(mashca.first.transf, col = "black", add = TRUE, pch = 21, cex = 0.2)
```

![Comparación de dos métodos para estimar áreas de vida, AKDE vs KDE de la primera mitad de datos del zorro andino  Mashca](docs/Figures/akde_vs_kde.png)

:::::: {.callout-note}
## Tarea
Investiga como realicé los gráficos de las áreas de vida y núcleo estimadas con AKDE y KDE.
:::

## _moveVis_

Una excelente forma de atraer la atención a tus proyectos y de alcanzar a una audiencia no académica, es mostrando como el animal se desplaza por su área de vida sobre mapas satelitales. Es con este proposito que fue diseñado [moveVis](https://movevis.org/index.html).

Ya que _moveVis_ ha pasado por algunos cambios que no han sido cargados al repositorio de CRAN y este ha sido removido, te recomiendo descargar la versión 0.10.5 [desde aquí](https://cran.r-project.org/src/contrib/Archive/moveVis/). Este modo de instalación también puede requerir en algunos casos que instales ciertas dependencias primero por lo que te recomiendo que los instales de esta manera.

```{r eval = FALSE}
library(pacman)
p_load(move, slippymath, magick, gifski, av, pbapply, moveVis)
```

Posteriormente, dirigete a la pestaña de _Packages_ en R Studio y elige la opción _Install_ en la esquina superior izquierda, y elige  instalar desde un archivo `.tar.gz`. Ahora selecciona el directorio donde lo descargaste e instalalo.

Accede a Movebank y descarga el set de datos que utilizarás hoy para animar los movimientos de un jaguar muestreado en la Caatinga en Brasil. Para esto, ingresa al mapa y escribe _"jaguar conservation in the caatinga biome"_, elige el estudio y haz click en "Open in studies page". Una vez aparezca la ventana que contiene todo acerca de este estudio, haz clic en "Download" y descarga este archivo en formato `csv`.

```{r include = FALSE}
library(moveVis)
library(move)

jaguars <- read_csv("docs/datos/Jaguars.csv")

# Selecciona a Courisco
courisco <- jaguars %>%
  filter(`individual-local-identifier` == "Courisco")
```

Debido a la gran cantidad de datos colectados en este estudio, selecciona solo a uno de los jaguares.
```{r eval = FALSE}
library(moveVis)
library(move)
library(tidyverse)

# Importa los datos
jaguars <- read_csv("Jaguars.csv")

# Selecciona a Courisco
courisco <- jaguars %>%
  filter(`individual-local-identifier` == "Courisco")
```

Ahora deberás llevar a cabo unos cuantos pasos para poder animar los movimientos de Courisco.

```{r}
# Transformar en objeto move y selecciona una porcion de los datos
courisco_move <- df2move(courisco[1:500,],
             proj = "+proj=longlat +south +zone=23 +ellps=WGS84",
             x = "location-long", y = "location-lat", time = "study-local-timestamp",
             track_id = "individual-local-identifier")
```

Es importante saber cual es la frecuencia de muestreo antes de obtener los frames de la animacion.

```{r}
# Visualizar frecuencia de muestreo
lagging <- timeLag(courisco_move,
                  unit="mins")
summary(lagging)

# Histograma de frecuencias de muestreo
hist(lagging, xlab= "mins", xlim= c(0, 400), breaks= 40)
```

Utiliza el número estimado anteriormente para alinear y estandarizar la frecuencia de muestreo sobre el objeto.

```{r}
courisco_move <- align_move(courisco_move, res = 60, unit = "mins")
```

Hay una variedad de mapas que puedes utilizar como fondo para tu animación.
```{r}
##ver mapas disponibles
get_maptypes()
```

También puedes indicar cuantos núcleos quieres utilizar para acelerar el proceso.

```{r eval = FALSE}
##utilizar nucleos
use_multicore()
use_multicore(n_cores = numero)
```

Finalmente, utiliza el objeto alineado y estandarizado para obtener los frames.
```{r eval = FALSE}
# Crear el objeto que alamcene los frames de la animacion
frames <- frames_spatial(courisco_move,
                         map_service = "osm",
                         map_type = "topographic",
                         path_legend = F,
                         path_colours = "black",
                         path_alpha = 0.8,
                         tail_size = 2,
                         tail_length = 5)

frames[[535]]
```

Puedes editar el formato de estas imágenes siguiendo el estilo de escritura de _ggplot2_.

```{r eval = FALSE}
frames_edit <- add_labels(frames, x = "Longitud" , y = "Latitud", verbose = TRUE) %>%
  add_scalebar(height = 0.015) %>%
  add_northarrow() %>%
  add_timestamps(courisco_move, type = "label") %>%
  add_gg(frames, gg = expr(labs(title = "Movimientos del jaguar Courisco en Brasil")))

frames_edit[[535]]
```

Puedes mirar varios formatos para exportar la animación utilizando `suggest_formats()`.

```{r eval = FALSE}
animate_frames(frames_edit, out_file = "courisco.gif",
              height = 500, width = 500, res = 82, overwrite = TRUE)
```

![Animación de 4500 puntos GPS del jaguar Courisco en la Caatinga brasilera ](docs/Figures/courisco.gif){#fig-courisco fig-alt="movimientos del jaguar Courisco en _movevis_"}

## Fuentes

Puedes utilizar el [manual de _adehabitatHR_](https://cran.r-project.org/web/packages/adehabitatHR/vignettes/adehabitatHR.pdf) para resolver algunos ejercicios de este tutorial, y puedes explorar más funciones que también ofrece este paquete