Map: point density

map_dens() creates a smoothed utilisation distribution (UD).

Usage

as.im.SpatRaster(.map)

as.owin.SpatRaster(.map, .im = NULL)

as.owin.sf(.poly, .bbox = sf::st_bbox(.poly), .invert = TRUE)

bw.h(X)

map_dens(
  .map,
  .owin = as.owin.SpatRaster(.map),
  .coord = NULL,
  .discretise = FALSE,
  .shortcut = list(),
  .sigma = bw.h,
  ...,
  .fterra = FALSE,
  .plot = TRUE,
  .tryCatch = TRUE,
  .verbose = getOption("patter.verbose")
)

Arguments

.map

A terra::SpatRaster that defines the grid on which the UD is represented. If .coord = NULL, .map also defines the points (and associated weights) that are smoothed (see .map_coord()). The coordinate reference system of .map must be planar and specified.

.im, .owin

A pixel image representation of .map (see as.im.SpatRaster() and spatstat.geom::im()) and an observation window (see as.owin.SpatRaster(), as.owin.sf() and spatstat.geom::owin()). If un-supplied, .owin is defined automatically from .map via as.owin.SpatRaster(), which uses as.im.SpatRaster() internally (see Details). For faster results, use a rectangular or polygon observation window (see as.owin.sf()).

.poly, .bbox, .invert

For as.owin.sf() to construct observation windows from sf objects.

• .poly is an sf polygon object;

• .bbox is the bounding of a simple feature (see sf::st_bbox());

• .invert is a logical variable that defines whether or not to invert .poly (e.g., to turn a terrestrial polygon into an aquatic polygon);

.coord

(optional) Coordinates for density estimation, provided in any format accepted by .map_coord(). Coordinates must be planar.

.discretise

If .coord is provided, .discretise is a logical variable that defines whether or not to discretise coordinates on .map (see .map_coord()).

.shortcut

(optional) A named list from a previous call to map_dens(). If supplied, the function short-cuts straight to smoothing (.owin, .coord and .discretise are silently unused).

.sigma, X

