Studying the complex ecological systems across both space and time is imperative in understanding the full dynamics of species. This vignette illustrates the construction of a spatiotemporal ISDM using PointedSDMs, using data of species Colinus virginianus across Alabama (United States of America). The first step in this vignette is to load the required packages:
library(PointedSDMs)
library(inlabru)
library(ggplot2)
library(spocc)
library(INLA)
library(dplyr)
library(sp)
library(sf)as well as define some objects required by the model to run.
proj <- "+proj=utm +zone=17 +datum=WGS84 +units=km"#"+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
AL <- USAboundaries::us_states(states = "Alabama", resolution = 'high')
AL <- as(AL, "sf")
AL <- st_transform(AL, proj)
mesh <- inla.mesh.2d(boundary = inla.sp2segment(AL[1]),
cutoff = 0.1 * 50,
max.edge = c(0.2, 0.8) * 70,
offset = c(0.1, 0.2) * 150,
crs = st_crs(proj))The first dataset we consider is obtained from the North American Breeding Bird Survey. This dataset may be loaded directly from the package, and contains observations of the species between 2015 and 2017. This dataset is treated as replicate present-absent, where every point is assumed to be a visited site (or route).
The second dataset considered is obtained via the citizen science program, eBird. These data are obtained via the R package, spocc using the script below, where a separate object of data points was created for each year to ensure that the number of records per year is equal.
eBird2015 <- spocc::occ(
query = 'Colinus virginianus',
from = 'gbif',
date = c("2015-01-01", "2015-12-31"),
geometry = st_bbox(st_transform(AL,
'+proj=longlat +datum=WGS84 +no_defs'))
)$gbif
eBird2016 <- spocc::occ(
query = 'Colinus virginianus',
from = 'gbif',
date = c("2016-01-01", "2016-12-31"),
geometry = st_bbox(st_transform(AL,
'+proj=longlat +datum=WGS84 +no_defs'))
)$gbif
eBird2017 <- spocc::occ(
query = 'Colinus virginianus',
from = 'gbif',
date = c("2017-01-01", "2017-12-31"),
geometry = st_bbox(st_transform(AL,
'+proj=longlat +datum=WGS84 +no_defs'))
)$gbif
eBird <- data.frame(eBird2015$data[[1]]) %>%
bind_rows(data.frame(eBird2016$data[[1]])) %>%
bind_rows(data.frame(eBird2017$data[[1]]))
eBird <- st_as_sf(x = eBird,
coords = c('longitude', 'latitude'),
crs = '+proj=longlat +datum=WGS84 +no_defs')
eBird$Year <- eBird$year
eBird <- st_transform(eBird, proj)
eBird <- eBird[unlist(st_intersects(AL, eBird)),]We then get onto the model description, which in this case includes a
shared spatial field between the two datasets. This shared spatial field
is characterized by an ar1 process. To add this structure into
the model, we specify the parameter temporalModel in the
function startISDM appropriately, Furthermore we specified
the hyper parameters for both the random field and the temporal effect
and the priors for the intercepts.
hyperParams <- list(model = 'ar1',
hyper = list(rho = list(prior = "pc.prec", param = c(0.9, 0.1))))
modelSetup <- startISDM(eBird, BBSColinusVirginianus,
temporalName = 'Year',
Boundary = AL,
Projection = proj, Mesh = mesh,
responsePA = 'NPres', trialsPA = 'Ntrials')
modelSetup$specifySpatial(sharedSpatial = TRUE, prior.sigma = c(1, 0.5),
prior.range = c(100, 0.5))
modelSetup$specifyRandom(temporalModel = hyperParams)
modelSetup$priorsFixed(Effect = 'intercept', mean.linear = 0, prec.linear = 0.001)The data is spread across the map like this
The components for this model look like this:
Next we run the model, using the function fitISDM. Due
to time considerations, inference for this model is not completed in the
vignette. However, both the data and script is provided for the user to
complete the analysis.
And finally create predictions for the three time periods, and plot them.