edwardlavender / patter

Particle filters, smoothers and sampling algorithms for animal movement modelling, with a focus on passive acoustic telemetry systems.
https://edwardlavender.github.io/patter/
Other
3 stars 1 forks source link
biologging biotelemetry electronic-tagging-and-tracking flapper-algorithms julia movement-modelling particle-filter particle-smoother passive-acoustic-telemetry patter r r-package simulation two-filter-smoother

patter: particle algorithms for animal movement

Particle filters, smoothers and sampling algorithms for animal movement modelling in R

Project Status: Active – The project has reached a stable, usable
state and is being actively
developed. Lifecycle:
experimental CRAN
status Coverage R-CMD-check

patter provides particle filtering, smoothing and sampling algorithms for animal movement modelling, with a focus on passive acoustic telemetry systems. This wraps and enhances a fast Julia backend (Patter.jl). The methodology enables the reconstruction of movement paths and patterns of space use. patter unifies a suite of methods formerly known as the flapper algorithms and supersedes the experimental flapper package (Lavender et al., 2023).

Note: patter is a new R package. Like all new packages, you should use it with a degree of caution. Please share feedback and issues.

Highlights

patter is designed to reconstruct movement paths and emergent patterns of space use from animal tracking data. A powerful, flexible, process-orientated, particle-based framework is used for this purpose. This framework unifies the flapper algorithms and provides important opportunities for development, which we exploit here.

The essential functions are pf_filter() and pf_smoother_*():

We hope to add backward sampling algorithms to the package in due course.

Evolution

patter evolved from the experimental flapper package, but is:

See NEWS for a summary of the evolution of flapper to patter.

At the time of writing (May 2024), patter is more streamlined than flapper and focuses on the implementation of fast particle-based algorithms for the reconstruction of movements and patterns of space use. Please get in touch if you would like to see additional functionality brought into patter.

Installation

Note: patter currently works on Windows and MacOS. On Windows, everything should work if you follow the instructions below. On MacOS, some additional set up (such as compiler configuration) may be required, depending on your set up. In our (limited) experience, patter installs but crashes on Debian/Ubuntu. This is due to a conflict between the GDAL/GEOS/PROJ libraries used by R and Julia (which we hope to solve in due course). Please let us know your experiences if you are using other Linux distributions. In case of issues, you should be able to use Patter.jl directly, which on some systems may be simpler than getting R and Julia to play together!

  1. Install R. This package requires R version ≥ 4.1 (but the most recent version is recommended). You can check your version from the R console using R.version.string.

  2. Install build packages. Package installation and configuration (may) require the devtools, pkgbuild and here packages. Install them with:

install.packages(c("devtools", "pkgbuild", "here"))
  1. Install system libraries.
  1. Install Julia. Julia is high-performance programming language that patter uses as a backend. If you do not have Julia installed on your system, you have three options.

Note: Install a recent Julia version. Julia v1.10.5 (long-term release) is recommended. Julia v1.11 is not currently supported by JuliaCall (see here).

A. Install Julia via R using JuliaCall:

# Install the {JuliaCall} package:
install.packages("JuliaCall")

# Use the development version if the CRAN version is unavailable:
devtools::install_github("Non-Contradiction/JuliaCall",
                         dependencies = TRUE)
# Install `Julia` via {JuliaCall}:
# * Set `JULIA_HOME` if Julia is not found (see `?julia_setup()`)
library(JuliaCall)
julia <- julia_setup(installJulia = TRUE)
# Validate the Julia installation:
# * TRUE: `Julia` is working! 
# * FALSE: `Julia` is not working (see below)! 
isTRUE(try(julia_eval('true'), silent = TRUE))

If julia_setup() fails with 'Julia is not found', you should tell R the location of the Julia binary via JULIA_HOME (see ?JuliaCall::julia_setup() and the JuliaCall README for troubleshooting and ways to get help).

B. If JuliaCall installation fails, download and install Julia from JuliaLang. On MacOS, make sure you download a version appropriate for your architecture (Intel or ARM). Then retry julia_setup():

# Set up Julia
# * Set `JULIA_HOME` if Julia is not found (see `?julia_setup()`)
julia <- julia_setup()
isTRUE(try(julia_eval('true'), silent = TRUE))

C. The other option is to use juliaup. Some users have found this easier on MacOS because you don’t have to worry about finding the right Julia installation for your architecture:

juliaup add 1.10.5  
juliaup default 1.10.5

Please share the solution(s) that worked for you and help to improve these instructions.

  1. Install patter. To install patter from the main branch, use:
devtools::install_github("edwardlavender/patter", 
                         dependencies = TRUE, 
                         build_vignettes = TRUE)

The dependencies = TRUE argument ensures that suggested packages are also installed, which are required for some functions and to build vignettes. This process may take several minutes. Set build_vignettes = FALSE for a faster installation.

To install patter from the development (dev) branch, if available, use:

devtools::install_github("edwardlavender/patter@dev", 
                         dependencies = TRUE, 
                         build_vignettes = TRUE)

This branch may include bug fixes and new features but should be used with caution.

We recommend using renv (or similar) and RStudio Projects to track the version of patter that you use in your projects. This will ensure that your code continues to work, even if we have to make breaking changes to patter as the package evolves in response to user feedback.

  1. Connect to Julia. At the start of every R session, you need to connect R to Julia (and patter to Patter.jl):
# Load & attach {patter}:
library(patter)

# Option (A): Connect to `Julia`: 
# * Set JULIA_HOME if 'Julia not found'
# * Set JULIA_PROJ to use a local Julia project (recommended)
# * Set JULIA_NUM_THREADS to exploit multi-threading (recommended)
# * See `julia_connect()` for guidance
julia <- julia_connect()

The first time you run julia_connect(), it will connect to Julia and install (and pre-compile) Patter.jl and the additional Julia dependencies. This may take a few minutes. Subsequent julia_connect() calls will be faster.

  1. Validate the RJulia connection. To validate that patter works on your system, run:
julia_validate()

This should return NULL, invisibly, in which case you are good to go. Otherwise, the function will return an error (or R may crash). Please report any issues you experience during this process.

Functionality

Vignettes

For an introduction to patter, use:

For a full list of all functions, see help(package = 'patter').

For a glossary of key arguments, see glossary.

Datasets

For example datasets from the Movement Ecology of Flapper Skate project (datasets-mefs), which inspired patter, see:

To validate new datasets for use with patter, see pat_setup_data() and/or the assemble_*() function documentation.

For example algorithm outputs (datasets-algorithms), see:

Set up Julia

To link patter and the Patter.jl Julia backend, use:

These functions should be run at the start of every R session.

Abstract Types

patter is based on three Abstract Types, defined in Julia:

Simulation

To simulate animal movement time series, see:

To evaluate model skill in reconstructing simulated patterns, see skill_*() functions:

Data exploration

For help with data acquisition, processing, checking and preliminary analyses, see the flapper package. This facilitates:

Please submit a feature request if you would like functions from flapper in patter.

Algorithms

The main thrust of patter is the provision of fast, integrated modelling workflow based on particle filtering for reconstructing animal movement paths and emergent patterns of space use from observational time series (with a focus on passive acoustic telemetry systems).

To assemble datasets for particle filtering, use assemble_*() functions:

Ancillary time series should be structured in the same way for inclusion in the particle filter.

To implement particle filtering (PF) routines, use:

These functions return pf_particles-class objects.

For convenience plotting functions, see:

For mapping utilisation distributions, use:

Options

For additional options in patter, see:

Usage

Set up

This is the basic patter workflow to reconstruct movement paths and patterns of space use from animal tracking data. First, we load some essential packages:

library(patter)
#> This is {patter} v.1.0.1. For an overview, see `?patter`. For support, contact edward.lavender@eawag.ch.
library(data.table)
#> Warning: package 'data.table' was built under R version 4.3.3
library(dtplyr)
library(dplyr, warn.conflicts = FALSE)
options(patter.verbose = FALSE)

Second, we connect R to Julia and set the seed in R and Julia to ensure reproducibility of our simulations:

julia_connect()
julia_validate()
set_seed()

Third, we define the properties of our study area; namely, a SpatRaster of our study area that defines the area within which movements are possible and the timeline over which we will model movements:

# Define map 
map <- dat_gebco()
set_map(map)

# Define timeline 
timeline <- seq(as.POSIXct("2016-03-17 01:50:00", tz = "UTC"),
                as.POSIXct("2016-03-18 01:48:00", tz = "UTC"), 
                by = "2 mins")

Movement

We will reconstruct the movements of a tagged flapper skate (Dipturus intermedius) within a study area off the west coast of Scotland, based on electronic tagging and tracking data. To do so, we need a model for the individual’s movements and a series of observation models that connect movements to observations. In this example, we are interested in the two-dimensional (x, y) location of our animal through time (that is, the animal’s ‘state’ is an object of type StateXY). The animal can move up to 750 m in two minutes, which is the resolution at which we will model movement, and we formulate a random walk model accordingly based on step lengths and turning angles:

# Define the animal's state:
state      <- "StateXY"

# Formulate a corresponding movement model:
mobility   <- 750.0
model_move <- move_xy(dbn_length = glue::glue("truncated(Gamma(1, 250.0), upper = {mobility})"),
                      dbn_angle = "Uniform(-pi, pi)")

# Visualise realisations of the movement model:
map |> 
  sim_path_walk(.timeline = timeline,
                .state = state,
                .model_move = model_move, 
                .n_path = 4L, .one_page = TRUE) |> 
  invisible()

Observations

We have collected acoustic and archival (depth) observations from tagged flapper skate. Let’s pull out the time series for a selected individual:

# Define acoustic detections
acc <-
  dat_acoustics |>
  filter(individual_id == 25L) |>
  mutate(individual_id = NULL) |>
  as.data.table()

# Define archival (depth) observations
arc <-
  dat_archival |>
  filter(individual_id == 25L) |>
  mutate(individual_id = NULL,
         depth_sigma = 50,
         depth_deep_eps = 30) |>
  rename(obs = depth) |>
  as.data.table()

Individual movements are connected to the observations by models of the observation process for each dataset. Without going into details, here we bundle together the observations with the parameters of the observation models:

model_1   <- "ModelObsAcousticLogisTrunc"
acoustics <- assemble_acoustics(.timeline = timeline,
                                .acoustics = acc,
                                .moorings = dat_moorings)

model_2  <- "ModelObsDepthNormalTrunc"
archival <- assemble_archival(.timeline = timeline,
                              .archival = arc)

Of course, you do not need acoustic and archival data to implement the algorithms (these are just the data we have collected from flapper skate)—other datasets can be used just as easily. To simulate observations instead, see sim_observations().

Particle filter

We are now in a position to run the particle filter. This runs a simulation forwards (or backwards) in time, sampling states (locations, termed ‘particles’) that are consistent with the movement model and the observations up to and including each time point. We end up with a time series (data.table) of particles that approximate the partial marginal distribution for the location of the animal, at each time step:

# List filter arguments
args <- list(.map = map,
             .timeline = timeline,
             .state = state,
             .xinit_pars = list(mobility = mobility),
             .yobs = list(acoustics, archival),
             .model_obs = c(model_1, model_2),
             .model_move = model_move,
             .n_record = 500L,
             .n_particle = 1e5L)

# Forward run
fwd <- do.call(pf_filter, args, quote = TRUE)
head(fwd$states)
#>    path_id timestep           timestamp map_value        x       y
#>      <int>    <int>              <POSc>     <num>    <num>   <num>
#> 1:       1        1 2016-03-17 01:50:00  59.76520 709142.1 6253007
#> 2:       1        2 2016-03-17 01:52:00  68.53316 709276.5 6253291
#> 3:       1        3 2016-03-17 01:54:00  45.86026 709476.1 6252964
#> 4:       1        4 2016-03-17 01:56:00  44.46762 709390.0 6252794
#> 5:       1        5 2016-03-17 01:58:00  60.64737 708976.4 6252849
#> 6:       1        6 2016-03-17 02:00:00  55.42853 709437.9 6253395

# Backward run
args$.direction <- "backward"
bwd <- do.call(pf_filter, args, quote = TRUE)

Particle smoother

Particle smoothers refine the outputs from the particle filter. Smoothed particles approximate the full marginal distribution for the location of the individual at each time step (accounting for all of the data before and after each step).

smo <- pf_smoother_two_filter(.n_particle = 100L, .n_sim = 100L)

Mapping

Particles can be used to reconstruct movement paths and patterns of space use. We can estimate a utilisation distribution from our particle samples as follows:

# Estimate UD
ud <- map_dens(.map = map,
               .coord = smo$states,
               sigma = spatstat.explore::bw.diggle)$ud
#> Observation window is gridded.

# Add home range
map_hr_home(ud, .add = TRUE)
mtext(side = 4, "Probability density", line = -3)

This basic workflow is highly customisable. You have the flexibility to define species-specific movement models, include any type of observational dataset and implement system-specific observation models. See the vignettes and function examples for further details and reach out with queries.

Resources

For full details on the methods, see the references below.

For further information of the patter package, see:

For further code examples, see:

Disclaimer and troubleshooting

patter is a new R package. All routines are experimental. Researchers interested in using the package are encouraged to get in touch while the methods and package remain at an early stage of evolution (edward.lavender@eawag.ch).

Citation

To cite patter in publications, please use:


Please note that patter is released with a Contributor Code of Conduct. By contributing to this project, you agree to abide by its terms.