.sigma is a numeric value or a function that specifies the smoothing bandwidth (passed to spatstat.explore::density.ppp()'s sigma argument). The default option is bw.h(), which sets the bandwidth based on combined variance of summarised coordinates, formatted as a point pattern X (see spatstat.geom::ppp), using the ad-hoc method (Worton, 1989). See spatstat functions (e.g., spatstat.explore::bw.diggle()) for more sophisticated methods.

...

Arguments for density estimation, passed to spatstat.explore::density.ppp(), such as sigma (i.e., the bandwidth). at and se are not permitted.

.fterra

A logical variable that defines whether or not to parallelise terra::resample().

.plot

A logical variable that defines whether or not to plot the output.

.tryCatch

A logical variable that controls error handling:

• If .tryCatch = FALSE, if density estimation fails with an error, the function fails with the same error.

• If .tryCatch = TRUE, if density estimation fails with an error, the function produces a warning with the error message and returns NULL.

.verbose

User output control (see patter-progress for supported options).

Value

The function returns a named list, with the following elements:

• x: a spatstat.geom::ppp object that defines points for density estimation;

• D: a spatstat.geom::im object of estimated intensities, from spatstat.explore::density.ppp();

• ud: a normalised terra::SpatRaster;

D and ud are NULL if spatstat.explore::density.ppp() fails and .tryCatch = TRUE.

Details

map_dens() smooths (a) a terra::SpatRaster or (b) a set of inputted coordinates.

.map_coord() (and .map_mark()) are used to define coordinates and weights:

• If .coords is NULL, .map cell coordinates are used for density estimation and cell values are used as weights.

• If coordinates are supplied, coordinates are optionally re-expressed on .map and then used for density estimation. This option is generally faster. Coordinate weights are defined by .map_mark().

Cell coordinates and associated weights are converted to a spatstat.geom::ppp() object, which is passed, alongside the observation window (.owin), to spatstat.explore::density.ppp() for the estimation. Weights must sum to one.

as.im.SpatRaster(), as.owin.SpatRaster() and as.owin.sf() are helper functions that convert a terra::SpatRaster to a pixel image and an observation window (see spatstat.geom::owin()). as.im.SpatRaster() is based on maptools::as.im.RasterLayer(). as.owin.SpatRaster() either defines a rectangular window, if there are no NAs on .map, or converts .map directly to an owin object. Gridded observation windows, especially if high resolution, considerably slow down density estimation and may exhaust vector memory. Use rectangular windows, or convert sf objects to polygon windows (via as.owin.sf()) if possible.

If .shortcut is supplied, the preceding steps can be skipped and the function short-cuts straight to smoothing. Use this option if the preceding steps are slow and you want to trial different smoothing options (such as sigma functions).

Coordinates and associated weights are smoothed via spatstat.explore::density.ppp() into an image. Pixel resolution and smoothing parameters such as bandwidth can be controlled via ... arguments which are passed directly to this function. The default bandwidth is set via bw.h() (see .sigma). The output is translated into a gridded probability density surface (on the geometry defined by .map). This process may use terra::resample(), which can be parallelised via .fterra (which controls the threads argument of that function).

This function replaces flapper::kud*() and flapper::pf_kud*() routines based on adehabitatHR (see here).

On Linux, these functions cannot be used within a Julia session.

References

Worton, B. J. (1989). Kernel Methods for Estimating the Utilization Distribution in Home-Range Studies. Ecology 70, 164–168. doi: 10.2307/1938423

See also

map_*() functions build maps of space use:

• map_pou() maps probability-of-use;

• map_dens() maps point density;

• map_hr_*() functions map home ranges;

All maps are represented as terra::SpatRasters.

To derive coordinates for mapping patterns of space use for tagged animals, see:

• coa() to calculate centre-of-activity;

• pf_filter() and associates to sample locations using particle filtering;

Examples

if (patter_run(.julia = FALSE, .geospatial = TRUE)) {

  library(data.table)
  library(dtplyr)
  library(dplyr, warn.conflicts = FALSE)
  library(spatstat.explore)

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

  #### Example (1): Use sample coordinates
  # Sample example coordinates
  coord <-
    map |>
    terra::spatSample(size = 100L, xy = TRUE, cell = TRUE, na.rm = TRUE) |>
    select("x", "y") |>
    as.data.table()
  # Use x, y coordinates
  map_dens(map, .coord = coord)
  # Other formats are acceptable
  map_dens(map, .coord = as.matrix(coord))
  map_dens(map, .coord = as.data.frame(coord))
  # `cell_x` and `cell_y` coordinates are acceptable
  map_dens(map, .coord = coord[, .(cell_x = x, cell_y = y)])
  # A SpatRaster is also acceptable
  pou <- map_pou(map, .coord = coord, .plot = FALSE)$ud
  map_dens(pou)

  #### Example (2): Use coordinates from `sim_path_walk()`
  coord <- dat_path()
  map_dens(map, .coord = coord)
  points(coord$x, coord$y, cex = 0.5)

  #### Example (3): Use coordinates from `coa()`
  # Use example dataset
  coord <- dat_coa()
  map_dens(map, .coord = coord)
  points(coord$x, coord$y, cex = 0.5)

  #### Example (4): Use a time series of coordinates from `pf_*()`
  # Use example dataset
  # * We use particles from the forward filter (`?pf_filter()`);
  # * Particles are equally weighted b/c re-sampling is implemented every time step;
  # * It is better to use outputs from the particle smoother;
  coord <- dat_pff()$states
  map_dens(map, .coord = coord)
  # points(coord$x, coord$y, cex = 0.5)

  #### Example (5): Control smoothing via `spatstat.explore::density.ppp()`
  # E.g. use fixed bandwidth:
  map_dens(map, .coord = coord, .sigma = 5)
  map_dens(map, .coord = coord, .sigma = 100)
  # E.g., perform automatic bandwidth selection using cross validation:
  if (FALSE) {
    # These examples are slow
    map_dens(map, .coord = coord, .sigma = bw.diggle) # 2 s
    map_dens(map, .coord = coord, .sigma = bw.scott)  # 1 s
    map_dens(map, .coord = coord, .sigma = bw.ppl)    # 65 s
    map_dens(map, .coord = coord, .sigma = bw.CvL)    # 25 s
  }
}
#> Loading required package: spatstat.data
#> Loading required package: spatstat.univar
#> spatstat.univar 3.1-2
#> Loading required package: spatstat.geom
#> spatstat.geom 3.3-6
#> 
#> Attaching package: ‘spatstat.geom’
#> The following object is masked from ‘package:data.table’:
#> 
#>     shift
#> Loading required package: spatstat.random
#> spatstat.random 3.3-3
#> Loading required package: nlme
#> 
#> Attaching package: ‘nlme’
#> The following object is masked from ‘package:dplyr’:
#> 
#>     collapse
#> spatstat.explore 3.4-2
#> `patter::map_dens()` called @ 2025-04-22 09:30:56... 
#> ... 09:30:56: Processing `.map`... 
#> ... 09:30:56: Building XYM... 
#> ... 09:30:56: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:56: Estimating density surface... 
#> ... 09:30:56: Scaling density surface... 

#> `patter::map_dens()` call ended @ 2025-04-22 09:30:56 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:56... 
#> ... 09:30:56: Processing `.map`... 
#> ... 09:30:56: Building XYM... 
#> ... 09:30:56: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:56: Estimating density surface... 
#> ... 09:30:56: Scaling density surface... 
#> `patter::map_dens()` call ended @ 2025-04-22 09:30:56 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:56... 
#> ... 09:30:56: Processing `.map`... 
#> ... 09:30:56: Building XYM... 
#> ... 09:30:56: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:56: Estimating density surface... 
#> ... 09:30:56: Scaling density surface... 
#> `patter::map_dens()` call ended @ 2025-04-22 09:30:57 (duration: ~1 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:57... 
#> ... 09:30:57: Processing `.map`... 
#> ... 09:30:57: Building XYM... 
#> ... 09:30:57: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:57: Estimating density surface... 
#> ... 09:30:57: Scaling density surface... 
#> `patter::map_dens()` call ended @ 2025-04-22 09:30:57 (duration: ~0 sec(s)). 
#> `patter::map_pou()` called @ 2025-04-22 09:30:57... 
#> ... Building XYM... 
#> ... Building SpatRaster... 
#> `patter::map_pou()` call ended @ 2025-04-22 09:30:57 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:57... 
#> ... 09:30:57: Processing `.map`... 
#> ... 09:30:57: Building XYM... 
#> ... 09:30:57: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:57: Estimating density surface... 
#> ... 09:30:57: Scaling density surface... 
#> `patter::map_dens()` call ended @ 2025-04-22 09:30:57 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:57... 
#> ... 09:30:57: Processing `.map`... 
#> ... 09:30:57: Building XYM... 
#> ... 09:30:57: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:57: Estimating density surface... 
#> ... 09:30:57: Scaling density surface... 

#> `patter::map_dens()` call ended @ 2025-04-22 09:30:57 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:57... 
#> ... 09:30:57: Processing `.map`... 
#> ... 09:30:57: Building XYM... 
#> ... 09:30:57: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:57: Estimating density surface... 
#> ... 09:30:57: Scaling density surface... 

#> `patter::map_dens()` call ended @ 2025-04-22 09:30:57 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:57... 
#> ... 09:30:57: Processing `.map`... 
#> ... 09:30:57: Building XYM... 
#> ... 09:30:57: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:57: Estimating density surface... 
#> ... 09:30:57: Scaling density surface... 

#> `patter::map_dens()` call ended @ 2025-04-22 09:30:57 (duration: ~0 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:57... 
#> ... 09:30:57: Processing `.map`... 
#> ... 09:30:57: Building XYM... 
#> ... 09:30:58: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:58: Estimating density surface... 
#> ... 09:30:58: Scaling density surface... 

#> `patter::map_dens()` call ended @ 2025-04-22 09:30:58 (duration: ~1 sec(s)). 
#> `patter::map_dens()` called @ 2025-04-22 09:30:58... 
#> ... 09:30:58: Processing `.map`... 
#> ... 09:30:58: Building XYM... 
#> ... 09:30:58: Defining `ppp` object... 
#> Observation window is gridded.
#> ... 09:30:58: Estimating density surface... 
#> ... 09:30:58: Scaling density surface... 

#> `patter::map_dens()` call ended @ 2025-04-22 09:30:58 (duration: ~0 sec(s)).