bma.cr {dga} | R Documentation |
This function averages over all decomposable graphical models for p lists.
bma.cr(Y, Nmissing, delta, graphs, logprior = NULL, log.prior.model.weights, normalize)
Y |
a |
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 |
graphs |
A pre-computed list of all decomposable graphical models for |
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). |
normalize |
Default is TRUE. This is mostly for debugging. Should we normalize everything to output posterior probabilities? |
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).
This function returns a matrix of weights, where rows correspond to models and columns correspond to values of Nmissing. Thus, the ij
th entry of the matrix is the posterior probability of the i
th model and the j
th entry of Nmissing. Row sums return posterior probabilities by model.Column sums return posterior probabilities by value of Nmissing.
This function is pretty robust relative to the more common log-linear model approach to capture-recapture, i.e. it will not issue a numerical warning even if there are no overlaps among the lists. Care should be taken that there is adequate list overlap.
James Johndrow james.johndrow@gmail.com and Kristian Lum (kl@hrdag.org)
Madigan, David, and Jeremy C. York. "Bayesian methods for estimation of the size of a closed population." Biometrika 84.1 (1997): 19-31.
#### 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)