Package 'secr'

Description

Functions to estimate the density and size of a spatially distributed animal population sampled with an array of passive detectors, such as traps, or by searching polygons or transects.

Details

Spatially explicit capture–recapture is a set of methods for studying marked animals distributed in space. Data comprise the locations of detectors (traps, searched areas, etc. described in an object of class ‘traps’), and the detection histories of individually marked animals. Individual histories are stored in an object of class ‘capthist’ that includes the relevant ‘traps’ object.

Models for population density (animals per hectare) and detection are defined in secr using symbolic formula notation. Density models may include spatial or temporal trend. Possible predictors for detection probability include both pre-defined variables (t, b, etc.) corresponding to ‘time’, ‘behaviour’ and other effects), and user-defined covariates of several kinds. Habitat is distinguished from nonhabitat with an object of class ‘mask’.

Models are fitted in secr by maximizing either the full likelihood or the likelihood conditional on the number of individuals observed ($n$). Conditional likelihood models are limited to homogeneous Poisson density, but allow continuous individual covariates for detection. A model fitted with secr.fit is an object of class secr. Generic methods (plot, print, summary, etc.) are provided for each object class.

A link at the bottom of each help page takes you to the help index. Several vignettes complement the help pages:

The datasets captdata and ovenbird include examples of fitted models. For models fitted to other datasets see secr-version4.pdf Appendix 2.

Add-on packages extend the capability of secr and are documented separately. secrlinear enables the estimation of linear density (e.g., animals per km) for populations in linear habitats such as stream networks (secrlinear-vignette.pdf). secrdesign enables the assessment of alternative study designs by Monte Carlo simulation; scenarios may differ in detector (trap) layout, sampling intensity, and other characteristics (secrdesign-vignette.pdf). ipsecr fits some awkward models (e.g., for single-catch traps) by simulation and inverse prediction (ipsecr-vignette.pdf). openCR fits open population models, both non-spatial and spatial (openCR-vignette.pdf).

The analyses in secr extend those available in the software Density (see www.otago.ac.nz/density/ for the most recent version of Density). Help is available on the ‘DENSITY | secr’ forum at www.phidot.org and the Google group secrgroup. Feedback on the software is also welcome, including suggestions for additional documentation or new features consistent with the overall design.

Author(s)

Murray Efford [email protected]

References

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

Borchers, D. L. and Fewster, R. M. (2016) Spatial capture–recapture models. Statistical Science 31, 219–232.

Efford, M. G. (2004) Density estimation in live-trapping studies. Oikos 106, 598–610.

Efford, M. G. (2011) Estimation of population density by spatially explicit capture–recapture with area searches. Ecology 92, 2202–2207.

Efford, M. G., Borchers D. L. and Byrom, A. E. (2009) Density estimation by spatially explicit capture-recapture: likelihood-based methods. In: D. L. Thomson, E. G. Cooch and M. J. Conroy (eds) Modeling Demographic Processes in Marked Populations. Springer, New York. Pp. 255–269.

Efford, M. G., Borchers D. L. and Mowat, G. (2013) Varying effort in capture–recapture studies. Methods in Ecology and Evolution 4, 629–636.

Efford, M. G., Dawson, D. K. and Borchers, D. L. (2009) Population density estimated from locations of individuals on a passive detector array. Ecology 90, 2676–2682.

Efford, M. G., Dawson, D. K. and Robbins C. S. (2004) DENSITY: software for analysing capture-recapture data from passive detector arrays. Animal Biodiversity and Conservation 27, 217–228.

Efford, M. G. and Fewster, R. M. (2013) Estimating population size by spatially explicit capture–recapture. Oikos 122, 918–928.

Efford, M. G. and Hunter, C. M. (2017) Spatial capture–mark–resight estimation of animal population density. Biometrics 74, 411–420.

Efford, M. G. and Mowat, G. (2014) Compensatory heterogeneity in capture–recapture data.Ecology 95, 1341–1348.

Royle, J. A., Chandler, R. B., Sollmann, R. and Gardner, B. (2014) Spatial capture–recapture. Academic Press.

Royle, J. A. and Gardner, B. (2011) Hierarchical spatial capture–recapture models for estimating density from trapping arrays. In: A.F. O'Connell, J.D. Nichols and K.U. Karanth (eds) Camera Traps in Animal Ecology: Methods and Analyses. Springer, Tokyo. Pp. 163–190.

See Also

read.capthist, secr.fit, traps, capthist, mask

Examples

## Not run: 

## generate some data & plot
detectors  <- make.grid (nx = 10, ny = 10, spacing = 20,
    detector = "multi")
plot(detectors, label = TRUE, border = 0, gridspace = 20)
detections <- sim.capthist (detectors, noccasions = 5,
    popn = list(D = 5, buffer = 100),
    detectpar = list(g0 = 0.2, sigma = 25))
session(detections) <- "Simulated data"
plot(detections, border = 20, tracks = TRUE, varycol = TRUE)

## generate habitat mask
mask <- make.mask (detectors, buffer = 100, nx = 48)

## fit model and display results
secr.model <- secr.fit (detections, model = g0~b, mask = mask)
secr.model


## End(Not run)

Description

Tools to construct spatial covariates for existing mask or traps objects from a spatial data source.

Usage

addCovariates(object, spatialdata, columns = NULL, strict = FALSE, replace = FALSE)

Arguments

Details

The goal is to obtain the value(s) of one or more spatial covariates for each point (i.e. row) in object. The procedure depends on the data source spatialdata, which may be either a spatial coverage (raster or polygon) or an object with covariate values at points (another mask or traps object). In the first case, an overlay operation is performed to find the pixel or polygon matching each point. In the second case, a search is conducted for the closest point in spatialdata.

If spatialdata is a character value then it is interpreted as the name of a polygon shape file (excluding ‘.shp’).

If spatialdata is a SpatialPolygonsDataFrame, SpatialGridDataFrame or 'sf' object from sf then it will be used in an overlay operation as described.

If package terra has been installed then spatialdata may also be a RasterLayer from package raster or SpatRaster from terra. If provided counts should be a single name that will be used for the values (otherwise 'raster' will be used).

If spatialdata is a mask or traps object then it is searched for the closest point to each point in object, and covariates are drawn from the corresponding rows in covariates(spatialdata). By default (strict = FALSE), values are returned even when the points lie outside any cell of the mask.

Value

An object of the same class as object with new or augmented covariates attribute. Column names and types are derived from the input.

Warning

Use of a SpatialGridDataFrame for spatialdata is untested.

See Also

make.mask, read.mask, read.traps

Examples

## In the Lake Station skink study (see ?skink), habitat covariates were
## measured only at trap sites. Here we extrapolate to a mask, taking
## values for each mask point from the nearest trap.

LSmask <- make.mask(LStraps, buffer = 30, type = "trapbuffer")
tempmask <- addCovariates(LSmask, LStraps)
## show first few lines
head(covariates(tempmask))

Description

Add sighting data on unmarked individuals and/or unidentified marked individuals to an existing capthist object.

Usage

addSightings(capthist, unmarked = NULL, nonID = NULL, uncertain = NULL, verify = TRUE,
    ...)

Arguments

Details

The capthist object for mark-resight analysis comprises distinct marking and sighting occasions, defined in the markocc attribute of traps(capthist). Add this attribute to traps(capthist) with markocc before using 'addSightings'. See also read.traps and read.capthist.

Mark-resight data may be binary (detector type ‘proximity’) or counts (detector types ‘count’, 'polygon' or 'transect'). The detector type is an attribute of traps(capthist). Values in unmarked and nonID should be whole numbers, and may be greater than 1 even for binary proximity detectors because multiple animals may be detected simultaneously at one place.

Arguments unmarked, nonID, uncertain provide data for attributes ‘Tu’, ‘Tm’, ‘Tn’ respectively. They may take several forms

• a single integer, the sum of all counts*

• a matrix of the count on each occasion at each detector (dimensions K x S, where K is the number of detectors and S is the total number of occasions). Columns corresponding to marking occasions should be all-zero.

• for multi-session data, a list with components as above

• a character value with the name of a text file containing the data; the file will be read with read.table. The ... argument allows some control over how the file is read. The data format comprises at least S+1 columns. The first is a session identifier used to split the file when the data span multiple sessions; it should be constant for a single-session capthist. The remaining S columns contain the counts for occasions 1:S, one row per detector. Further columns may be present; they are ignored at present.

* although this is convenient, the full matrix of counts provides more flexibility (e.g., when you wish to subset by occasion), and enables modelling of variation across detectors and occasions.

Value

A capthist object with the same structure as the input, but with new sighting-related attributes Tu (sightings of unmarked animals) and/or Tm (unidentified sightings of marked animals). Input values, including NULL, overwrite existing values.

Warning

** Mark-resight data formats and models are experimental and subject to change **

See Also

markocc, read.capthist, read.traps, sim.resight, Tm, Tu, Tn, secr-markresight.pdf

Examples

## construct capthist object MRCH from text files provided in 
## 'extdata' folder, assigning attribute 'markocc' and add unmarked
## and marked sightings from respective textfiles

datadir <- system.file("extdata", package = "secr")
captfile <- paste0(datadir, '/MRCHcapt.txt')
trapfile <- paste0(datadir, '/MRCHtrap.txt')
Tufile <- paste0(datadir, '/Tu.txt')
Tmfile <- paste0(datadir, '/Tm.txt')

MRCH <- read.capthist(captfile, trapfile, detector = c("multi", 
    rep("proximity",4)), markocc = c(1,0,0,0,0))
MRCH1 <- addSightings(MRCH, Tufile, Tmfile)

## alternatively (ignoring marked, not ID sightings)

MRCH <- read.capthist(captfile, trapfile, detector = c("multi", 
    rep("proximity",4)), markocc = c(1,0,0,0,0))
Tu <- read.table(Tufile)[,-1]  # drop session column
MRCH2 <- addSightings(MRCH, unmarked = Tu)
summary(MRCH2)

Description

Animal locations determined by radiotelemetry can be used to augment capture–recapture data. The procedure in secr is first to form a capthist object containing the telemetry data and then to combine this with true capture–recapture data (e.g. detections from hair-snag DNA) in another capthist object. secr.fit automatically detects the telemetry data in the new object.

Usage

addTelemetry (detectionCH, telemetryCH, type = c('concurrent','dependent','independent'), 
   collapsetelemetry = TRUE, verify = TRUE) 

xy2CH (CH, inflation = 1e-08)

telemetrytype (object) <- value

telemetrytype (object, ...)

Arguments

Details

It is assumed that a number of animals have been radiotagged, and their telemetry data (xy-coordinates) have been input to telemetryCH, perhaps using read.capthist with detector = "telemetryonly" and fmt = "XY", or with read.telemetry.

A new capthist object is built comprising all the detection histories in detectionCH, plus empty (all-zero) histories for every telemetered animal not in detectionCH. Telemetry is associated with new sampling occasions and a new detector (nominally at the same point as the first in detectionCH). The number of telemetry fixes of each animal is recorded in the relevant cell of the new capthist object (CH[i, s, K+1] for animal i and occasion s if there were K detectors in detectionCH).

