dga (version 2.0.1)

bma.cr: Bayesiam Model Averaging for Capture-Recapture

Description

This function averages over all decomposable graphical models for p lists.

Usage

bma.cr(
  Y,
  Nmissing,
  delta,
  graphs,
  logprior = NULL,
  log.prior.model.weights = NULL
)

Arguments

Y

a 2^p array of list intersection counts. See details.

Nmissing

A vector of all possible values for the number of individuals that appear on no list.

delta

The hyper-parameter for the hyper-Dirichlet prior distribution on list intersection probabilities. A smaller value indicates fewer prior observations per cell. A suggested default is 2^-p

graphs

A pre-computed list of all decomposable graphical models for p lists. These should be loaded using data(graphsp); see example. Currently, this package includes a list of graphs for three, four, or five lists.

logprior

The log of the prior probability of each value in Nmissing. If left blank, this will default to the -log(Nmissing).

log.prior.model.weights

Prior weights on the graphs. This should be a vector of length length(graphs).

Value

This function returns a matrix of weights, where rows correspond to models and columns correspond to values of Nmissing. Thus, the ijth entry of the matrix is the posterior probability of the ith model and the jth entry of Nmissing. Row sums return posterior probabilities by model.Column sums return posterior probabilities by value of Nmissing.

Details

This is the main function in this package. It performs capture-recapture (or multiple systems estimation) using Bayesian model averaging as outlined in Madigan and York (1997).

Y can be created by the array() command from a vector that is ordered lexigraphically by the cell names, e.g., c(x000, x001, x010, x011, x100, x101, x110, x111).

References

Madigan, David, and Jeremy C. York. "Bayesian methods for estimation of the size of a closed population." Biometrika 84.1 (1997): 19-31.

Examples

# NOT RUN {
#### 5 list example from M & Y ##########
delta <- .5
Y <- c(0, 27, 37, 19, 4, 4, 1, 1, 97, 22, 37, 25, 2, 1, 3, 5,
       83, 36, 34, 18, 3, 5, 0, 2, 30, 5, 23, 8, 0, 3, 0, 2)
Y <- array(Y, dim = c(2, 2, 2, 2, 2))
Nmissing <- 1:300
N <- Nmissing + sum(Y)
data(graphs5)
weights <- bma.cr(Y, Nmissing, delta, graphs5)
plotPosteriorN(weights, N)

##### 3 list example from M & Y #######
Y <- c(0, 60, 49, 4, 247, 112, 142, 12)
Y <- array(Y, dim = c(2, 2, 2))

delta <- 1
a <- 13.14
b <- 55.17


Nmissing <- 1:300
N <- Nmissing + sum(Y)

logprior <- N * log(b) - (N + a) * log(1 + b) + lgamma(N + a) - lgamma(N + 1) - lgamma(a)

data(graphs3)
weights <- bma.cr(Y, Nmissing, delta, graphs3, logprior)
plotPosteriorN(weights, N)
# }