unmarked (version 0.9-1)

pcountOpen: Fit the open N-mixture model of Dail and Madsen

Description

Fit the model of Dail and Madsen, which is a generalized form of the Royle (2004) N-mixture model for open populations.

Usage

pcountOpen(lambdaformula, gammaformula, omegaformula, pformula, data, 
	mixture = c("P", "NB"), K, dynamics=c("constant", "autoreg", "notrend"),
    fix=c("none", "gamma", "omega"), starts, method = "BFGS", se = TRUE, ...)

Arguments

lambdaformula

Right-hand sided formula for initial abundance

gammaformula

Right-hand sided formula for recruitment rate

omegaformula

Right-hand sided formula for apparent survival probability

pformula

Right-hand sided formula for detection probability

data

An object of class unmarkedFramePCO. See details

mixture

character specifying mixture: either "P" or "NB" for the Poisson and negative binomial distributions.

K

Integer defining upper bound of discrete integration. This should be higher than the maximum observed count and high enough that it does not affect the parameter estimates. However, the higher the value the slower the compuatation.

dynamics

Character string describing the type of population dynamics. "constant" indicates that there is no relationship between omega and gamma. "autoreg" is an auto-regressive model in which recruitment is modeled as gamma*N[i,t-1]. "notrend" model gamma as lamb

fix

If "omega", omega is fixed at 1. If "gamma", gamma is fixed at 0.

starts

vector of starting values

method

Optimization method used by optim.

se

logical specifying whether or not to compute standard errors.

...

additional arguments to be passed to optim.

Value

• An object of class unmarkedFitPCO.

Warning

This function can be extremely slow, especially if there are covariates of gamma or omega. Consider testing the timing on a small subset of the data, perhaps with se=FALSE. Finding the lowest value of K that does not affect estimates will also help with speed.

Details

This model generalizes the Royle (2004) N-mixture model by relaxing the closure assumption. The model includes two additional parameters: gamma, the recruitment rate (births and immigrations), and omega, the apparent survival rate (deaths and emigrations). Estimates of population size at each time period can be derived from these parameters, and thus so can trend estimates. The latent initial abundance distribution, $f(N | \mathbf{\theta})$ can be set as either a Poisson or a negative binomial random variable, depending on the setting of the mixture argument. mixture = "P" or mixture = "NB" select the Poisson or negative binomial distribution respectively. The mean of $N_i$ is $\lambda_i$. If $N_i \sim NB$, then an additional parameter, $\alpha$, describes dispersion (lower $\alpha$ implies higher variance). The latent abundance state following the initial sampling period arises from a Markovian process in which survivors are modeled as $S_it \sim Binomial(N_it-1, omega_it)$, and recruits follow $G_it \sim Poisson(gamma_it)$. The detection process is modeled as binomial: $y_{it} \sim Binomial(N_it, p_it)$. Covariates of $\lambda_i$ use the log link and covariates of $p_it$ use the logit link.

References

Royle, J. A. (2004) N-Mixture Models for Estimating Population Size from Spatially Replicated Counts. Biometrics 60, pp. 108--105. Dail, D. and L. Madsen (In press) Models for Estimating Abundance from Repeated Counts of an Open Metapopulation. Biometrics.

See Also

pcount, unmarkedFramePCO

Examples

## Simulation 
## No covariates, constant time intervals between primary periods, and 
## no secondary sampling periods

set.seed(3)
M <- 50
T <- 5
lambda <- 1
gamma <- 0.5
omega <- 0.8
p <- 0.7
y <- N <- matrix(NA, M, T)
S <- G <- matrix(NA, M, T-1)
N[,1] <- rpois(M, lambda)
for(t in 1:(T-1)) {
	S[,t] <- rbinom(M, N[,t], omega)
	G[,t] <- rpois(M, gamma)
	N[,t+1] <- S[,t] + G[,t]
	}
y[] <- rbinom(M*T, N, p)

                            
# Prepare data                               
umf <- unmarkedFramePCO(y = y, numPrimary=T)

summary(umf)


# Fit model and backtransform
(m1 <- pcountOpen(~1, ~1, ~1, ~1, umf, K=10))

(lam <- coef(backTransform(m1, "lambda"))) # or
lam <- exp(coef(m1, type="lambda"))
gam <- exp(coef(m1, type="gamma"))
om <- plogis(coef(m1, type="omega"))
p <- plogis(coef(m1, type="det"))


# Calculate population size at each time period
N.hat <- matrix(NA, M, T)
N.hat[,1] <- lam
for(t in 2:T) {
    N.hat[,t] <- om*N.hat[,t-1] + gam
    }

rbind(N=colSums(N), N.hat=colSums(N.hat))