The new sampling occasion(s) are assigned the detector type ‘telemetry’ in the traps attribute of the output capthist object, and the traps attribute telemetrytype is set to the value provided. The telemetry type may be “independent” (no matching of individuals in captured and telemetered samples), “dependent” (telemetered animals are a subset of captured animals) or “concurrent” (histories may be capture-only, telemetry-only or both capture and telemetry).

The telemetry locations are carried over from telemetryCH as attribute ‘xylist’ (each component of xylist holds the coordinates of one animal; use telemetryxy to extract).

The default behaviour of 'addTelemetry' is to automatically collapse all telemetry occasions into one. This is computationally more efficient than the alternative, but closes off some possible models.

xy2CH partly reverses addTelemetry: the location information in the telemetryxy attribute is converted back to a capthist with detector type ‘telemetry’.

Value

A single-session capthist object with the same detector type as detectionCH, but possibly with empty rows and an ‘telemetryxy’ attribute.

Note

Telemetry provides independent data on the location and presence of a sample of animals. These animals may be missed in the main sampling that gives rise to detectionCH i.e., they may have all-zero detection histories.

The ‘telemetry’ detector type is used for telemetry occasions in a combined dataset.

See Also

capthist, make.telemetry, read.telemetry, telemetryxy telemetered

Examples

## Not run: 

# Generate some detection and telemetry data, combine them using
# addTelemetry, and perform analyses

# detectors
te <- make.telemetry()
tr <- make.grid(detector = "proximity")

# simulated population and 50% telemetry sample
totalpop <- sim.popn(tr, D = 20, buffer = 100)
tepop <- subset(totalpop, runif(nrow(totalpop)) < 0.5)

# simulated detection histories and telemetry
# the original animalID (renumber = FALSE) are needed for matching
trCH <- sim.capthist(tr,  popn = totalpop, renumber = FALSE, detectfn = "HHN")
teCH <- sim.capthist(te, popn = tepop, renumber=FALSE, detectfn = "HHN",
    detectpar = list(lambda0 = 3, sigma = 25))

combinedCH <- addTelemetry(trCH, teCH)

# summarise and display
summary(combinedCH)
plot(combinedCH, border = 150)
ncapt <- apply(combinedCH,1,sum)
points(totalpop[row.names(combinedCH)[ncapt==0],], pch = 1)
points(totalpop[row.names(combinedCH)[ncapt>0],], pch = 16)

# for later comparison of precision we must fix the habitat mask
mask <- make.mask(tr, buffer = 100)
fit.tr <- secr.fit(trCH, mask = mask, CL = TRUE, detectfn = "HHN")  ## trapping alone
fit.te <- secr.fit(teCH, mask = mask, CL = TRUE, start = log(20),   ## telemetry alone
    detectfn = "HHN") 
fit2   <- secr.fit(combinedCH, mask = mask, CL = TRUE,              ## combined
    detectfn = "HHN")                                 

# improved precision when focus on realised population
# (compare CVD)
derived(fit.tr, distribution = "binomial")
derived(fit2, distribution = "binomial")


# may also use CL = FALSE
secr.fit(combinedCH, CL = FALSE, detectfn = "HHN", trace = FALSE)

## End(Not run)

Description

Terse report on the fit of one or more spatially explicit capture–recapture models. Models with smaller values of AIC (Akaike's Information Criterion) are preferred. Extraction ([) and logLik methods are included.

Usage

## S3 method for class 'secr'
AIC(object, ..., sort = TRUE, k = 2, dmax = 10, criterion = c("AIC","AICc"), chat = NULL)
## S3 method for class 'secrlist'
AIC(object, ..., sort = TRUE, k = 2, dmax = 10, criterion = c("AIC","AICc"), chat = NULL)
## S3 method for class 'secr'
logLik(object, ...)
secrlist(..., names = NULL)
## S3 method for class 'secrlist'
x[i]

Arguments

Details

Models to be compared must have been fitted to the same data and use the same likelihood method (full vs conditional). From version 4.1 a warning is issued if AICcompatible reveals a problem.

AIC is given by

$\mbox{AIC} = -2\log(L(\hat{\theta})) + 2K$

where $K$ is the number of "beta" parameters estimated.

AIC with small sample adjustment is given by

$\mbox{AIC}_c = -2\log(L(\hat{\theta})) + 2K + \frac{2K(K+1)}{n-K-1}.$

The sample size $n$ is the number of individuals observed at least once (i.e. the number of rows in capthist).

Model weights are calculated as

$w_i = \frac{\exp(-\Delta_i / 2),}{ \sum{\exp(-\Delta_i / 2)}}$

where $\Delta$ refers to differences in AIC or AICc depending on the argument ‘criterion’ (see Notes).

Models for which delta > dmax are given a weight of zero and are excluded from the summation. Model weights may be used to form model-averaged estimates of real or beta parameters with modelAverage (see also Buckland et al. 1997, Burnham and Anderson 2002).

The argument k is included for consistency with the generic method AIC.

secrlist forms a list of fitted models (an object of class ‘secrlist’) from the fitted models in .... Arguments may include secrlists. If secr components are named the model names will be retained unless ‘names’ is specified. (see Examples).

If chat ($\hat c$) is provided then quasi-AIC values are computed (secr >= 4.6.0):

$\mbox{QAIC} = -2\log(L(\hat{\theta}))/ \hat c + 2K.$

Value

A data frame with one row per model. By default, rows are sorted by ascending 'criterion' (default AIC see Notes).

And depending on criterion:

or

logLik.secr returns an object of class ‘logLik’ that has attribute df (degrees of freedom = number of estimated parameters).

If the variance inflation factor 'chat' is provided then outputs AIC, AICc etc. are replaced by the corresponding quasi-AIC values labelled QAIC, QAICc etc.

Note

It is not be meaningful to compare models by AIC if they relate to different data (see AICcompatible).

Specifically:

• an ‘secrlist’ generated and saved to file by mask.check may be supplied as the object argument of AIC.secrlist, but the results are not informative

• models fitted by the conditional likelihood (CL = TRUE) and full likelihood (CL = FALSE) methods cannot be compared

• hybrid mixture models (using hcov argument of secr.fit) should not be compared with other models

• grouped models (using groups argument of secr.fit) should not be compared with other models

• multi-session models should not be compared with single-session models based on the same data.

A likelihood-ratio test (LR.test) is a more direct way to compare two models.

The issue of goodness-of-fit and possible adjustment of AIC for overdispersion has yet to be addressed (cf QAIC in MARK).

The user may select between AIC and AICc for comparing models. AICc is widely used, but AIC may be better for model averaging even when samples are small (Turek and Fletcher 2012; Fletcher 2019, p. 60). The default was changed to AIC in version 5.0.0).

References

Buckland S. T., Burnham K. P. and Augustin, N. H. (1997) Model selection: an integral part of inference. Biometrics 53, 603–618.

Burnham, K. P. and Anderson, D. R. (2002) Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach. Second edition. New York: Springer-Verlag.

Fletcher, D. (2019) Model averaging. SpringerBriefs in Statistics. Berlin: Springer-Verlag.

Hurvich, C. M. and Tsai, C. L. (1989) Regression and time series model selection in small samples. Biometrika 76, 297–307.

Turek, D. and Fletcher, D. (2012) Model-averaged Wald confidence intervals. Computational statistics and data analysis 56, 2809–2815.

See Also

AICcompatible, modelAverage, AIC, secr.fit, print.secr, score.test, LR.test, deviance.secr

Examples

## Compare two models fitted previously
## secrdemo.0 is a null model
## secrdemo.b has a learned trap response

AIC(secrdemo.0, secrdemo.b)

## Form secrlist and pass to AIC.secr
temp <- secrlist(null = secrdemo.0, learnedresponse = secrdemo.b)
AIC(temp)

Description

Determine whether models can be compared by AIC. Incompatibility may be due to difference in the data or the specifications of the groups, hcov or binomN arguments to secr.fit,

Usage

## S3 method for class 'secr'
AICcompatible(object, ...)
## S3 method for class 'secrlist'
AICcompatible(object, ...)

Arguments

Details

The capthist objects are checked for strict identity with the function identical.

All elements in the output must be TRUE for valid AIC comparison or model averaging using AIC or AICc.

Value

Named logical vector with elements ‘data’, ‘CL’, ‘groups’, ‘hcov’ and ‘binomN’.

See Also

AIC.secr, modelAverage

Examples

AICcompatible(secrdemo.0, secrdemo.CL)

## Not run: 

## A common application of AICcompatible() is to determine 
## the compatibility of models fitted with and without the 
## fastproximity option.

ovenCHp1 <- reduce(ovenCHp, by = 'all', outputdetector = 'count')
ob1 <- secr.fit(ovenCHp, buffer = 300, details = list(fastproximity = TRUE))
ob2 <- secr.fit(ovenCHp1, buffer = 300, details = list(fastproximity = FALSE))
ob3 <- secr.fit(ovenCHp1, buffer = 300, details = list(fastproximity = FALSE), binomN = 1)
AICcompatible(ob1,ob2)
AICcompatible(ob1,ob3)


## End(Not run)

Description

Method for generic as.data.frame function that partially reverses make.capthist.

Usage

## S3 method for class 'capthist'
as.data.frame(x, row.names = NULL, optional = FALSE, covariates = FALSE, 
                                    fmt = c("trapID", "XY"), ...)
## S3 method for class 'traps'
as.data.frame(x, row.names = NULL, optional = FALSE, usage = FALSE, 
                              covariates = FALSE, ...)
## S3 method for class 'capthist'
as.array(x, ...)

Arguments

Details

By default individual covariates are not exported. When exported they are repeated for each detection of an individual.

Value

A data frame or list of data frames (in the case of a multisession input).

For capthist objects –

The core columns are (Session, ID, Occasion, TrapID) or (Session, ID, Occasion, x, y), depending on the value of fmt. Additional columns for covariates and signal strength (detector ‘signal’) are appended to the right.

For traps objects –

The core columns are (x, y). Usage columns are named u1, u2, ..., uS where S is the number of occasions.

The as.array method for capthist objects returns an object with the same dimensions and dimnames but different class, or a list of such objects in the case of multisession input.

Examples

as.data.frame (captdata)
  as.data.frame (traps(captdata))

Description

This function is used primarily for plotting covariates, for which the plot.mask function has greater functionality than plot.traps. It also generates pretty maps of grid cells.

Usage

as.mask(x)

Arguments

Details

A mask derived by coercion with as.mask may behave unpredictably e.g., in secr.fit.

Value

If x is a single-session traps object –

an object of class c("mask", "data.frame")

If x is a multi-session traps object –

an object of class c("mask", "list"), for which each component is a single-session mask.

See Also

make.mask, plot.mask, mask, traps

Examples

plot(as.mask(traps(captdata)), dots = FALSE, meshcol = "black")
plot(traps(captdata), add = TRUE)

Description

This function converts a spatstat "ppp" object (Baddeley et al. 2015), making it easier to use the simulation capability of spatstat in secr.

Usage

as.popn(x)

Arguments

Details

Not all attributes are carried over.

Value

An object of class c("popn", "data.frame") with attribute "boundingbox". The attribute "Lambda" (spatstat class "im") is also carried over if present (used for the intensity surface of LGCP simulations).

References

