The `closedpCI.t`

and `closedpCI.0`

functions fit a loglinear model specified by the user
and compute a confidence interval for the abundance estimation. For a normal heterogeneous model,
a log-transformed confidence interval (Chao 1987) is produced.
For any other model, the multinomial profile likelihood confidence interval (Cormack 1992) is produced.

The model is identified with the argument `m`

or `mX`

.
For heterogeneous models, the form of the heterogeneity is specified with the arguments
`h`

and `h.control`

. If `h`

is given with `mX`

, heterogeneity is added in `mX`

.

These functions extend `closedp.t`

and `closedp.0`

as
they broaden the range of models one can fit and they compute confidence intervals.
Unlike the `closedp`

functions, it fits only one model at a time.

```
closedpCI.t(X, dfreq=FALSE, m=c("M0","Mt","Mh","Mth"), mX=NULL,
h=NULL, h.control=list(), mname=NULL, alpha=0.05,
fmaxSupCL=3, …)
```closedpCI.0(X, dfreq=FALSE, dtype=c("hist","nbcap"), t=NULL, t0=NULL,
m=c("M0","Mh"), mX=NULL, h=NULL, h.control=list(),
mname=NULL, alpha=0.05, fmaxSupCL=3, …)

plotCI(x.closedpCI, main="Profile Likelihood Confidence Interval", …)
# S3 method for closedpCI
print(x, …)

# S3 method for closedpCI
boxplot(x, main="Boxplots of Pearson Residuals", …)

# S3 method for closedpCI
plot(x, main="Scatterplot of Pearson Residuals", …)

X

The matrix of the observed capture histories (see `Rcapture-package`

for a
description of the accepted formats).

dfreq

A logical. By default `FALSE`

, which means that `X`

has one row per unit.
If `TRUE`

, it indicates that the matrix `X`

contains frequencies in its last column.

dtype

A characters string, either `"hist"`

or `"nbcap"`

, to specify the type of data.
`"hist"`

, the default, means that `X`

contains complete observed capture histories.
`"nbcap"`

means that `X`

contains numbers of captures (see `Rcapture-package`

for details on data formats).

t

Requested only if `dtype="nbcap"`

. A numeric specifying the total number of
capture occasions in the experiment. For `closedpCI.0`

, the value `t=Inf`

is accepted. It indicates that captures occur in continuous time
(see `Rcapture-package`

).

t0

A numeric. Models are fitted considering only the frequencies of units captured 1 to `t0`

times.
By default, if `t`

is not equal to `Inf`

, `t0=t`

. When `t=Inf`

, the default value
of `t0`

is the larger number of captures observed.

m

A character string indicating the model to fit. For `closedpCI.0`

it can be either "M0"=M0 model
or "Mh"=Mh model. For `closedpCI.t`

it can also be "Mt"=Mt model or "Mth"=Mth model.

mX

The design matrix of the loglinear model. By default, the design matrix is built based
on the `m`

argument. If a `mX`

argument is given, it must be a matrix (or
an object that can be coerced to a matrix by `as.matrix`

). The order of the capture
histories in the `mX`

matrix must be as defined in the `histpos.t`

or
`histpos.0`

function.
For the `closedpCI.t`

function only, `mX`

can also be an object of class "`formula`

" or a hierarchical loglinear model name. The details of a formula specification are given under **Details** and the details of a hierarchical loglinear model name specification are given in the `closedpMS.t`

documentation.
An intercept is always added to the model. Therefore, `mX`

must not contain a column of ones if it is a matrix. Also, if the argument `h`

is not `NULL`

, one or more columns for heterogeneity are added to the design matrix.

h

A character string ("LB", "Chao", "Poisson", "Darroch", "Gamma" or "Normal") or a
numerical `R`

function specifying the form of the column(s) for heterogeneity
in the design matrix.
"LB" and "Chao" represents Chao's lower bound model (the default),
"Poisson" represents the function \(f(k)=theta^k-1\),
"Darroch" represents the function \(f(k)=k^2/2\),
and "Gamma" represents the function \(f(k)=-\log(theta + k) + \log(theta)\)
where \(k\) is the number of captures and \(theta\) is a parameter specified with
the `h.control`

element `theta`

.
"Normal" represents heterogeneity
modeled with a normal mixing distribution, as presented in Rivest (2011).
If an `R`

function is given, it must be the implementation of any convex
mathematical function \(f(k)\) (it must have only one argument).

h.control

A list of elements to control the heterogeneous part of the model, if any. For a Poisson or Gamma heterogeneous model:

`theta`

The value of the parameter \(theta\) in \(f(k)=theta^k-1\) for the Poisson model (the default value is 2) and in \(f(k)=-\log(theta + k) + \log(theta)\) for the Gamma model (the default value is 3.5).

For a Chao's lower bound heterogeneous model:

`neg`

If this option is set to

`TRUE`

(the default), negative eta parameters in Chao's lower bound models are set to zero.

For a Normal heterogeneous model:

mname

A character string specifying the name of the customized model. By default, it is derived from the arguments specifying the model.

