Introduction

Example: White-bearded wildebeest (2010-2012)

Today, we will use part of a dataset collected during Jared Stabach’s field research in Kenya from 2010-2012. This dataset is freely available on Movebank and consist of GPS tracking data from 12 adult white-bearded wildebeest (Connochaetes taurinus). The goal of this exercise is to demonstrate how to work with animal movement data and combine it with various types of remote sensing data.

In this exercise, you will:

  • Import and clean an animal movement trajectory
  • Annotate animal movement data using existing environmental data from remote sensing, including land cover, elevation, and NDVI (Normalized Difference Vegetation Index)
  • Annotate animal movement with customized environmental information. In this practice, we will calculate the distance to the nearest urban area and various terrain indices as an example.

References

You can refer to these data in your research using the following citations:

Stabach JA, Hughey LF, Crego RD, Fleming CH, Hopcraft JGC, Leimgruber P, Morrison TA, Ogutu JO, Reid RS, Worden JS, Boone RB. 2022. Increasing anthropogenic disturbance restricts wildebeest movement across East African grazing systems. Frontiers in Ecology and Evolution. doi.org/10.3389/fevo.2022.846171

Stabach JA, Hughey LF, Reid RS, Worden JS, Leimgruber P, Boone RB. 2020. Data from: Comparison of movement strategies of three populations of white-bearded wildebeest. Movebank Data Repository. doi:10.5441/001/1.h0t27719


Setup

We have already covered data cleaning in previous sessions of our workshop. Therefore, we will start by loading a pre-processed white-beared wildebeest dataset. This will allow us to focus on the integration and analysis of remote sensing data.

Load Libraries

First, we need to load the required libearies and clear any existing data from R’s memory. This ensures that we start with a clean workspace and have all the necessary tools for our analysis.

# Clear all objects from the workspace to ensure a clean environment
rm(list=ls())

# Load required libraries
library(sf)         # Simple features for R
library(rgdal)      # Geospatial data abstraction library
library(sp)         # Spatial data classes
library(rgeos)      # Interface to geometry engine
library(dplyr)      # Data manipulation
library(raster)     # Raster data manipulation
library(viridis)    # Color scales
library(ggplot2)    # Data visualization
library(mapview)    # Interactive maps
library(terra)      # Raster and vector data manipulation
library(lubridate)  # Date and time manipulation

##########################################
#Version 2: only use sf and terra package

# Load required libraries
#library(sf)         # Simple features for R
#library(dplyr)      # Data manipulation
#library(terra)      # Raster and vector data manipulation
#library(ggplot2)    # Data visualization
#library(mapview)    # Interactive maps
#library(lubridate)  # Date and time manipulation

Animal Data

Next, we will load the animal movement data, filter out any rows with missing coordinates and set the coordinate reference system (CRS) to UTM zone 37S. This ensures that all spatial data is in the same projection, making it easier to integrate and analyze.

# Load the pre-processed animal movement data from an RDS file
WB <- readRDS("wildebeest_3hr_adehabitat.rds")

# Check the basic structure and class of the data to understand its format
class(WB) 
## [1] "data.frame"
str(WB)
## 'data.frame':    52847 obs. of  17 variables:
##  $ x        : num  281442 282340 281930 282040 282187 ...
##  $ y        : num  9808391 9809132 9808660 9808734 9808773 ...
##  $ date     : POSIXct, format: "2010-10-19 03:00:00" "2010-10-19 06:00:00" ...
##  $ dx       : num  898.3 -410.7 109.9 147.7 -23.7 ...
##  $ dy       : num  740.4 -471.6 74.1 38.3 97.9 ...
##  $ dist     : num  1164 625 133 153 101 ...
##  $ dt       : num  10800 10800 10800 10800 10800 10800 10800 10800 10800 10800 ...
##  $ R2n      : num  0 1355089 310028 474624 700726 ...
##  $ abs.angle: num  0.689 -2.287 0.593 0.254 1.808 ...
##  $ rel.angle: num  NA -2.977 2.88 -0.339 1.555 ...
##  $ id       : Factor w/ 12 levels "Karbolo","Kikaya",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ burst    : chr  "Karbolo" "Karbolo" "Karbolo" "Karbolo" ...
##  $ sex      : Factor w/ 2 levels "1","2": 2 2 2 2 2 2 2 2 2 2 ...
##  $ DOP      : num  3.4 2.8 1.8 1.6 1.8 2.2 4.4 3.6 1.6 3.2 ...
##  $ fixType  : num  3 3 3 3 3 3 3 3 3 3 ...
##  $ temp     : num  17 17 25 30 28 25 22 18 15 17 ...
##  $ speed    : num  0.10779 0.0579 0.01227 0.01413 0.00933 ...
# Filter out rows with missing coordinates to ensure all data points are valid
WB <- WB %>%
  filter(!is.na(x) & !is.na(y))

