Introduction

Spatial capture–recapture estimates population density from a sample of \(n\) marked individuals detected at an array of detectors. Density is treated as the intensity of a point process where each point is the activity centre of an individual. Overdispersion of \(n\) relative to the expected variance (Poisson or binomial) leads to poor coverage of confidence intervals for density. This document demonstrates the use of secrRFS to simulate overdispersion of \(n\) from a known detection model and generating process for the AC.

The generating process is assumed to be a Cox process i.e. an inhomogeneous Poisson process (IHPP) in which the intensity surface varies stochastically between realisations (Mller and Waagepetersen 2004). The Thomas cluster process is treated as a Cox process with intensity surface defined by each realisation of parent points. The conditional case (fixed number in known area \(A\)) is not handled.

The available generating processes are

Generating process	Source function	secrRFS function	Parameters
Binary habitat mosaic	secr::randomHabitat	randomDensity	D, p, A, rescale
Log-gaussian random field	spatstat.random::rLGCP	randomGaussian	D, var, scale
Thomas cluster process	spatstat.random::rThomas	randomParents	D, mu, scale

secrRFS requires version 4.6.5 or later of package secr (Efford 2025). Package spatstat is required for the log-gaussian and Thomas process options.

Theory

We define the local density \(D_a\) as \[ D_a \equiv \frac{\int D(\mathbf{x}) \, p_\cdot (\mathbf{x} | \theta) \, d\mathbf{x}}{\int p_\cdot (\mathbf{x} | \theta) \, d\mathbf{x}} \] where \(D(\mathbf{x})\) is the intensity at \(\mathbf{x}\) of one realisation and \(p_\cdot(\mathbf{x}|\theta)\) is the probability that an AC at \(\mathbf{x}\) will be detected at least once, given vector of detection parameters \(\theta\).

Values of \(D_a\) are simulated by generating multiple realisations of \(D(\mathbf{x})\).

Overdispersion is estimated as \(c = 1 + a(\theta) \mbox{var}(D_a) / \mbox E(D_a).\)

Simple example

See here for an overview of ‘secr’.

library(secrRFS)
setNumThreads(2)   # adjust to number of available cores

detectpar <- list(lambda0 = 0.5, sigma = 1)
grid144 <- make.grid(12,12, detector = 'proximity', spacing = 2.0)
grid144mask <- make.mask(grid144, spacing = 0.5, buffer = 4)
D <- 256/maskarea(grid144mask)

parm1 <- list(D = D, A = 0.5, p = 0.5)
RFS(
    randomfn   = randomDensity,  
    parm       = parm1,
    nrepl      = 1000, 
    traps      = grid144, 
    mask       = grid144mask,
    detectpar  = detectpar,
    noccasions = 5)

## [1] 6.403554

The result is a value for \(c\) that may be used as a variance inflation factor.

Generating functions

Parameters are provided in the ‘parm’ argument whose interpretation is specific to each generating function. In each case ‘the ’parm’ includes the global density D.

`randomDensity`

This function generates an intensity surface as a random binary habitat mosaic, following Saura and Martinez-Millan (2000) and the randomHabitat function in secr.

Parameter ‘A’ (range 0–1) controls the habitat proportion, designated \(f\) by Efford and Fletcher (unpubl.).

Parameter ‘p’ controls fragmentation of the habitat patches. This is non-intuitive. Cluster size increases up to the ‘percolation threshold’ (about p \(\approx 0.59\) in the default case; Saura and Martinez-Millan 2000 p.664).

Parameter ‘rescale’ (default TRUE) determines whether within-patch density is scale by 1/A to maintain constant global expected density.

Di <- randomDensity(grid144mask,  parm = parm1, plt = TRUE, border = 1, 
    breaks = c(0,10,10000), col=c('grey','white'), legend = FALSE)
pop <- sim.popn(D = Di, grid144mask, model2D = 'IHP')
plot(pop, add = TRUE, frame = FALSE)

Fig. 1. Realisation of random habitat mosaic from randomDensity (white) with one point simulation.

`randomGaussian`

This function generates an intensity surface as a log-gaussian random field using rLGCP in the spatstat package.

Parameter ‘var’ is the variance of in situ log intensity. Mean log intensity (\(\mu\)) is computed as log(D) - var/2.