alpha

A confidence interval with confidence level 1-`alpha`

is constructed. The value of alpha
must be between 0 and 1; the default is 0.05.

fmaxSupCL

A numeric (by default 3). The upper end point of the interval to be searched by
`uniroot`

to find the upper bound of the multinomial profile likelihood
confidence interval (Cormack 1992) is defined by
`fmaxSupCL`

\(\times \hat{N}\). In
this product, `fmaxSupCL`

is a factor multiplying \(\hat{N}\), the population
size estimated by Poisson regression. If the upper bound obtained is equal
the upper end point `fmaxSupCL`

\(\times \hat{N}\), the
element `SupCL`

of the output value `CI`

will begin with the symbol
\(>\). In that case, one could try to increase the value of `fmaxSupCL`

to find the upper bound of the multinomial profile likelihood confidence interval
instead of a minimum value for this upper bound.

…

Further arguments to be passed to `glm`

, `optim`

,
`print.default`

, `plot.default`

or `boxplot.default`

.

x.closedpCI

An object, produced by a `closedpCI`

function, to produce a plot
of the multinomial profile likelihood for \(N\).

main

A main title for the plot

x

An object, produced by a `closedpCI`

function, to print or to plot.

The number of captured units

The total number of capture occasions in the data matrix `X`

.
When captures occur in continuous time (input argument `t=Inf`

),
this output value is the larger number of captures observed.

For `closedpCI.0`

only: the value of the argument `t0`

used in the computations.

A table containing, for the fitted model:

`abundance`

: the estimated population size,

`stderr`

: the standard error of the estimated population size,

`deviance`

: the model's deviance,

`df`

: the number of degrees of freedom,

`AIC`

: the Akaike's information criterion,

`BIC`

: the bayesian information criterion,

`infoFit`

: a numerical code giving information about error or warnings encountered when fitting the model (see

`Rcapture-package`

for details).

Not produced for normal heterogeneous models (`h="Normal"`

): The asymptotic bias of the estimated population size.

The 'glm' object obtained from fitting the model except for
normal heterogeneous models (`h="Normal"`

). These models are not fitted
with `glm`

. For them, `fit`

is a list with the following elements:

`parameters`

the matrix of parameters (loglinear coefficients + sigma parameter) estimates with their standard errors,

`varcov`

the estimated variance-covariance matrix of the estimated parameters,

`y`

the y vector used to fit the model,

`fitted.values`

the model fitted values,

`initparam`

the initial values for the parameters (loglinear coefficients + sigma parameter) used by

`optim`

,`optim`

the output produced by

`optim`

.

A vector of character strings. If the `glm`

or `optim`

(for normal heterogeneous models) function generates
one or more warnings when fitting the model, a copy of these warnings are
stored in `fit.warn`

. `NULL`

if no warnings occured.

For Chao's lower bound models only: the position of the eta parameters set to zero in the loglinear parameter vector, if any.

Not produced for normal heterogeneous models (`h="Normal"`

): A table containing
the abundance estimation and its multinomial profile likelihood confidence interval.
The last column of the table, named `infoCI`

, contains a numerical code giving information about error or warnings encountered when calculating the confidence interval. Here is a description of the meaning of this numerical code:

- 0
no error or warning occured;

- -1
an error occured while calculating the multinomial abundance estimation;

- 4
a warning occured while calculating the multinomial abundance estimation;

- 5
a warning occured while calculating the lower bound of the multinomial profile likelihood confidence interval;

- 6
a warning occured while calculating the upper bound of the multinomial profile likelihood confidence interval.

Not produced for normal heterogeneous models (`h="Normal"`

): If
an error occured while calculating the multinomial abundance, a copy of the error message is
stored in `CI.err`

. `NULL`

if no error occured.

Not produced for normal heterogeneous models (`h="Normal"`

): If
one or more warnings occur while calculating the multinomial abundance estimation or the profile likelihood confidence interval, a copy of these warnings are
stored in `CI.warn`

. `NULL`

if no warnings occured.

1-the confidence level of the interval.

Not produced for normal heterogeneous models (`h="Normal"`

):
The x-coordinates for `plot.closedpCI.t`

.

Not produced for normal heterogeneous models (`h="Normal"`

):
The y-coordinates for `plot.closedpCI.t`

.

The `closedpCI.t`

function fits models using the frequencies of the observable capture histories (vector of size \(2^t-1\)), whereas `closedpCI.0`

uses the number of units capture i times, for \(i=1,\ldots,t\) (vector of size \(t\)). Thus, `closedpCI.0`

can be used with data sets larger than those for `closedpCI.t`

, but it cannot fit models with a temporal effect. See `Rcapture-package`

for more details about the distinction between `.t`

and `.0`

functions.

These functions do not work for closed population models featuring a behavioral effect, such as Mb and Mbh. The abundance estimation is calculated as the number of captured units plus the exponential of the Poisson regression intercept. However, models with a behavioral effect can by fitted with `closedp.t`

(Mb and Mbh), `closedp.Mtb`

and `closedp.bc`