Baddeley, A., Rubak, E., and Turner, R. 2015. Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press, London. ISBN 9781482210200, https://www.routledge.com/Spatial-Point-Patterns-Methodology-and-Applications-with-R/Baddeley-Rubak-Turner/p/book/9781482210200/.

See Also

sim.popn, popn

Description

Find plausible initial parameter values for secr.fit. A simple SECR model is fitted by a fast ad hoc method.

Usage

autoini(capthist, mask, detectfn = 0, thin = 0.2, tol = 0.001, 
    binomN = 1, adjustg0 = TRUE, adjustsigma = 1.2, ignoreusage = FALSE, 
    ncores = NULL)

Arguments

Details

Plausible starting values are needed to avoid numerical problems when fitting SECR models. Actual models to be fitted will usually have more than the three basic parameters output by autoini; other initial values can usually be set to zero for secr.fit. If the algorithm encounters problems obtaining a value for g0, the default value of 0.1 is returned.

Only the halfnormal detection function is currently available in autoini (cf other options in e.g. detectfn and sim.capthist).

autoini implements a modified version of the algorithm proposed by Efford et al. (2004). In outline, the algorithm is

1. Find value of sigma that predicts the 2-D dispersion of individual locations (see RPSV).

2. Find value of g0 that, with sigma, predicts the observed mean number of captures per individual (by algorithm of Efford et al. (2009, Appendix 2))

3. Compute the effective sampling area from g0, sigma, using thinned mask (see esa)

4. Compute D = $n$/esa(g0, sigma), where $n$ is the number of individuals detected

Here ‘find’ means solve numerically for zero difference between the observed and predicted values, using uniroot.

Halfnormal sigma is estimated with RPSV(capthist, CC = TRUE). The factor adjustsigma is applied as a crude correction for truncation of movements at the edge of the detector array.

If RPSV cannot be computed the algorithm tries to use observed mean recapture distance $\bar{d}$. Computation of $\bar{d}$ fails if there no recaptures, and all returned values are NA.

If the mask has more than 100 points then a proportion 1–thin of points are discarded at random to speed execution.

The argument tol is passed to uniroot. It may be a vector of two values, the first for g0 and the second for sigma.

If traps(capthist) has a usage attribute (defining effort on each occasion at each detector) then the value of g0 is divided by the mean of the non-zero elements of usage. This adjustment is not precise.

If adjustg0 is TRUE then an adjustment is made to g0 depending on the value of binomN. For Poisson counts (binomN = 0) the adjustment is linear on effort (adjusted.g0 = g0 / usage). Otherwise, the adjustment is on the hazard scale (adjusted.g0 = 1 – (1 – g0) ^ (1 / (usage x binomN))). An arithmetic average is taken over all non-zero usage values (i.e. over used detectors and times). If usage is not specified it is taken to be 1.0.

Setting ncores = NULL uses the existing value from the environment variable RCPP_PARALLEL_NUM_THREADS (see setNumThreads).

Value

A list of parameter values :

Note

autoini always uses the Euclidean distance between detectors and mask points.

You may get this message from secr.fit: “'autoini' failed to find g0; setting initial g0 = 0.1”. If the fitted model looks OK (reasonable estimates, non-missing SE) there is no reason to worry about the starting values. If you get this message and model fitting fails then supply your own values in the start argument of secr.fit.

References

Efford, M. G., Dawson, D. K. and Robbins C. S. (2004) DENSITY: software for analysing capture–recapture data from passive detector arrays. Animal Biodiversity and Conservation 27, 217–228.

Efford, M. G., Dawson, D. K. and Borchers, D. L. (2009) Population density estimated from locations of individuals on a passive detector array. Ecology 90, 2676–2682.

See Also

capthist, mask, secr.fit, dbar

Examples

## Not run: 

demotraps <- make.grid()
demomask <- make.mask(demotraps)
demoCH <- sim.capthist (demotraps, popn = list(D = 5, buffer = 100), seed = 321)
autoini (demoCH, demomask)


## End(Not run)

Description

Forms a new covariate, replacing values of an old covariate by the central value of equal-width bins.

Usage

binCovariate(object, covname, width)

Arguments

Details

The name of the new covariate is paste0(covname, width).

Fails if covariate not found or is not numeric or there is already a covariate with the new name.

Multi-session objects are handled appropriately.

Value

Object of the same class as the input with new covariate.

See Also

covariates, skink

Examples

# bin values of skink snout-vent length (mm)
infraCH <- binCovariate (infraCH, 'SVL', 5)
table(covariates(infraCH[[1]])$SVL5)

# bin values of trap covariate 'HtBrack' (height of bracken, cm)
traps(infraCH) <- binCovariate(traps(infraCH), "HtBrack", 20)
table(covariates(traps(infraCH)[[1]])$HtBrack20)

Description

American black bears Ursus americanus were surveyed with baited hair snags in the Great Smoky Mountains National Park, Tennessee, in the summer of 2003.

Usage

blackbearCH
GSM
blackbear.0
blackbear.h2bk

Details

American black bears Ursus americanus were surveyed in the Tennessee sector of Great Smoky Mountains National Park over 9 June–15 August 2003. Baited hair snags (barbed wire enclosures) were operated for 10 weeks at 65 sites, about 1 km apart and mostly close to trails. Bait consisted of bakery products in a small waxed-paper bag. Raspberry extract was used as a scent lure.

Genotyping and non-spatial capture-recapture analysis of a data subset were described by Settlage et al. (2008). The sex of each genotyped bear was determined subsequently and some additional samples were included (J. Laufenberg pers. comm. 2012-05-09).

The dataset is a single-session capthist object with binary proximity detector type. Snags were visited weekly, so there were 10 occasions in the raw data. A single covariate ‘sex’ was recorded for each individual.

The dataset comprises 282 detections of 81 females and 58 males. Female 15 apparently made a long movement (17 km) between occasions 1 and 3.

The hair snag array sampled less than 20% of the area of the park. The unforested area outside the park on the northwestern boundary of the study area was not considered to be black bear habitat (F. van Manen pers. comm. 2012-05-18) and should be excluded in analyses. The approximate boundary of the park is included as a shapefile ‘GSMboundary.shp’ in the ‘extdata’ folder of the package and as the sf sfc_POLYGON object GSM. The latter may be used in make.mask (see Examples).

Two models (blackbear.0 and blackbear.h2bk) were fitted as shown in the Examples.

Source

The data were provided by Jared Laufenberg, Frank van Manen and Joe Clark.

References

Settlage, K. E., Van Manen, F. T., Clark, J. D., and King, T. L. (2008) Challenges of DNA-based mark–recapture studies of American black bears. Journal of Wildlife Management 72, 1035–1042.

Examples

summary(blackbearCH)



# GSM is the approximate boundary of Great Smoky Mountains National Park
# Make a habitat mask restricted to the park

tr <- traps(blackbearCH)
msk <- make.mask(tr, buffer = 6000, type = 'trapbuffer', poly = GSM)

# Plot

plot(GSM)
plot(msk, add = TRUE)
plot(blackbearCH, tracks = TRUE, add = TRUE)
plot(tr, add = TRUE)

# Fit models
# suppress fastproximity to allow learned response

setNumThreads()   # as appropriate

# null model
blackbear.0 <- secr.fit(blackbearCH,  detectfn = 'EX', hcov = 'sex', 
    mask = msk, details = list(fastproximity = FALSE), trace = FALSE)

# sex differences and site-specific behavioural response
blackbear.h2bk <- secr.fit(blackbearCH,  detectfn = 'EX', hcov = 'sex', 
    model = list(g0~bk+h2, sigma~h2), mask = msk, 
    details = list(fastproximity = FALSE), trace = FALSE)
    
AIC(blackbear.0, blackbear.h2bk)
summary(blackbear.h2bk)   

# How many if we extrapolate to GSM NP?
region.N(blackbear.h2bk, region = GSM)

Description

Convert data between ‘capthist’ and BUGS input format.

Usage

read.DA(DAlist, detector = "polygonX", units = 1, session = 1,
    Y = "Y", xcoord = "U1", ycoord = "U2", xmin = "Xl",
    xmax = "Xu", ymin = "Yl", ymax = "Yu", buffer = "delta",
    verify = TRUE)

write.DA(capthist, buffer, nzeros = 200, units = 1)

Arguments

Details

Data for OpenBUGS or WinBUGS called from R using the package R2WinBUGS (Sturtz et al. 2005) take the form of an R list.

These functions are limited at present to binary data from a square quadrat such as used by Royle and Young (2008). Marques et al. (2011) provide an R function create.data() for generating simulated datasets of this sort (see sim.capthist for equivalent functionality).

When reading BUGS data –

The character values Y, xcoord, ycoord, xmin etc. are used to locate the data within DAlist, allowing for variation in the input names.

The number of sampling occasions is taken from the number of columns in Y. Each value in Y should be 0 or 1. Coordinates may be missing

A numeric value for buffer is the distance (in the original units) by which the limits Xl, Xu etc. should be shrunk to give the actual plot limits. If buffer is character then a component of DAlist contains the required numeric value.

Coordinates in the output will be multiplied by the scalar units.

Augmentation rows corresponding to ‘all-zero’ detection histories in Y, xcoord, and ycoord are discarded.

When writing BUGS data –

Null (all-zero) detection histories are added to the matrix of detection histories Y, and missing (NA) rows are added to the coordinate matrices xcoord and ycoord.

Coordinates in the output will be divided by the scalar units.

Value

For read.DA, an object of class ‘capthist’.

For write.DA, a list with the components

U1 and U2 are ‘NA’ where animal was not detected.

References

Marques, T. A., Thomas, L. and Royle, J. A. (2011) A hierarchical model for spatial capture–recapture data: Comment. Ecology 92, 526–528.

Royle, J. A. and Young, K. V. (2008) A hierarchical model for spatial capture–recapture data. Ecology 89, 2281–2289.

Sturtz, S., Ligges, U. and Gelman, A. (2005) R2WinBUGS: a package for running WinBUGS from R. Journal of Statistical Software 12, 1–16.

See Also

hornedlizardCH, verify, capthist

Examples

write.DA (hornedlizardCH, buffer = 100, units = 100)

## In this example, the input uses Xl, Xu etc.
## for the limits of the plot itself, so buffer = 0.
## Input is in hundreds of metres.
## First, obtain the list lzdata
olddir <- setwd (system.file("extdata", package="secr"))
source ("lizarddata.R")
setwd(olddir)
str(lzdata)
## Now convert to capthist
tempcapt <- read.DA(lzdata, Y = "H", xcoord = "X",
    ycoord = "Y", buffer = 0, units = 100)

## Not run: 

plot(tempcapt)
secr.fit(tempcapt, trace = FALSE)
## etc.


## End(Not run)

Description

A capthist object encapsulates all data needed by secr.fit, except for the optional habitat mask.

Details

An object of class capthist holds spatial capture histories, detector (trap) locations, individual covariates and other data needed for a spatially explicit capture-recapture analysis with secr.fit.

A capthist is primarily an array of values with dim(capthist) = c(nc, noccasions, ntraps) where nc is the number of detected individuals. Values maybe binary ({–1, 0, 1}) or integer depending on the detector type.

Deaths during the experiment are represented as negative values.