Parameter ‘scale’ is the spatial scale of the exponential covariance function.

parm2 <- list(D = D, var = 0.25, scale = 5)
Di <- randomGaussian(grid144mask,  parm = parm2, plt = TRUE, border = 1)

## Loading required namespace: spatstat

pop <- sim.popn(D = Di, grid144mask, model2D = 'IHP')
plot(pop, add = TRUE, frame = FALSE)

Fig. 2. Realisation of log-gaussian density surface from randomGaussian with one point simulation.

`randomParents`

This function generates an intensity surface from the distribution of parents for one realisation of a Thomas process.

Parameter ‘mu’ is the expected number of offspring per parent.

Parameter ‘scale’ is the spatial scale of the circular bivariate normal scatter of offspring around the parent.

The surface is derived from the spatstat function rThomas with argument ‘saveLambda = TRUE’.

parm3 <- list(D = D, mu = 8, scale = 2)
Di <- randomParents(grid144mask,  parm = parm3, plt = TRUE, parentcex = 0.6,
                    border = 5)
pop <- sim.popn(D = Di, grid144mask, model2D = 'IHP')
plot(pop, add = TRUE, frame = FALSE)

Fig. 3. Realisation of density surface from randomParents with one point simulation. Parent locations are shown as filled dots.

Parameter combinations

The function carray makes it easy to see how varying parameters of the Cox process affect overdispersion.

parmlevels <- list(D = D, mu = 2^(0:5), scale = c(1e-04, 1, 2, 4, 8))
RFSargs <- list (
    randomfn   = randomParents, 
    nrepl      = 1000, 
    traps      = grid144, 
    mask       = grid144mask, 
    detectfn   = 'HHN', 
    detectpar  = detectpar, 
    noccasions = 5)
ca <- carray (parmlevels, RFSargs)

## Completed in 11.36 minutes at 10:43:55 27 May 2026

round(ca,1)

##     scale
## mu   1e-04    1    2    4    8
##   1    1.9  1.9  1.8  1.7  1.4
##   2    2.8  2.7  2.7  2.3  1.9
##   4    4.7  4.7  4.3  3.6  2.7
##   8    8.2  8.3  7.5  6.3  4.5
##   16  15.5 15.3 14.3 11.5  7.9
##   32  32.1 30.1 25.8 22.7 13.6

Scaling to constant expected number

The parameter list for each generating function includes the component maskscale. When this is set to TRUE the intensity surface is scaled to provide a constant expected number within the area of mask. The default is maskscale = FALSE. For randomParents the implied effect is to scale the expected number of offspring per parent \(\mu\) to keep \(\mbox{E}[N(A)]\) constant.

Scaling can have a large effect on \(c\). Here is an example with randomDensity (not run).

parmlevels <- list(D = D, p = c(0.25,0.5), A = c(0.0625,0.125,0.25,0.5,0.75,1), 
                   maskscale = c(FALSE,TRUE))
RFSargs <- list (
    randomfn   = randomDensity, 
    nrepl      = 1000, 
    traps      = grid144, 
    mask       = grid144mask, 
    detectfn   = 'HHN', 
    detectpar  = detectpar, 
    noccasions = 5)
ca <- carray (parmlevels, RFSargs)
round(ca,1)

References

Baddeley, A. and Turner, R. (2005) spatstat: An R package for analyzing spatial point patterns. Journal of Statistical Software 12, 1–42. DOI: 10.18637/jss.v012.i06

Borchers, D. L. and Efford, M. G. (2008) Spatially explicit maximum likelihood methods for capture–recapture studies. Biometrics 64, 377–385.

Efford, M. G. (2023). secr: Spatially explicit capture–recapture models. R package version 4.6.5. https://CRAN.R-project.org/package=secr/

Efford, M. G. and Fletcher, D. (2024) The effect of spatial overdispersion on confidence intervals for population density estimated by spatially explicit capture–recapture. bioRxiv https://doi.org/10.1101/2024.03.12.584742.

Mller, J., and Waagepetersen, R. P. (2004) Statistical inference and simulation for spatial point processes. Chapman & Hall/CRC.

Saura, S. and Martinez-Millan, J. (2000) Landscape patterns simulation with a modified random clusters method. Landscape Ecology, 15, 661–678.

secrRFS - SECR random field simulator