gdistsamp.RdExtends the distance sampling model of Royle et al. (2004) to estimate the probability of being available for detection. Also allows abundance to be modeled using the negative binomial distribution.
gdistsamp(lambdaformula, phiformula, pformula, data, keyfun = c("halfnorm", "exp", "hazard", "uniform"), output = c("abund", "density"), unitsOut = c("ha", "kmsq"), mixture = c("P", "NB"), K, starts, method = "BFGS", se = TRUE, engine=c("C","R"), rel.tol=1e-4, threads=1, ...)
| lambdaformula | A right-hand side formula describing the abundance covariates. |
|---|---|
| phiformula | A right-hand side formula describing the availability covariates. |
| pformula | A right-hand side formula describing the detection function covariates. |
| data | An object of class |
| keyfun | One of the following detection functions: "halfnorm", "hazard", "exp", or "uniform." See details. |
| output | Model either "density" or "abund" |
| unitsOut | Units of density. Either "ha" or "kmsq" for hectares and square kilometers, respectively. |
| mixture | Either "P" or "NB" for the Poisson and negative binomial models of abundance. |
| K | An integer value specifying the upper bound used in the integration. |
| starts | A numeric vector of starting values for the model parameters. |
| method | Optimization method used by |
| se | logical specifying whether or not to compute standard errors. |
| engine | Either "C" to use fast C++ code or "R" to use native R code during the optimization. |
| rel.tol | relative accuracy for the integration of the detection function. See integrate. You might try adjusting this if you get an error message related to the integral. Alternatively, try providing different starting values. |
| threads | Set the number of threads to use for optimization in C++, if
OpenMP is available on your system. Increasing the number of threads
may speed up optimization in some cases by running the likelihood
calculation in parallel. If |
| ... | Additional arguments to optim, such as lower and upper bounds |
This model extends the model of Royle et al. (2004) by estimating the
probability of being available for detection \(\phi\). This
effectively relaxes the assumption that \(g(0)=1\). In other words,
inividuals at a distance of 0 are not assumed to be detected with
certainty. To estimate this additional parameter, replicate distance
sampling data must be collected at each transect. Thus the data are
collected at i = 1, 2, ..., R transects on t = 1, 2, ..., T
occassions. As with the model of Royle et al. (2004), the detections
must be binned into distance classes. These data must be formatted in
a matrix with R rows, and JT columns where J is the number of distance
classses. See unmarkedFrameGDS for more information.
If you aren't interested in estimating phi, but you want to use the negative binomial distribution, simply set numPrimary=1 when formatting the data.
An object of class unmarkedFitGDS.
Royle, J. A., D. K. Dawson, and S. Bates. 2004. Modeling abundance effects in distance sampling. Ecology 85:1591-1597.
Chandler, R. B, J. A. Royle, and D. I. King. 2011. Inference about density and temporary emigration in unmarked populations. Ecology 92:1429--1435.
Richard Chandler rbchan@uga.edu
You cannot use obsCovs, but you can use yearlySiteCovs (a confusing name since this model isn't for multi-year data. It's just a hold-over from the colext methods of formatting data upon which it is based.)
# Simulate some line-transect data set.seed(36837) R <- 50 # number of transects T <- 5 # number of replicates strip.width <- 50 transect.length <- 100 breaks <- seq(0, 50, by=10) lambda <- 5 # Abundance phi <- 0.6 # Availability sigma <- 30 # Half-normal shape parameter J <- length(breaks)-1 y <- array(0, c(R, J, T)) for(i in 1:R) { M <- rpois(1, lambda) # Individuals within the 1-ha strip for(t in 1:T) { # Distances from point d <- runif(M, 0, strip.width) # Detection process if(length(d)) { cp <- phi*exp(-d^2 / (2 * sigma^2)) # half-normal w/ g(0)<1 d <- d[rbinom(length(d), 1, cp) == 1] y[i,,t] <- table(cut(d, breaks, include.lowest=TRUE)) } } } y <- matrix(y, nrow=R) # convert array to matrix # Organize data umf <- unmarkedFrameGDS(y = y, survey="line", unitsIn="m", dist.breaks=breaks, tlength=rep(transect.length, R), numPrimary=T) summary(umf)#> unmarkedFrame Object #> #> 50 sites #> Maximum number of observations per site: 25 #> Mean number of observations per site: 25 #> Number of primary survey periods: 5 #> Number of secondary survey periods: 1 #> Sites with at least one detection: 50 #> #> Tabulation of y observations: #> 0 1 2 3 4 #> 838 334 65 12 1#> #> Call: #> gdistsamp(lambdaformula = ~1, phiformula = ~1, pformula = ~1, #> data = umf, output = "density", K = 50) #> #> Abundance (log-scale): #> Estimate SE z P(>|z|) #> 1.53 0.125 12.2 1.89e-34 #> #> Availability (logit-scale): #> Estimate SE z P(>|z|) #> 0.519 0.318 1.63 0.103 #> #> Detection (log-scale): #> Estimate SE z P(>|z|) #> 3.44 0.0706 48.8 0 #> #> AIC: 1812.616 #> Number of sites: 50 #> optim convergence code: 0 #> optim iterations: 38 #> Bootstrap iterations: 0 #>#> Backtransformed linear combination(s) of Abundance estimate(s) #> #> Estimate SE LinComb (Intercept) #> 4.6 0.574 1.53 1 #> #> Transformation: exp#> Backtransformed linear combination(s) of Availability estimate(s) #> #> Estimate SE LinComb (Intercept) #> 0.627 0.0745 0.519 1 #> #> Transformation: logistic#> Backtransformed linear combination(s) of Detection estimate(s) #> #> Estimate SE LinComb (Intercept) #> 31.3 2.21 3.44 1 #> #> Transformation: exp