Ancillary data are retained as attributes of a capthist object as follows:

• traps — object of class traps (required)

• session — session identifier (required)

• covariates — dataframe of individual covariates (optional)

• cutval — threshold of signal strength for detection (‘signal’ only)

• signalframe — signal strength values etc., one row per detection (‘signal’ only)

• detectedXY — dataframe of coordinates for location within polygon (‘polygon’-like detectors only)

• xylist — coordinates of telemetered animals

• Tu — detectors x occasions matrix of sightings of unmarked animals

• Tm — detectors x occasions matrix of sightings of marked but unidentified animals

• Tn — detectors x occasions matrix of sightings with unknown mark status

read.capthist is adequate for most data input. Alternatively, the parts of a capthist object can be assembled with the function make.capthist. Use sim.capthist for Monte Carlo simulation (simple models only). Methods are provided to display and manipulate capthist objects (print, summary, plot, rbind, subset, reduce) and to extract and replace attributes (covariates, traps, xy).

A multi-session capthist object is a list in which each component is a capthist for a single session. The list maybe derived directly from multi-session input in Density format, or by combining existing capthist objects with MS.capthist.

Note

Early versions of secr (before 3.0) used an individual x occasion matrix for data from single-catch and multi-catch traps, instead of a 3-D array. Entries in the matrix corresponded to trap numbers. The function updateCH converts the old format.

References

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

Efford, M. G., Borchers D. L. and Byrom, A. E. (2009) Density estimation by spatially explicit capture-recapture: likelihood-based methods. In: D. L. Thomson, E. G. Cooch and M. J. Conroy (eds) Modeling Demographic Processes in Marked Populations. Springer, New York. Pp. 255–269.

See Also

traps, secr.fit, read.capthist, make.capthist, sim.capthist, subset.capthist, rbind.capthist, MS.capthist, reduce.capthist, mask

Description

Extract parts of an object of class ‘capthist’.

Usage

animalID(object, names = TRUE, sortorder = c("snk", "ksn"))
occasion(object, sortorder = c("snk", "ksn"))
trap(object, names = TRUE, sortorder = c("snk", "ksn"))
alive(object, sortorder = c("snk", "ksn"))
alongtransect(object, tol = 0.01)
xy(object)
xy(object) <- value
telemetryxy(object, includeNULL = FALSE)
telemetryxy(object) <- value
telemetered(object)

Arguments

Details

These functions extract data on detections, ignoring occasions when an animal was not detected. By default, detections are ordered by occasion, animalID and trap (sortorder = "snk"). The alternative is to order by trap, occasion and animalID (sortorder = "ksn"). (‘n’, ‘s’ and ‘k’ are the indices used internally for animals, occasions and traps respectively).

For historical reasons, "ksn" is used for locations within polygons and similar (xy).

trap returns polygon or transect numbers if traps(object) has detector type ‘polygon’ or ‘transect’.

alongtransect returns the distance of each detection from the start of the transect with which it is associated.

Replacement values must precisely match object in number of detections and in their order. xy<- expects a dataframe of x and y coordinates for points of detection within a ‘polygon’ or ‘transect’ detector. telemetryxy<- expects a list of dataframes, one per telemetered animal.

Value

For animalID and trap a vector of numeric or character values, one per detection.

For alive a vector of logical values, one per detection.

For occasion, a vector of numeric values, one per detection.

For xy, a dataframe with one row per detection and columns ‘x’ and ‘y’.

If object has multiple sessions, the result is a list with one component per session.

See Also

capthist, polyID, signalmatrix

Examples

## `captdata' is a demonstration dataset
animalID(captdata)

temp <- sim.capthist(popn = list(D = 1), make.grid(detector
    = "count"))
cbind(ID = as.numeric(animalID(temp)), occ = occasion(temp),
    trap = trap(temp))

Description

Activity centres may be clumped (overdispersed) relative to a Poisson distribution, the model used in secr.fit (Borchers and Efford 2008). This can cause the sampling variance of density estimates to be understated. One solution currently under investigation is to apply a variance inflation factor, a measure of overdispersion, based on the number of individuals detected at each detector (Bischof et al. 2020).

Functions described here compute the observed (nk) or expected (Enk) number of individuals detected at each detector and use that to compute Fletcher's $\hat c$ estimate of overdispersion for use as a variance inflation factor.

Enk uses exact formulae for 'multi', 'proximity' and 'count' detector types. Other types may be simulated by setting a positive value for 'nrepl', which should be large (e.g., nrepl = 10000).

adjustVarD adjusts the SE and confidence limits of density estimates using Fletcher's $\hat c$. The implementation is limited to simple detection models (see Warnings).

See Cooch and White (2022) for an introduction to measurement of overdispersion in capture–recapture. The focus here is on overdispersion of activity centres relative to a Poisson distribution, rather than on non-independence in the spatial detection process.

Usage

nk(capthist)

Enk(D, mask, traps, detectfn = 0, detectpar = list(g0 = 0.2,
    sigma = 25, z = 1), noccasions = NULL, binomN = NULL,
    userdist = NULL, ncores = NULL, nrepl = NULL) 

chat.nk(object, nsim = NULL, ...)

adjustVarD(object, chatmin = 1, alpha = 0.05, chat = NULL)

Arguments

Details

If traps has a usage attribute then noccasions is set accordingly; otherwise it must be provided.

The environment variable RCPP_PARALLEL_NUM_THREADS determines the number of parallel threads. It is set to the value of ncores, unless that is NULL (see setNumThreads).

A conventional variance inflation factor due to Wedderburn (1974) is $\hat c_X = X^2/(K-p)$ where $K$ is the number of detectors, $p$ is the number of estimated parameters, and

$X^2 = \sum_k (n_k - E (n_k))^2/ E(n_k).$

Fletcher's $\hat c$ is an improvement on $\hat c_X$ that is less affected by small expected counts. It is defined by

$\hat c = c_X / (1+ \bar s),$

where $\bar s = \sum_k s_k / K$ and $s_k = (n_k - E(n_k)) / E(n_k)$.

chat.nk may be used to simulate $\hat c$ values under the given model (set nsim > 0). The ... argument may include 'ncores = x' (x>1) to specify parallel processing of simulations - the speed up is large on unix-like machines for which the cluster type of makeCluster is "FORK" rather than "PSOCK". If 'ncores' is not provided then the value returned by setNumThreads() is used.

No adjustment is made by adjustVarD when $\hat c$ is less than the minimum. adjustVarD by default computes Fletcher's ‘chat’ using chat.nk, but a value may be provided.

If chat has been computed separately and provided in the argument of that name, adjustVarD also accepts a single dataframe as the argument ‘object’; the dataframe should have row ‘D’ and columns ‘link’, ‘estimate’, ‘SE.estimate’ as in the output from predict.secr.

Value

For nk, a vector of observed counts, one for each detector in traps(capthist).

For Enk, a vector of expected counts, one for each detector in traps.

For chat.nk, usually a list comprising –

There are two variations –

If ‘verbose = FALSE’ then only the numeric value of $\hat c$ is returned (a vector of 2 values if ‘type = "both"’).

If chat.nk is called with ‘nsim > 0’ then the output is a list comprising –

For adjustVarD, a dataframe with one row for each session, based on predict.secr or derived.secr, with extra column ‘c-hat’.

Warning

These functions are experimental in secr 4.6, and do not work with polygon-like and single-catch detectors. No allowance is made for modelled variation in detection parameters with respect to occasion, detector or animal; this includes mixture models (e.g., g0~h2).

Versions before 4.5.11 did not correctly compute expected counts for multi-catch detectors.

Furthermore, we doubt that the adjustment actually solves the problem of overdispersion (Efford and Fletcher unpubl.).

References

Bischof, R., P. Dupont, C. Milleret, J. Chipperfield, and J. A. Royle. 2020. Consequences of ignoring group association in spatial capture–recapture analysis. Wildlife Biology wlb.00649. DOI 10.2981/wlb.00649

Cooch, E. and White, G. (eds) (2022) Program MARK: A Gentle Introduction. 22nd edition. Most recent edition available online at www.phidot.org/software/mark/docs/book/.

Fletcher, D. (2012) Estimating overdispersion when fitting a generalized linear model to sparse data. Biometrika 99, 230–237.

Wedderburn, R. W. M. (1974) Quasi-likelihood functions, generalized linear models, and the Gauss-Newton method. Biometrika 61, 439–47.

See Also

secr, make.mask, Detection functions, Fletcher.chat

Examples

temptrap <- make.grid()
  msk <- make.mask(temptrap)
  ## expected number of individuals per detector (multi-catch) 
  Enk (D = 5, msk, temptrap, detectpar = list(g0 = 0.2, sigma = 25),
    noccasions = 5)



# useful plotting function for simulated chat (nsim>0)
plotchat <- function(chat, head = '', breaks = seq(0.5,2.5,0.05)) {
    hist(chat$sim.chat, xlim = range(breaks), main = head, xlab = 'c-hat',
        breaks = breaks, cex.main = 1, yaxs = 'i')
    abline(v = chat$chat, lwd = 1.5, col = 'blue')
}

Description

Functions to answer the question "what radius is expected to include proportion p of points from a circular bivariate distribution corresponding to a given detection function", and the reverse. These functions may be used to relate the scale parameter(s) of a detection function (e.g., $\sigma$) to home-range area (specifically, the area within an activity contour for the corresponding simple home-range model) (see Note).

WARNING: the default behaviour of these functions changed in version 2.6.0. Integration is now performed on the cumulative hazard (exposure) scale for all functions unless hazard = FALSE. Results will differ.

Usage

circular.r (p = 0.95, detectfn = 0, sigma = 1, detectpar = NULL, hazard
= TRUE, upper = Inf, ...) 

circular.p (r = 1, detectfn = 0, sigma = 1, detectpar = NULL, hazard
= TRUE, upper = Inf, ...)

Arguments

Details

circular.r is the quantile function of the specified circular bivariate distribution (analogous to qnorm, for example). The quantity calculated by circular.r is sometimes called 'circular error probable' (see Note).

For detection functions with two parameters (intercept and scale) it is enough to provide sigma. Otherwise, detectpar should be a named list including parameter values for the requested detection function (g0 may be omitted, and order does not matter).

Detection functions in secr are expressed in terms of the decline in probability of detection with distance $g(d)$, and both circular.r and circular.p integrate this function by default. Rather than integrating $g(d)$ itself, it may be more appropriate to integrate $g(d)$ transformed to a hazard i.e. $1 - log(-g(d))$. This is selected with hazard = TRUE.

Integration may also fail with the message “roundoff error is detected in the extrapolation table”. Setting upper to a large number less than infinity sometimes corrects this.

Value

Vector of values for the required radii or probabilities.

Note

The term ‘circular error probable’ has a military origin. It is commonly used for GPS accuracy with the default probability level set to 0.5 (i.e. half of locations are further than CEP from the true location). A circular bivariate normal distriubution is commonly assumed for the circular error probable; this is equivalent to setting detectfn = "halfnormal".

Closed-form expressions are used for the normal and uniform cases; in the circular bivariate normal case, the relationship is $r = \sigma \sqrt{-2\mbox{ln}(1-p)}$. Otherwise, the probability is computed numerically by integrating the radial distribution. Numerical integration is not foolproof, so check suspicious or extreme values.

