anipaths: spline-based interpolation

Henry Scharf

2024-02-01

Animating animal trajectories

The package anipaths contains a collection of telemetry observations for turkey vultures originally analyzed in:

Dodge S, Bohrer G, Bildstein K, Davidson SC, Weinzierl R, Mechard MJ, Barber D, Kays R, Brandes D, Han J (2014) Environmental drivers of variability in the movement ecology of turkey vultures (Cathartes aura) in North and South America. Philosophical Transactions of the Royal Society B 20130195.

To animate the locations, we first need to create a time stamp variable of class numeric or POSIX. One advantage to using the POSIX class is that we can specify the gaps in the interpolation (delta.t) using convenient character strings like "hour" or "week".

library(anipaths)
vultures$POSIX <- as.POSIXct(vultures$timestamp, tz = "UTC")
vultures_paths <- vultures[format(vultures$POSIX, "%Y") == 2009, ] ## limit attention to 2009
animate_paths(paths = vultures_paths, 
              delta.t = "day",
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier")

Using background maps

Often it may be useful to add a map of the relevant study area. There are lots of ways to incorporate a map. The simplest way is to set background to TRUE, in which case anipaths will do the best it can to select a map based on the data. In the next example, we’ve changed the time step to help the animations load a little faster.

library(ggmap)
animate_paths(paths = vultures_paths, 
              delta.t = 2 * 24 * 60 * 60, ## number of seconds in two days
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier", 
              background = TRUE, img.name = "background_TRUE")

You can also give a long/lat location, zoom level (3-21; see ?ggmap::get_map()), and maptype (satellite, terrain, hybrid) to be passed to ggmap::get_map(), and anipaths will make a background for you.

As shown below, the value of delta.t can be specified as a numeric, which will be interpreted in whatever units used by as.numeric(paths['Time.name']) (in our case, this is seconds).

vultures_paths <- vultures[format(vultures$POSIX, "%Y") == 2009:2010, ]
background <- list(center = c(-90, 10), 
                   zoom = 3, 
                   maptype = "satellite")
animate_paths(paths = vultures_paths, 
              delta.t = 3 * 24 * 60 * 60, ## number of seconds in three days
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier", 
              background = background, img.name = "background_google")

You can also supply your own background image. The projection and units should match the data.

background <- geodata::world(path = ".")
sf::st_crs(background)$proj4string ## matches default projection in animate_paths()
animate_paths(paths = vultures_paths, 
              delta.t = "week",
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier", 
              background = background, img.name = "background_user")
system("rm -r gadm") ## remove geodata map from machine

The function animate_paths() must project the data to match Google’s map tiles, so if your data aren’t in long/lat format, make sure you update the paths.proj variable with the correct string.

Adding individual-level covariate information

If you have a covariate of interest for each individual, animate_paths() can display that information as a colored ring around each individual. We don’t have any natural individual-level covariates available for the vultures, so we’ll make one up for the purposed of demonstration. The following code assigns each individual a random interval for each of three behaviors: exploratory, directed, and stationary.

behaviors <- c("exploratory", "directed", "stationary")
set.seed(1)
vultures_paths$behavior <- 
  unlist(sapply(unique(vultures_paths$individual.local.identifier), function(id){
    v_id <- vultures_paths[vultures_paths$individual.local.identifier == id, ]
    switches <- c(0, sort(sample(1:nrow(v_id), 2)), nrow(v_id))
    rep(behaviors[sample(1:3, 3)], diff(switches))
  }))

The covariates have now been appended to the paths data frame. We can let animate_paths() know we would like it to display this information by setting the argument covariate to match the name of the appropriate column in paths (i.e., behavior). The default colors are a gray scale. Any collection of colors can be provided (e.g., covariate.colors = viridis::viridis(3)) which will be turned into a palette to match the support of the covariate.

delta.t <- "day"
background <- geodata::world(path = ".")
sf::st_crs(background)$proj4string ## matches default projection in animate_paths()
animate_paths(paths = vultures_paths, 
              delta.t = delta.t,
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              covariate = "behavior", covariate.colors = viridis::viridis(3),
              ID.name = "individual.local.identifier",
              background = background, img.name = "covariates")
system("rm -r gadm") ## remove geodata map from machine

Note: if a covariate is represented numerically but should be treated as a factor, changing the class using as.factor() should prevent animate_paths() from interpolating to “intermediate” values (e.g., if male and female are coded as 0 and 1, defining the column of the data frame as a factor will constrain individuals to only take on values of 0 and 1).

Output formats

The default format for displaying the animation is an index.html file. This is great for broad functionality, and the generated images may be used as the basis of a variety of other animation formats (e.g., the “animate” package in LaTeX). If you have ffmpeg installed on your system, you can also set method = "mp4" to create a stand-alone video that is easy to share with others. For more information, see https://www.ffmpeg.org/.

Checking the interpolation

As a way to check that anipaths is producing reasonable interpolations of the telemetry observations, a generic plot() function is provided that takes an argument of class paths_animation produced by calling animate_paths() with return.paths = TRUE. See example(plot.paths_animation). Adjusting the parameters of the interpolation can be done by modifying the s_args argument.

The following code demonstrates the use of s_args by intentionally over-smoothing the trajectories.

interp <- animate_paths(paths = vultures_paths, 
                        delta.t = "day",
                        coord = c("location.long", "location.lat"),
                        Time.name = "POSIX",
                        ID.name = "individual.local.identifier",
                        s_args = rep(list(list(k = 10)), 10),
                        return.paths = T)
plot(interp, i = 2)

This piece of code shows the number of observations for each individual vulture, in the same order they are handled by animate_paths() (i.e., the same order as unique(vultures_paths$individual.local.identifier)). The wide variability in duration and number of observations required individual-specific numbers of knots in the GAM interpolation. In this case, the number of observations divided by 4 is used as the number of knots for each individual, with a maximum of 306 knots because each frame in the animation corresponds to a single day between March 1 and December 31.

obs_counts <- merge(data.frame(individual.local.identifier = unique(vultures_paths$individual.local.identifier)), 
                    aggregate(timestamp ~ individual.local.identifier, data = vultures_paths, FUN = length),
                    by = "individual.local.identifier", sort = F)
obs_counts
s_args <- lapply(obs_counts$timestamp, function(x) c(k = floor(min(x / 4, 306))))
interp <- animate_paths(paths = vultures_paths, 
                        delta.t = "day",
                        coord = c("location.long", "location.lat"),
                        Time.name = "POSIX",
                        ID.name = "individual.local.identifier",
                        s_args = s_args,
                        return.paths = T)
plot(interp, i = 2)

As of version 0.10.2, the default value for the number of knots also works well for this example.

interp <- animate_paths(paths = vultures_paths, 
                        delta.t = "day",
                        coord = c("location.long", "location.lat"),
                        Time.name = "POSIX",
                        ID.name = "individual.local.identifier",
                        return.paths = T, verbose = T)
plot(interp, i = 2)

Additional features

Showing complete trajectories

Set whole.path = TRUE. It is not required, but probably clearer visually, if we also set tail.length = 0.

animate_paths(paths = vultures_paths, 
              delta.t = delta.t,
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier", 
              background = background, 
              whole.path = TRUE, tail.length = 0, img.name = "whole_traj")

Dimming selected individuals

vultures_paths <- vultures[format(vultures$POSIX, "%Y") == 2009 & 
                             vultures$location.lat > 32, ]
animate_paths(paths = vultures_paths, 
              delta.t = delta.t,
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier", 
              background = background, dimmed = c(1, 3, 5), img.name = "dim")

Showing relational information (beta)

It is also possible to visualize dynamic pair-wise relational information with anipaths. This information needs to be coded as an array of adjacency matrices (binary and weighted edges allowed, no support yet for directed information). If the network information is on a different time scale than the position information, it will be interpolated using smoothing splines to match.

We don’t have existing relationships among the vultures to display, so we’ll make some up for the purposes of demonstration.

vultures_paths <- vultures[format(vultures$POSIX, "%Y") == 2009, ]
set.seed(1)
n_indiv <- length(unique(vultures_paths$individual.local.identifier))
change_pts <- 5
network.times <- seq(min(vultures_paths$POSIX), max(vultures_paths$POSIX), l = change_pts + 2)
network <- array(NA, dim = c(n_indiv, n_indiv, length(network.times)))
for(time_i in 2:length(network.times)){
  network_mat <- matrix(sample(1:0, n_indiv^2, prob = c(0.1, 0.9), replace = T), 
                        n_indiv, n_indiv)
  network_mat[lower.tri(network_mat)] <- t(network_mat)[lower.tri(network_mat)]
  diag(network_mat) <- 1
  network[, , (time_i - 1):time_i] <- network_mat
}
delta.t <- 3600*24*2
animate_paths(paths = vultures_paths, 
              delta.t = delta.t,
              coord = c("location.long", "location.lat"),
              Time.name = "POSIX",
              ID.name = "individual.local.identifier",
              background = background, network = network, network.times = network.times,
              img.name = "network")
system("rm -r js; rm -r css; rm -r images; rm index.html")