#Convert the data frame to an sf object and set the CRS to UTM zone 37 S
WB_utm <- WB %>%
  st_as_sf(coords = c("x", "y"), crs = 32737)

Remote Sensing Data

To analyze the environment where the animal move, we will use remote sensing data. We can download this data from NASA AppEEARS platform, which requires creating an area sample to upload. Note that you need to register for an Earthdata account and log in to the website.

Create An Area Sample

Firstly, we need to create a minimum convex polygon (MCP) that encompasses the area where the animals move. This MCP will serve as the area sample to upload. The projection must be WGS 84.

# Convert the projection to WGS 84 (EPSG:4326)
WB_84 <- st_transform(WB_utm, crs = 4326)

# Create a minimum convex polygon (MCP) to cover the moving areas of all individuals
WB_84_mcp <- st_convex_hull(st_union(WB_84))

# Save the MCP as a shapefile
st_write(WB_84_mcp, "Boundary.shp")
## Writing layer `Boundary' to data source `Boundary.shp' using driver `ESRI Shapefile'
## Writing 1 features with 0 fields and geometry type Polygon.
# Load the MCP shapefile to confirm is was saved correctly
MCP=st_read("Boundary.shp")
## Reading layer `Boundary' from data source 
##   `/Users/Qianr/Documents/GitHub/Movement_workshop/Movement_workshop/animove-2023-remote-sensing-main/Boundary.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 1 field
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: 36.78754 ymin: -1.865096 xmax: 37.11519 ymax: -1.338212
## Geodetic CRS:  WGS 84

Download Data

Next, we need to combine all the “Boundary” files into a single zip file and upload it to NASA’s AppEEARS platform.


Step 1

To start, navigate to the AppEEARS website, click “Extract” in the top left corner, and select “Area”. On the next page, choose “Start a new request”.

Fig 1.
Fig 1.

Step 2

From the top to the bottom, set a name for your project, upload the zip file that contains Boundary.shp, Boundary.dbf, Boundary.prj, and Boundary.shx files (the shapefile documents we just created), and select the start date as 10-15-2010, and the end date as 12-04-2012.

Fig 2.
Fig 2.

Step 3

Next, we will select the layers we need for our analysis:

  1. Landcover Type
  2. Elevation Data
  3. NDVI (Normalized Difference Vegetation Index)
    • Type “MOD13Q1.061” or “Terra MODIS Vegetation Indices (NDVI & EVI)” - “_250M_16_days_NDVI”

Fig 3.
Fig 3.
Fig 4.
Fig 4.

Fig 5.

Step 4

After selecting the environmental layers, we need to choose the format and projection for the output. For the format, choose “GeoTiff”, and for the projection, choose “Geographic”, which is WGS 84. Finnaly, click “Submit”.

Fig 6.
Fig 6.

Step 5

Once submitted, you will receive a confirmation email indicating that the data is ready for download. After downloading, create three folders (Landcover, Elevation and NDVI) to store the data downloaded from NASA AppEEARS, as shown in the following picture.

Fig 7.
Fig 7.

Land Cover Data

We used the Terra and Aqua combined Moderate Resolution Imaging Spectroradiometer (MODIS) Land Cover Type (MCD12Q1) Version 6.1 data product. It utilized data from both Terra and Aqua satellites and provide global land cover types at yearly intervals (2001-2022). The land cover types are derived from the International Geosphere-Biosphere Programme (IGBP). And the spatial resolution is 500 meters.

In this section, we will visualize the land cover types and calculate the distance from each GPS location to the nearest urban areas.

Visualize Land Cover Types

We will load, stack, rename each layers of the land cover data according to the Annual International Geosphere-Biosphere Programme (IGBP) classification, and then visualize the land cover types for 2010.