When circular.r is used with the default sigma = 1, the result may be interpreted as the factor by which sigma needs to be inflated to include the desired proportion of activity (e.g., 2.45 sigma for 95% of points from a circular bivariate normal distribution fitted on the hazard scale (detectfn = 14) OR 2.24 sigma on the probability scale (detectfn = 0)).

References

Calhoun, J. B. and Casby, J. U. (1958) Calculation of home range and density of small mammals. Public Health Monograph No. 55. United States Government Printing Office.

Johnson, R. A. and Wichern, D. W. (1982) Applied multivariate statistical analysis. Prentice-Hall, Englewood Cliffs, New Jersey, USA.

See Also

detectfn, detectfnplot

Examples

## Calhoun and Casby (1958) p 3.
## give p = 0.3940, 0.8645, 0.9888
circular.p(1:3, hazard = FALSE)

## halfnormal, hazard-rate and exponential
circular.r ()
circular.r (detectfn = "HR", detectpar = list(sigma = 1, z = 4))
circular.r (detectfn = "EX")
circular.r (detectfn = "HHN")
circular.r (detectfn = "HHR", detectpar = list(sigma = 1, z = 4))
circular.r (detectfn = "HEX")

plot(seq(0, 5, 0.05), circular.p(r = seq(0, 5, 0.05)),
    type = "l", xlab = "Radius (multiples of sigma)", ylab = "Probability")
lines(seq(0, 5, 0.05), circular.p(r = seq(0, 5, 0.05), detectfn = 2),
    type = "l", col = "red")
lines(seq(0, 5, 0.05), circular.p(r = seq(0, 5, 0.05), detectfn = 1,
    detectpar = list(sigma = 1,z = 4)), type = "l", col = "blue")
abline (h = 0.95, lty = 2)

legend (2.8, 0.3, legend = c("halfnormal","hazard-rate, z = 4", "exponential"),
    col = c("black","blue","red"), lty = rep(1,3))

## in this example, a more interesting comparison would use
## sigma = 0.58 for the exponential curve.

Description

Clone rows of an object a constant or random number of times

Usage

## Default S3 method:
clone(object, type, ...)
  ## S3 method for class 'popn'
clone(object, type, ...)
  ## S3 method for class 'capthist'
clone(object, type, ...)

Arguments

Details

The ... argument specifies the number of times each row should be repeated. For random distributions (Poisson or negative binomial) ... provides the required parameter values: lambda for Poisson, size, prob or size, mu for negative binomial.

One application is to derive a population of cues from a popn object, where each animal in the original popn generates a number of cues from the same point.

Cloning a capthist object replicates whole detection histories. Individual covariates and detection-specific attributes (e.g., signal strength or xy location in polygon) are also replicated. Cloned data from single-catch traps will cause verify() to fail, but a model may still be fitted in secr.fit by overriding the check with verify = FALSE.

Value

Object of same class as object but with varying number of rows. For clone.popn and capthist an attribute ‘freq’ is set, a vector of length equal to the original number of rows giving the number of repeats (including zeros).

If popn or capthist is a multi-session object the returned value will be a multi-session object of the same length.

See Also

sim.popn

Examples

## population of animals at 1 / hectare generates random
## Poisson number of cues, lambda = 5
mics4 <- make.grid( nx = 2, ny = 2, spacing = 44, detector = "signal")
pop <- sim.popn (D = 1, core = mics4, buffer = 300, nsessions = 6)
pop <- clone (pop, "poisson", 5)
attr(pop[[1]],"freq")

clone(captdata, "poisson", 3)

# To avoid losing any individuals use zero-truncated Poisson
# First find lambda of truncated Poisson with given mean
getlambda <- function (target) {
    fn <- function(x) x / (1-exp(-x)) - target
    uniroot(interval = c(1e-8, target), f = fn)$root
}
clone(captdata, "truncatedpoisson", getlambda(3))

Description

Estimate N, the size of a closed population, by several conventional non-spatial capture–recapture methods.

Usage

closedN(object, estimator = NULL, level = 0.95, maxN = 1e+07,
    dmax = 10 )

Arguments

Details

Data are provided as spatial capture histories, but the spatial information (trapping locations) is ignored.

AIC-based model selection is available for the maximum-likelihood estimators null, zippin, darroch, h2, and betabinomial.

Model weights are calculated as

$w_i = \frac{\exp(-\Delta_i / 2)}{ \sum{\exp(-\Delta_i / 2)}}$

Models for which dAICc > dmax are given a weight of zero and are excluded from the summation, as are non-likelihood models.

Computation of null, zippin and darroch estimates differs slightly from Otis et al. (1978) in that the likelihood is maximized over real values of N between Mt1 and maxN, whereas Otis et al. considered only integer values.

Asymmetric confidence intervals are obtained in the same way for all estimators, using a log transformation of $\hat{N}-Mt1$ following Burnham et al. (1987), Chao (1987) and Rexstad and Burnham (1991).

The available estimators are

Value

A dataframe with one row per estimator and columns

Warning

If your data are from spatial sampling (e.g. grid trapping) it is recommended that you do not use these methods to estimate population size (see Efford and Fewster 2013). Instead, fit a spatial model and estimate population size with region.N.

Note

Prof. Anne Chao generously allowed me to adapt her code for the variance of the ‘chao.th1’ and ‘chao.th2’ estimators.

Chao's estimators have been subject to various improvements not included here (e.g., Chao et al. 2016).

References

Burnham, K. P. and Overton, W. S. (1978) Estimating the size of a closed population when capture probabilities vary among animals. Biometrika 65, 625–633.

Chao, A. (1987) Estimating the population size for capture–recapture data with unequal catchability. Biometrics 43, 783–791.

Chao, A., Ma, K. H., Hsieh, T. C. and Chiu, Chun-Huo (2016) SpadeR: Species-Richness Prediction and Diversity Estimation with R. R package version 0.1.1. https://CRAN.R-project.org/package=SpadeR

Dorazio, R. M. and Royle, J. A. (2003) Mixture models for estimating the size of a closed population when capture rates vary among individuals. Biometrics 59, 351–364.

Efford, M. G. and Fewster, R. M. (2013) Estimating population size by spatially explicit capture–recapture. Oikos 122, 918–928.

Hurvich, C. M. and Tsai, C. L. (1989) Regression and time series model selection in small samples. Biometrika 76, 297–307.

Lee, S.-M. and Chao, A. (1994) Estimating population size via sample coverage for closed capture-recapture models. Biometrics 50, 88–97.

Otis, D. L., Burnham, K. P., White, G. C. and Anderson, D. R. (1978) Statistical inference from capture data on closed animal populations. Wildlife Monographs 62, 1–135.

Pledger, S. (2000) Unified maximum likelihood estimates for closed capture-recapture models using mixtures. Biometrics 56, 434–442.

Rexstad, E. and Burnham, K. (1991) User's guide for interactive program CAPTURE. Colorado Cooperative Fish and Wildlife Research Unit, Fort Collins, Colorado, USA.

See Also

capthist, closure.test, region.N

Examples

closedN(deermouse.ESG)

Description

Perform tests to determine whether a population sampled by capture-recapture is closed to gains and losses over the period of sampling.

Usage

closure.test(object, SB = FALSE, min.expected = 2)

Arguments

Details

The test of Stanley and Burnham in part uses a sum over 2x2 contingency tables; any table with a cell whose expected count is less than min.expected is dropped from the sum. The default value of 2 is that used by CloseTest (Stanley and Richards 2005, T. Stanley pers. comm.; see also Stanley and Burnham 1999 p. 203).

Value

In the case of a single-session capthist object, either a vector with the statistic (z-value) and p-value for the test of Otis et al. (1978 p. 120) or a list whose components are data frames with the statistics and p-values for various tests and test components as follows –

Check the original papers for an explanation of the components of the Stanley and Burnham test.

In the case of a multi-session object, a list with one component (as above) for each session.

Note

No omnibus test exists for closure: the existing tests may indicate nonclosure even when a population is closed if other effects such as trap response are present (see White et al. 1982 pp 96–97). The test of Stanley and Burnham is sensitive to individual heterogeneity which is inevitable in most spatial sampling, and it should not in general be used for this sort of data.

References

Otis, D. L., Burnham, K. P., White, G. C. and Anderson, D. R. (1978) Statistical inference from capture data on closed animal populations. Wildlife Monographs 62, 1–135.

Stanley, T. R. and Burnham, K. P. (1999) A closure test for time-specific capture–recapture data. Environmental and Ecological Statistics 6, 197–209.

Stanley, T. R. and Richards, J. D. (2005) A program for testing capture–recapture data for closure. Wildlife Society Bulletin 33, 782–785.

White, G. C., Anderson, D. R., Burnham, K. P. and Otis, D. L. (1982) Capture-recapture and removal methods for sampling closed populations. Los Alamos National Laboratory, Los Alamos, New Mexico.

See Also

capthist

Examples

closure.test(captdata)

Description

Clusters are uniform groups of detectors. Use these functions to extract or replace cluster information of a traps object, or extract cluster information for each detection in a capthist object.

Usage

clusterID(object)
clusterID(object) <- value
clustertrap(object)
clustertrap(object) <- value

Arguments

Details

Easy access to attributes used to define compound designs, those in which a detector array comprises several similar subunits (‘clusters’). ‘clusterID’ identifies the detectors belonging to each cluster, and ‘clustertrap’ is a numeric index used to relate matching detectors in different clusters.

For replacement (‘traps’ only), the number of rows of value must match exactly the number of detectors in object.

‘clusterID’ and ‘clustertrap’ are assigned automatically by trap.builder.

Value

Factor (clusterID) or integer-valued vector (clustertrap).

clusterID(object) may be NULL.

See Also

traps, trap.builder, mash, derivedCluster, cluster.counts, cluster.centres

Examples

## 25 4-detector clusters
mini <- make.grid(nx = 2, ny = 2)
tempgrid <- trap.builder (cluster = mini , method = "all",
    frame = expand.grid(x = seq(100, 500, 100), y = seq(100,
    500, 100)))
clusterID(tempgrid)
clustertrap(tempgrid)

tempCH <- sim.capthist(tempgrid)
table(clusterID(tempCH)) ## detections per cluster
cluster.counts(tempCH)   ## distinct individuals

Description

Extract coefficients (estimated beta parameters) from a spatially explicit capture–recapture model.

Usage

## S3 method for class 'secr'
 coef(object, alpha = 0.05, ...)

Arguments

Value

A data frame with one row per beta parameter and columns for the coefficient, SE(coefficient), asymptotic lower and upper 100(1–alpha) confidence limits.

See Also

secr.fit, esaPlot

Examples

## load & extract coefficients of previously fitted null model
coef(secrdemo.0)

Description

Estimates from one or more openCR models are formed into an array.

Usage

## S3 method for class 'secr'
collate(object, ..., realnames = NULL, betanames = NULL, newdata = NULL, 
    alpha = 0.05, perm = 1:4, fields = 1:4)

## S3 method for class 'ipsecr'
collate(object, ..., realnames = NULL, betanames = NULL, newdata = NULL, 
    alpha = 0.05, perm = 1:4, fields = 1:4)