.

CHAO'S LOWER BOUND MODELS

Chao's lower bound models estimate a lower bound for the abundance. Rivest (2011)
presents a generalized loglinear model underlying this estimator. To test whether a
certain model for heterogeneity is adequate,
one can conduct a likelihood ratio test by subtracting the deviance of a Chao's lower
bound model to the deviance of the heterogeneous model under study. The two models should
have the same `mX`

argument.
Under the null hypothesis of equivalence between the two models, the difference of deviances
follows a chi-square distribution with degrees of freedom equal to the difference between
the models' degrees of freedom.

A Chao's lower bound model contains \(t-2\) parameters, called
eta parameters, for the heterogeneity. These parameters should theoretically be greater
or equal to zero (see Rivest and Baillargeon (2007)). When the element `neg`

of
the argument `h.control`

is set to `TRUE`

(the default), negative eta parameters are
set to zero (to do so, columns are removed from the design matrix of the model). Degrees
of freedom of Chao's lower bound model increase when eta parameters are set to zero.

ARGUMENT `mX`

: FORMULA SPECIFICATION
For the `closedpCI.t`

function, `mX`

can be an object of class "`formula`

". The only accepted variables
in this formula are `c1`

to `ct`

. The variable `ci`

represents
a capture indicator (1 for a capture, 0 otherwise) for the \(i\)th capture occasions.
Also, the formula must not contain a response variable since
it is only used to construct the design matrix of the model.
For example, if `t=3`

, the Mt model is fitted if
`mX = ~ .`

or `mX = ~ c1 + c2 + c3`

. The symbol `.`

in this formula
is a shortcut for `c1`

+ `c2`

+ ... + `ct`

. Formula `mX`

arguments
facilitate the addition of interactions between capture occasions in the model. For
example, if `t=3`

, the Mt model with an interaction between the first and
the second capture occasion is fitted if `mX = ~ . + c1:c2`

.
See `formula`

for more details of allowed formulae.

PLOT METHODS AND FUNCTIONS

The `boxplot.closedpCI`

function produces a boxplot of the Pearson residuals of the customized model.

The `plot.closedpCI`

function traces the scatterplot of the Pearson residuals in terms of \(f_i\)
(number of units captured i times) for the customized model.

The `plotCI`

function produces a plot of the multinomial profile likelihood for \(N\).
The value of N maximizing the profile likelihood and the bounds of the confidence interval are identified.

Baillargeon, S. and Rivest, L.P. (2007) Rcapture: Loglinear models for capture-recapture in R. *Journal of Statistical Software*, **19**(5), http://www.jstatsoft.org/v19/i05.

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

Cormack, R. M. (1992) Interval estimation for mark-recapture studies of closed populations. *Biometrics*, **48**, 567--576.

Rivest, L.P. (2011) A lower bound model for multiple record systems estimation with heterogeneous catchability.
*The International Journal of Biostatistics*, **7**(1), Article 23.

Rivest, L.P. and Baillargeon, S. (2007) Applications and extensions of Chao's moment estimator for the size of a closed population. *Biometrics*, **63**(4), 999--1006.

```
# NOT RUN {
# hare data set
CI <- closedpCI.t(hare, m = "Mth", h = "Poisson", h.control = list(theta = 2))
CI
plotCI(CI)
# HIV data set
mat <- histpos.t(4)
mX2 <- cbind(mat, mat[, 1] * mat[ ,2])
closedpCI.t(HIV, dfreq = TRUE, mX = mX2, mname = "Mt interaction 1,2")
# which can be obtained more conveniently with
closedpCI.t(HIV, dfreq = TRUE, mX = ~ . + c1:c2, mname = "Mt interaction 1,2")
# BBS2001 data set
CI0 <- closedpCI.0(BBS2001, dfreq = TRUE, dtype = "nbcap", t = 50, t0 = 20,
m = "Mh", h = "Gamma", h.control = list(theta = 3.5))
CI0
plot(CI0)
plotCI(CI0)
### As an alternative to a gamma model, one can fit a negative Poisson model.
### It is appropriate in experiments where very small capture probabilities
### are likely. It can lead to very large estimators of abundance.
# Third primary period of mvole data set
period3 <- mvole[, 11:15]
psi <- function(x) { 0.5^x - 1 }
closedpCI.t(period3, m = "Mh", h = psi)
### Example of normal heterogeneous models
### diabetes data of Bruno et al. (1994)
histpos <- histpos.t(4)
diabetes <- cbind(histpos, c(58,157,18,104,46,650,12,709,14,20,7,74,8,182,10))
# chosen interaction set I in Rivest (2011)
closedpCI.t(X = diabetes, dfreq = TRUE, mX = ~ . + c1:c3 + c2:c4 + c3:c4,
h = "Normal", mname = "Mth normal with I")
### Example of captures in continuous time
# Illegal immigrants data set
closedpCI.0(ill, dtype = "nbcap", dfreq = TRUE, t = Inf, m = "Mh", h = "LB")
# }
```

Run the code above in your browser using DataCamp Workspace