## S3 method for class 'secrlist'
collate(object, ..., realnames = NULL, betanames = NULL, newdata = NULL, 
    alpha = 0.05, perm = 1:4, fields = 1:4)

Arguments

Details

collate extracts parameter estimates from a set of fitted secr model objects.

fields may be used to select a subset of summary fields ("estimate","SE.estimate","lcl","ucl") by name or number.

Value

A 4-dimensional array of model-specific parameter estimates. By default, the dimensions correspond respectively to

• rows in newdata (usually sessions),

• models,

• statistic fields (estimate, SE.estimate, lcl, ucl), and

• parameters ("phi", "sigma" etc.).

It often helps to reorder the dimensions with the perm argument.

See Also

modelAverage, secr.fit

Examples

collate (secrdemo.0, secrdemo.b, perm = c(4,2,3,1))[,,1,]

Description

Compute profile likelihood confidence intervals for ‘beta’ or ‘real’ parameters of a spatially explicit capture-recapture model,

Usage

## S3 method for class 'secr'
 confint(object, parm, level = 0.95, newdata = NULL,
tracelevel = 1, tol = 0.0001, bounds = NULL, ncores = NULL, ...)

Arguments

Details

If parm is numeric its elements are interpreted as the indices of ‘beta’ parameters; character values are interpreted as ‘real’ parameters. Different methods are used for beta parameters and real parameters. Limits for the $j$-th beta parameter are found by a numerical search for the value satisfying $-2(l_j(\beta_j) - l) = q$, where $l$ is the maximized log likelihood, $l_j(\beta_j)$ is the maximized profile log likelihood with $\beta_j$ fixed, and $q$ is the $100(1-\alpha)$ quantile of the $\chi^2$ distribution with one degree of freedom. Limits for real parameters use the method of Lagrange multipliers (Fletcher and Faddy 2007), except that limits for constant real parameters are backtransformed from the limits for the relevant beta parameter.

If bounds is provided it should be a 2-vector or matrix of 2 columns and length(parm) rows.

Setting ncores = NULL uses the existing value from the environment variable RCPP_PARALLEL_NUM_THREADS (see setNumThreads).

Value

A matrix with one row for each parameter in parm, and columns giving the lower (lcl) and upper (ucl) 100*level

Note

Calculation may take a long time, so probably you will do it only after selecting a final model.

The R function uniroot is used to search for the roots of $-2(l_j(\beta_j) - l) = q$ within a suitable interval. The interval is anchored at one end by the MLE, and at the other end by the MLE inflated by a small multiple of the asymptotic standard error (1, 2, 4 or 8 SE are tried in turn, using the smallest for which the interval includes a valid solution).

A more efficient algorithm was proposed by Venzon and Moolgavkar (1988); it has yet to be implemented in secr, but see plkhci in the package Bhat for another R implementation.

References

Evans, M. A., Kim, H.-M. and O'Brien, T. E. (1996) An application of profile-likelihood based confidence interval to capture–recapture estimators. Journal of Agricultural, Biological and Experimental Statistics 1, 131–140.

Fletcher, D. and Faddy, M. (2007) Confidence intervals for expected abundance of rare species. Journal of Agricultural, Biological and Experimental Statistics 12, 315–324.

Venzon, D. J. and Moolgavkar, S. H. (1988) A method for computing profile-likelihood-based confidence intervals. Applied Statistics 37, 87–94.

Examples

## Not run: 

## Limits for the constant real parameter "D"
confint(secrdemo.0, "D")   


## End(Not run)

Description

Display contours of the net probability of detection p.(X), or the area within a specified distance of detectors. bufferContour adds a conventional ‘boundary strip’ to a detector (trap) array, where buffer equals the strip width.

Usage

pdotContour(traps, border = NULL, nx = 64, detectfn = 0,
    detectpar = list(g0 = 0.2, sigma = 25, z = 1), noccasions = NULL,
    binomN = NULL, levels = seq(0.1, 0.9, 0.1), poly =
    NULL, poly.habitat = TRUE, plt = TRUE, add = FALSE, fill = NULL, ...)

bufferContour(traps, buffer, nx = 64, convex = FALSE, ntheta = 100,
     plt = TRUE, add = FALSE, poly = NULL, poly.habitat = TRUE, 
     fill = NULL, ...)

Arguments

Details

pdotContour constructs a rectangular mask and applies pdot to compute the p.(X) at each mask point.

If convex = FALSE, bufferContour constructs a mask and contours the points on the basis of distance to the nearest detector at the levels given in buffer.

If convex = TRUE, bufferContour constructs a set of potential vertices by adding points on a circle of radius = buffer to each detector location; the desired contour is the convex hull of these points (this algorithm derives from Efford, 2012).

If traps has a usage attribute then noccasions is set accordingly; otherwise it must be provided.

If traps is for multiple sessions then detectpar should be a list of the same length, one component per session, and noccasions may be a numeric vector of the same length.

Increase nx for smoother lines, at the expense of speed.

Value

Coordinates of the plotted contours are returned as a list with one component per polygon. The list is returned invisibly if plt = TRUE.

For multi-session input (traps) the value is a list of such lists, one per session.

Note

The precision (smoothness) of the fitted line in bufferContour is controlled by ntheta rather than nx when convex = TRUE.

To suppress contour labels, include the argument drawlabels = FALSE (this will be passed via ... to contour). Other useful arguments of contour are col (colour of contour lines) and lwd (line width).

You may wish to consider function st_buffer in package sf as an alternative to bufferContour.

bufferContour failed with multi-session traps before secr 2.8.0.

References

Efford, M. G. (2012) DENSITY 5.0: software for spatially explicit capture–recapture. Department of Mathematics and Statistics, University of Otago, Dunedin, New Zealand https://www.otago.ac.nz/density/.

See Also

pdot, make.mask

Examples

possumtraps <- traps(possumCH)

## convex and concave buffers
plot(possumtraps, border = 270)
bufferContour(possumtraps, buffer = 100, add = TRUE, col = "blue")
bufferContour(possumtraps, buffer = 100, convex = TRUE, add = TRUE)

## areas
buff.concave <- bufferContour(possumtraps, buffer = 100,
    plt = FALSE)
buff.convex  <- bufferContour(possumtraps, buffer = 100,
    plt = FALSE, convex = TRUE)
sum (sapply(buff.concave, polyarea)) ## sum over parts
sapply(buff.convex, polyarea)

## effect of nx on area
buff.concave2 <- bufferContour(possumtraps, buffer = 100,
    nx = 128, plt = FALSE)
sum (sapply(buff.concave2, polyarea))

## Not run: 

plot(possumtraps, border = 270)
pdotContour(possumtraps, detectfn = 0, nx = 128, detectpar =
    detectpar(possum.model.0), levels = c(0.1, 0.01, 0.001),
    noccasions = 5, add = TRUE)

## clipping to polygon
olddir <- setwd(system.file("extdata", package = "secr"))
possumtraps <- traps(possumCH)
possumarea <- read.table("possumarea.txt", header = TRUE)
par(xpd = TRUE, mar = c(1,6,6,6))
plot(possumtraps, border = 400, gridlines = FALSE)
pdotContour(possumtraps, detectfn = 0, nx = 256, detectpar =
    detectpar(possum.model.0), levels = c(0.1, 0.01, 0.001),
    noccasions = 5, add = TRUE, poly = possumarea, col = "blue")
lines(possumarea)
setwd(olddir)
par(xpd = FALSE, mar = c(5,4,4,2) + 0.1)    ## reset to default


## End(Not run)

Description

Extract or replace covariates

Usage

covariates(object, ...)
covariates(object) <- value

Arguments

Details

For replacement, the number of rows of value must match exactly the number of rows in object.

Value

covariates(object) returns the dataframe of covariates associated with object. covariates(object) may be NULL.

Individual covariates are stored in the ‘covariates’ attribute of a capthist object.

Covariates used for modelling density are stored in the ‘covariates’ attribute of a mask object.

Detector covariates may vary between sampling occasions. In this case, columns in the detector covariates data.frame are associated with particular times; the matching is controlled by the timevaryingcov attribute.

See Also

timevaryingcov

Examples

## detector covariates
temptrap <- make.grid(nx = 6, ny = 8)
covariates (temptrap) <- data.frame(halfnhalf = 
    factor(rep(c("left","right"),c(24,24))) )
summary(covariates(temptrap))

Description

The coefficient of variation of effective sampling area predicts the bias in estimated density (Efford and Mowat 2014). These functions assist its calculation from fitted finite mixture models.

Usage

CV(x, p, na.rm = FALSE)
CVa0(object, ...)
CVa(object, sessnum = 1, ...)

Arguments

Details

CV computes the coefficient of variation of x. If p is provided then the distribution is assumed to be discrete, with support x and class membership probabilities p (scaled automatically to sum to 1.0).

CVa computes CV($a$) where $a$ is the effective sampling area of Borchers and Efford (2008).

CVa0 computes CV(a0) where a0 is the single-detector sampling area defined as $a_0 = 2 \pi \lambda_0 \sigma^2$ (Efford and Mowat 2014); a0 is a convenient surrogate for a, the effective sampling area. CV(a0) uses either the fitted MLE of a0 (if the a0 parameterization has been used), or a0 computed from the estimates of lambda0 and sigma.

CVa and CVa0 do not work for models with individual covariates.

Value

Numeric

Note

Do not confuse the function CVa with the estimated relative standard error of the estimate of a from derived, also labelled CVa in the output. The relative standard error RSE is often labelled CV in the literature on capture–recapture, but this can cause unnecessary confusion. See also RSE.

References

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

Efford, M. G. and Mowat, G. (2014) Compensatory heterogeneity in capture–recapture data. Ecology 95, 1341–1348.

See Also

CVpdot, derived, details, RSE

Examples

## Not run: 

## housemouse model
morning <- subset(housemouse, occ = c(1,3,5,7,9))
msk <- make.mask((traps(morning)), nx = 32) 
morning.h2   <- secr.fit(morning, buffer = 20, model = list(g0~h2), mask = msk, 
    trace = FALSE)
CVa0(morning.h2 )


## End(Not run)

Description

Internal function used by secr.fit, confint.secr, and score.test.

Usage

D.designdata (mask, Dmodel, grouplevels, sessionlevels, sessioncov =
NULL, meanSD = NULL)

Arguments

Details

This is an internal secr function that you are unlikely ever to use. Unlike secr.design.MS, this function does not call model.matrix.

Value

Dataframe with one row for each combination of mask point, group and session. Conceptually, we use a 3-D rectangular array with enough rows to accommodate the largest mask, so some rows in the output may merely hold space to enable easy indexing. The dataframe has an attribute ‘dimD’ that gives the relevant dimensions: attr(dframe, "dimD") = c(nmask, ngrp, R), where nmask is the number of mask points, ngrp is the number of groups, and R is the number of sessions. Columns correspond to predictor variables in Dmodel.

The number of valid rows (points in each session-specific mask) is stored in the attribute ‘validMaskRows’.

For a single-session mask, meanSD is a 2 x 2 matrix of mean and SD (rows) for x- and y-coordinates. For a multi-session mask, a list of such objects. Ordinarily these values are from the meanSD attribute of the mask, but they must be specified when applying a new mask to an existing model.

See Also

secr.design.MS

Description

Data of V. H. Reid from live trapping of deermice (Peromyscus maniculatus) at two sites in Colorado, USA.

Usage

deermouse.ESG
deermouse.WSG

Details

Two datasets of V. H. Reid were described by Otis et al. (1978) and distributed with their CAPTURE software (now available from https://eesc.usgs.gov/mbr/software/capture.shtml). They have been used in several other papers on closed population methods (e.g., Huggins 1991, Stanley and Richards 2005). This description is based on pages 32 and 87–93 of Otis et al. (1978).

Both datasets are from studies in Rio Blanco County, Colorado, in the summer of 1975. Trapping was for 6 consecutive nights. Traps were arranged in a 9 x 11 grid and spaced 50 feet (15.2 m) apart.

The first dataset was described by Otis et al. (1978: 32) as from 'a drainage bottom of sagebrush, gambel oak, and serviceberry with pinyon pine and juniper on the uplands'. By matching with the ‘examples’ file of CAPTURE this was from East Stuart Gulch (ESG).

The second dataset (Otis et al. 1978: 87) was from Wet Swizer Creek or Gulch (WSG) in August 1975. No specific vegetation description is given for this site, but it is stated that Sherman traps were used and trapping was done twice daily.

Two minor inconsistencies should be noted. Although Otis et al. (1978) said they used data from morning trap clearances, the capture histories in ‘examples’ from CAPTURE include a ‘pm’ tag on each record. We assume the error was in the text description, as their numerical results can be reproduced from the data file. Huggins (1991) reproduced the East Stuart Gulch dataset (omitting spatial data that were not relevant to his method), but omitted two capture histories.

The data are provided as two single-session capthist objects ‘deermouse.ESG’ and ‘deermouse.WSG’. Each has a dataframe of individual covariates, but the fields differ between the two study areas. The individual covariates of deermouse.ESG are sex (factor levels ‘f’, ‘m’), age class (factor levels ‘y’, ‘sa’, ‘a’) and body weight in grams. The individual covariates of deermouse.WSG are sex (factor levels ‘f’,‘m’) and age class (factor levels ‘j’, ‘y’, ‘sa’, ‘a’) (no data on body weight). The aging criteria used by Reid are not recorded.

The datasets were originally in the CAPTURE ‘xy complete’ format which for each detection gives the ‘column’ and ‘row’ numbers of the trap (e.g. ‘ 9 5’ for a capture in the trap at position (x=9, y=5) on the grid). Trap identifiers have been recoded as strings with no spaces by inserting zeros (e.g. ‘905’ in this example).

Sherman traps are designed to capture one animal at a time, but the data include double captures (1 at ESG and 8 at WSG – see Examples). The true detector type therefore falls between ‘single’ and ‘multi’. Detector type is set to ‘multi’ in the distributed data objects.

Source

File ‘examples’ distributed with program CAPTURE.

References

Huggins, R. M. (1991) Some practical aspects of a conditional likelihood approach to capture experiments. Biometrics 47, 725–732.

Otis, D. L., Burnham, K. P., White, G. C. and Anderson, D. R. (1978) Statistical inference from capture data on closed animal populations. Wildlife Monographs 62, 1–135.

Stanley, T. R. and Richards, J. D. (2005) A program for testing capture–recapture data for closure. Wildlife Society Bulletin 33, 782–785.

See Also

closure.test

Examples

par(mfrow = c(1,2), mar = c(1,1,4,1))
plot(deermouse.ESG, title = "Peromyscus data from East Stuart Gulch",
    border = 10, gridlines = FALSE, tracks = TRUE)
plot(deermouse.WSG, title = "Peromyscus data from Wet Swizer Gulch",
    border = 10, gridlines = FALSE, tracks = TRUE)

closure.test(deermouse.ESG, SB = TRUE)

## reveal multiple captures
table(trap(deermouse.ESG), occasion(deermouse.ESG))
table(trap(deermouse.WSG), occasion(deermouse.WSG))

Description

Mask points may be removed by one of three methods: clicking on points, clicking on vertices to define a polygon from which points will be removed, or specifying a polygon to which the mask will be clipped.

Usage

deleteMaskPoints(mask, onebyone = TRUE, add = FALSE, poly = NULL,
 poly.habitat = FALSE, ...)

Arguments

Details

The default method (onebyone = TRUE, poly = NULL) is to click on each point to be removed. The nearest mask point will be selected.

Setting onebyone = FALSE allows the user to click on the vertices of a polygon within which all points are to be removed (the default) or retained (poly.habitat = TRUE). Vertices need not coincide with mask points.

Defining poly here is equivalent to calling make.mask with poly defined. poly one of the several types described in boundarytoSF. Whether poly represents habitat or non-habitat is toggled with poly.habitat – the default here differs from make.mask.

Value

A mask object, usually with fewer points than the input mask.

See Also

make.mask, subset.mask

Examples

if (interactive()) {
    mask0 <- make.mask (traps(captdata))
    ## Method 1 - click on each point to remove
    mask1 <- deleteMaskPoints (mask0)
    ## Method 2 - click on vertices of removal polygon
    mask2 <- deleteMaskPoints (mask0, onebyone = FALSE)
    ## Method 3 - predefined removal polygon
    plot(captdata)
    poly1 <- locator(5)
    mask3 <- deleteMaskPoints (mask0, poly = poly1)
}

Description

Compute derived parameters of spatially explicit capture-recapture model. Density is a derived parameter when a model is fitted by maximizing the conditional likelihood. So also is the effective sampling area (in the sense of Borchers and Efford 2008).

Usage

derived(object, ...)

## S3 method for class 'secr'
derived(object, sessnum = NULL, groups = NULL, alpha = 0.05, 
    se.esa = FALSE, se.D = TRUE, loginterval = TRUE, distribution = NULL, 
    ncores = NULL, bycluster = FALSE, ...)

## S3 method for class 'secrlist'
derived(object, sessnum = NULL, groups = NULL, alpha = 0.05, 
    se.esa = FALSE, se.D = TRUE, loginterval = TRUE, distribution = NULL, 
    ncores = NULL, bycluster = FALSE, ...)

## S3 method for class 'secr'
esa(object, sessnum = 1, beta = NULL, real = NULL, noccasions = NULL, 
    ncores = NULL, ...)

Arguments

Details

The derived estimate of density is a Horvitz-Thompson-like estimate:

$\hat{D} = \sum\limits _{i=1}^{n} {a_i (\hat{\theta})^{-1}}$

where $a_i (\hat{\theta})$ is the estimate of effective sampling area for animal $i$ with detection parameter vector $\theta$.

A non-null value of the argument distribution overrides the value in object$details. The sampling variance of $\hat{D}$ from secr.fit by default is spatially unconditional (distribution = "Poisson"). For sampling variance conditional on the population of the habitat mask (and therefore dependent on the mask area), specify distribution = "binomial". The equation for the conditional variance includes a factor $(1 - a/A)$ that disappears in the unconditional (Poisson) variance (Borchers and Efford 2007). Thus the conditional variance is always less than the unconditional variance. The unconditional variance may in turn be an overestimate or (more likely) an underestimate if the true spatial variance is non-Poisson.

Derived parameters may be estimated for population subclasses (groups) defined by the user with the groups argument. Each named factor in groups should appear in the covariates dataframe of object$capthist (or each of its components, in the case of a multi-session dataset).

esa is used by derived to compute individual-specific effective sampling areas:

$a_i (\hat{\theta}) = \int _A \: p.(\mathbf{X};\mathbf{z}_i, \mathbf{\hat{\theta}}) \; \mathrm{d} \mathbf{X}$

where $p.(\mathbf{X})$ is the probability an individual at X is detected at least once and the $\mathbf{z}_i$ are optional individual covariates. Integration is over the area $A$ of the habitat mask.

The argument noccasions may be used to vary the number of sampling occasions; it works only when detection parameters are constant across individuals and across time.

Setting ncores = NULL uses the existing value from the environment variable RCPP_PARALLEL_NUM_THREADS (see setNumThreads).

The effective sampling area ‘esa’ ($a(\hat{\theta})$) reported by derived is equal to the harmonic mean of the $a_i (\hat{\theta})$ (arithmetic mean prior to version 1.5). The sampling variance of $a(\hat{\theta})$ is estimated by

$\widehat{\mbox{var}}(a(\hat{\theta})) = \hat{G}_\theta^T \hat{V}_\theta \hat{G}_\theta,$

where $\hat{V}$ is the asymptotic estimate of the variance-covariance matrix of the beta detection parameters ($\theta$) and $\hat{G}$ is a numerical estimate of the gradient of $a(\theta)$ with respect to $\theta$, evaluated at $\hat{\theta}$.

A 100(1–alpha)% asymptotic confidence interval is reported for density. By default, this is asymmetric about the estimate because the variance is computed by backtransforming from the log scale. You may also choose a symmetric interval (variance obtained on natural scale).

The vector of detection parameters for esa may be specified via beta or real, with the former taking precedence. If neither is provided then the fitted values in object$fit$par are used. Specifying real parameter values bypasses the various linear predictors. Strictly, the ‘real’ parameters are for a naive capture (animal not detected previously).

The computation of sampling variances is relatively slow and may be suppressed with se.esa and se.D as desired.

For computing derived across multiple models in parallel see par.derived.

From secr 5.0.0 the ... argument may be used to control the step size (.relStep) used by fdHess when estimating gradients for SE(D) and SE(esa).

Value

Dataframe with one row for each derived parameter (‘esa’, ‘D’) and columns as below

For a multi-session or multi-group analysis the value is a list with one component for each session and group.

The result will also be a list if object is an ‘secrlist’.

Warning

derived() may be applied to detection models fitted by maximizing the full likelihood (CL = FALSE). However, models using g (groups) will not be handled correctly.

Note

Before version 2.1, the output table had columns for ‘varcomp1’ (the variance in $\hat{D}$ due to variation in $n$, i.e., Huggins' $s^2$), and ‘varcomp2’ (the variance in $\hat{D}$ due to uncertainty in estimates of detection parameters).

These quantities are related to CVn and CVa as follows:

$\mbox{CVn} = \sqrt{ \mbox{varcomp1} } / \hat{D}$

$\mbox{CVa} = \sqrt{ \mbox{varcomp2} } / \hat{D}$

References

Borchers, D. L. and Efford, M. G. (2007) Supplements to Biometrics paper. Available online at https://www.otago.ac.nz/density/.

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

Huggins, R. M. (1989) On the statistical analysis of capture experiments. Biometrika 76, 133–140.

See Also

predict.secr, print.secr, secr.fit, empirical.varD

Examples

## Not run: 
## extract derived parameters from a model fitted previously
## by maximizing the conditional likelihood 
derived (secrdemo.CL)

## what happens when sampling variance is conditional on mask N?
derived(secrdemo.CL, distribution = "binomial")
## fitted g0, sigma
esa(secrdemo.CL)
## force different g0, sigma
esa(secrdemo.CL, real = c(0.2, 25))

## End(Not run)

Description

The function secr.fit allows many options. Some of these are used infrequently and have been bundled as a single argument details to simplify the documentation. They are described here.

Detail components

details$autoini specifies the session number from which to compute starting values (multi-session data only; default 1). From 4.0.0, the character value ‘all’ first forms a single-session capthist using join(); this may be slow or not work at all (especially with telemetry data).

details$centred = TRUE causes coordinates of both traps and mask to be centred on the centroid of the traps, computed separately for each session in the case of multi-session data. This may be necessary to overcome numerical problems when x- or y-coordinates are large numbers. The default is not to centre coordinates.

details$chat optionally specifies the overdispersion of unmarked sightings Tu and unidentified marked sightings Tm. It is used only for mark-resight models, and is usually computed within secr.fit (details$nsim > 0), but may be provided by the user. For a single session 'chat' is a vector of length 2; for multiple sessions it is a 2-column matrix.

details$chatonly = TRUE used with details$nsim > 0 causes the overdispersion statistics for sighting counts Tu and Tm to be estimated and returned as a vector or 2-column matrix (multi-session models), with no further model fitting.

details$contrasts may be used to specify the coding of factor predictors. The value should be suitable for the 'contrasts.arg' argument of model.matrix. See ‘Trend across sessions’ in secr-multisession.pdf for an example.

details$convexpolygon may be set to FALSE for searches of non-convex polygons. This is slower than the default which requires poygons to be convex east-west (secr-polygondetectors.pdf).

details$debug is an integer code used to control the printing of intermediate values (1,2) and to switch on the R code browser (3). In ordinary use it should not be changed from the default (0).

details$Dfn is a function for reparameterizing density models; this is set internally when Dlambda = TRUE. Exotic variations may be specified directly by the user when Dlambda = FALSE. The defaults (Dfn = NULL, Dlambda = FALSE) leave the original density model unchanged. Note there is no connection to userDfn (except that the two are incompatible).

Dlambda if TRUE causes reparameterization of density as the session-on-session finite rate of increase $lambda$. Details at (secr-trend.pdf).

details$distribution specifies the distribution of the number of individuals detected $n$; this may be conditional on the number in the masked area ("binomial") or unconditional ("poisson"). distribution affects the sampling variance of the estimated density. The default is "poisson". The component ‘distribution’ may also take a numeric value larger than nrow(capthist), rather than "binomial" or "poisson". The likelihood then treats n as a binomial draw from a superpopulation of this size, with consequences for the variance of density estimates. This can help to reconcile MLE with Bayesian estimates using data augmentation.

details$fastproximity controls special handling of data from binary proximity and count detectors. If TRUE and other conditions are met (no temporal variation or groups) then a multi-occasion capthist is automatically reduced to a count for a single occasion and further compressed by storing only non-zero counts, which can greatly speed up computation of the likelihood (default TRUE).

details$fixedbeta may be used to fix values of beta parameters. It should be a numeric vector of length equal to the total number of beta parameters (coefficients) in the model. Parameters to be estimated are indicated by NA. Other elements should be valid values on the link scale and will be substituted during likelihood maximisation. Check the order of beta parameters in a previously fitted model.

details$grain sets the grain argument for multithreading in RcppParallel parallelFor (default 1). details$grain = 0 suppresses multithreading (equivalent to ncores = 1).

details$hessian is a character string controlling the computation of the Hessian matrix from which variances and covariances are obtained. Options are "none" (no variances), "auto" (the default) or "fdhess" (use the function fdHess in nlme). If "auto" then the Hessian from the optimisation function is used. See also method = "none" below.

details$ignoreusage = TRUE causes the function to ignore usage (varying effort) information in the traps component. The default (details$ignoreusage = FALSE) is to include usage in the model.

details$intwidth2 controls the half-width of the interval searched by optimise() for the maximum likelihood when there is a single parameter. Default 0.8 sets the search interval to $(0.2s, 1.8s)$ where $s$ is the ‘start’ value.

details$knownmarks = FALSE causes secr.fit to fit a zero-truncated sightings-only model that implicitly estimates the number of marked individuals, rather than inferring it from the number of rows in the capthist object.

details$LLonly = TRUE causes the function to returns a single evaluation of the log likelihood at the ‘start’ values.

details$maxdistance sets a limit to the centroid-to-mask distances considered. The centroid is the geometric mean of detection locations for each individual. If no limit is specified then summation is over all mask points. Specifying maxdistance can speed up computation; it is up to the user to select a limit that is large enough not to affect the likelihood ($5\sigma$?).

details$miscparm (default NULL) is an optional numeric vector of starting values for additional parameters used in a user-supplied distance function (see ‘userdist’ below). If the vector has a names attribute then the names will be used for the corresponding coefficients (‘beta’ parameters) which will otherwise be named ‘miscparm1’, miscparm2' etc. These parameters are constant across each model and do not appear in the model formula, but are estimated along with other coefficients when the likelihood is maximised. Any transformation (link function) etc. is handled by the user in the userdist function. The coefficients appear in the output from coef.secr and vcov.secr, but not predict.secr.

details$newdetector specifies a detector type to use for this fit, replacing the previous detector(traps(capthist)). The value may be a vector (one value per occasion) or for multi-session data, a list of vectors. A scalar value (e.g. "proximity") is otherwise used for all occasions and sessions. The true detector type is usually known and will be specified in the 'traps' attribute; newdetector is useful in simulation studies that examine the effect of misspecification. The capthist component of the output from secr.fit has the new type.

details$nsim specifies the number of replicate simulations to perform to estimate the overdispersion statistics for the sighting counts Tu and Tm. See also details$chat and details$chatonly.

details$param chooses between various parameterisations of the SECR model. The default details$param = 0 is the formulation in Borchers and Efford (2008) and later papers.

details$param = 1 was once used to select the Gardner & Royle parameterisation of the detection model (p0, $\sigma$; Gardner et al. 2009) when the detector type is ‘multi’. This parameterisation was discontinued in 2.10.0.

details$param = 2 selects parameterisation in terms of ($esa(g_0, \sigma)$, $\sigma$) (Efford and Mowat 2014).

details$param = 3 selects parameterisation in terms of ($a_0(\lambda_0, \sigma)$, $\sigma$) (Efford and Mowat 2014). This parameterization is used automatically if a0 appears in the model (e.g., a0 ~ 1).

details$param = 4 selects parameterisation of sigma in terms of the coefficient sigmak and constant c (sigma = sigmak / D^0.5 + c) (Efford et al. 2016). If c is not included explicitly in the model (e.g., c ~ 1) then it is set to zero. This parameterization is used automatically if sigmak appears in the model (e.g., sigmak ~ 1)

details$param = 5 combines parameterisations (3) and (4) (first compute sigma from D, then compute lambda0 from sigma).

details$relativeD fits a density model conditional on $n$ that describes relative density instead of absolute density. This describes the distribution of tagged animals.

details$savecall determines whether the full call to secr.fit is saved in the output object. The default is TRUE except when called by list.secr.fit as names in the call are then evaluated, causing the output to become unwieldy.

details$splitmarked determines whether the home range centre of marked animals is allowed to move between the marking and sighting phases of a spatial capture–mark–resight study. The default is to assume a common home-range centre (splitmarked = FALSE).

details$telemetrytype determines how telemetry data in the attribute ‘xylist’ are treated. ‘none’ causes the xylist data to be ignored. ‘dependent’ uses information on the sampling distribution of each home-range centre in the SECR likelihood. ‘concurrent’ does that and more: it splits capthist according to telemetry status and appends all-zero histories to the telemetry part for any animals present in xylist. The default is ‘concurrent’.

details$usecov selects the mask covariate to be used for normalization. NULL limits denominator for normalization to distinguishing habitat from non-habitat.

details$userDfn is a user-provided function for modelling a density surface. See secr-densitysurfaces.pdf

details$userdist is either a function to compute non-Euclidean distances between detectors and mask points, or a pre-computed matrix of such distances. The first two arguments of the function should be 2-column matrices of x-y coordinates (respectively $k$ detectors and $m$ mask points). The third argument is a habitat mask that defines a non-Euclidean habitat geometry (a linear geometry is described in documentation for the package ‘secrlinear’). The matrix returned by the function must have exactly $k$ rows and $m$ columns. When called with no arguments the function should return a character vector of names for the required covariates of ‘mask’, possibly including the dynamically computed density 'D' and a parameter ‘noneuc’ that will be fitted. A slightly expanded account is at userdist, and full documentation is in the separate document secr-noneuclidean.pdf.

**Do not use ‘userdist’ for polygon or transect detectors**

References

Efford, M. G., Dawson, D. K., Jhala, Y. V. and Qureshi, Q. (2016) Density-dependent home-range size revealed by spatially explicit capture–recapture. Ecography 39, 676–688.

Efford, M. G. and Mowat, G. (2014) Compensatory heterogeneity in capture–recapture data.Ecology 95, 1341–1348.

Gardner, B., Royle, J. A. and Wegan, M. T. (2009) Hierarchical models for estimating density from DNA mark-recapture studies. Ecology 90, 1106–1115.

Royle, J. A., Chandler, R. B., Sun, C. C. and Fuller, A. K. (2013) Integrating resource selection information with spatial capture–recapture. Methods in Ecology and Evolution 4, 520–530.

See Also

secr.fit , userdist

Examples

## Not run: 

## Demo of miscparm and userdist
## We fix the usual 'sigma' parameter and estimate the same 
## quantity as miscparm[1]. Differences in CI reflect the implied use 
## of the identity link for miscparm[1]. 

mydistfn3 <- function (xy1,xy2, mask) {
    if (missing(xy1)) return(character(0))
    xy1 <- as.matrix(xy1)
    xy2 <- as.matrix(xy2)
    miscparm <- attr(mask, 'miscparm')
    distmat <- edist(xy1,xy2) / miscparm[1]
    distmat
}

fit0 <- secr.fit (captdata)
fit <- secr.fit (captdata, fixed = list(sigma=1), details = 
    list(miscparm = c(sig = 20), userdist = mydistfn3))    
predict(fit0)
coef(fit)


## End(Not run)

Description

A detection function relates the probability of detection $g$ or the expected number of detections $\lambda$ for an animal to the distance of a detector from a point usually thought of as its home-range centre. In secr only simple 2- or 3-parameter functions are used. Each type of function is identified by its number or by a 2–3 letter code (version $\ge$ 2.6.0; see below). In most cases the name may also be used (as a quoted string).

Choice of detection function is usually not critical, and the default ‘HN’ is usually adequate.

Functions (14)–(20) are parameterised in terms of the expected number of detections $\lambda$, or cumulative hazard, rather than probability. ‘Exposure’ (e.g. Royle and Gardner 2011) is another term for cumulative hazard. This parameterisation is natural for the ‘count’ detector type or if the function is to be interpreted as a distribution of activity (home range). When one of the functions (14)–(19) is used to describe detection probability (i.e., for the binary detectors ‘single’, ‘multi’,‘proximity’,‘polygonX’ or ‘transectX’), the expected number of detections is internally transformed to a binomial probability using $g(d) = 1-\exp(-\lambda(d))$.

The hazard halfnormal (14) is similar to the halfnormal exposure function used by Royle and Gardner (2011) except they omit the factor of 2 on $\sigma^2$, which leads to estimates of $\sigma$ that are larger by a factor of sqrt(2). The hazard exponential (16) is identical to their exponential function.