Learn R Programming

MCMCpack (version 0.4-5)

MCMCdynamicEI: Markov chain Monte Carlo for Quinn's Dynamic Ecological Inference Model

Description

MCMCdynamicEI is used to fit Quinn's dynamic ecological inference model for partially observed 2 x 2 contingency tables.

Usage

MCMCdynamicEI(r0, r1, c0, c1, burnin=5000, mcmc=200000, thin=200,
              tune=2.65316, verbose=FALSE, seed=0, W=0, nu0=1,
              delta0=0.03, nu1=1, delta1=0.03, ...)

Arguments

r0
$(ntables \times 1)$ vector of row sums from row 0.
r1
$(ntables \times 1)$ vector of row sums from row 1.
c0
$(ntables \times 1)$ vector of column sums from column 0.
c1
$(ntables \times 1)$ vector of column sums from column 1.
burnin
The number of burn-in scans for the sampler.
mcmc
The number of mcmc scans to be saved.
thin
The thinning interval used in the simulation. The number of mcmc iterations must be divisible by this value.
tune
Tuning parameter for the Metropolis-Hasting sampling.
verbose
A switch which determines whether or not the progress of the sampler is printed to the screen. Information is printed if TRUE.
seed
The seed for the random number generator. The code uses the Mersenne Twister, which requires an integer as an input. If nothing is provided, the Scythe default seed is used.
W
Weight (not precision) matrix structuring the temporal dependence among elements of $\theta_{0}$ and $\theta_{1}$. The default value of 0 will construct a weight matrix that corresponds to random walk priors for $\theta_{0}$ and
nu0
Shape parameter for the inverse-gamma prior on the $\sigma^2_0$ parameter.
delta0
Scale parameter for the inverse-gamma prior on the $\sigma^2_0$ parameter.
nu1
Shape parameter for the inverse-gamma prior on the $\sigma^2_1$ parameter.
delta1
Scale parameter for the inverse-gamma prior on the $\sigma^2_1$ parameter.
...
further arguments to be passed

Value

  • An mcmc object that contains the posterior density sample. This object can be summarized by functions provided by the coda package.

Details

Consider the following partially observed 2 by 2 contingency table for unit $t$ where $t=1,\ldots,ntables$: llll{ | $Y=0$ | $Y=1$ | - - - - - - - - - - - - - - - - - - - - $X=0$ | $Y_{0t}$ | | $r_{0t}$ - - - - - - - - - - - - - - - - - - - - $X=1$ | $Y_{1t}$ | | $r_{1t}$ - - - - - - - - - - - - - - - - - - - - | $c_{0t}$ | $c_{1t}$ | $N_t$ } Where $r_{0t}$, $r_{1t}$, $c_{0t}$, $c_{1t}$, and $N_t$ are non-negative integers that are observed. The interior cell entries are not observed. It is assumed that $Y_{0t}|r_{0t} \sim \mathcal{B}inomial(r_{0t}, p_{0t})$ and $Y_{1t}|r_{1t} \sim \mathcal{B}inomial(r_{1t}, p_{1t})$. Let $\theta_{0t} = log(p_{0t}/(1-p_{0t}))$, and $\theta_{1t} = log(p_{1t}/(1-p_{1t}))$.

The following prior distributions are assumed: $$p(\theta_0|\sigma^2_0) \propto \sigma_0^{-ntables} \exp \left(-\frac{1}{2\sigma^2_0} \theta'_{0} P \theta_{0}\right)$$ and $$p(\theta_1|\sigma^2_1) \propto \sigma_1^{-ntables} \exp \left(-\frac{1}{2\sigma^2_1} \theta'_{1} P \theta_{1}\right)$$ where $P_{ts}$ = $-W_{ts}$ for $t$ not equal to $s$ and $P_{tt}$ = $\sum_{s \ne t}W_{ts}$. The $\theta_{0t}$ is assumed to be a priori independent of $\theta_{1t}$ for all t. In addition, the following hyperpriors are assumed: $\sigma^2_0 \sim \mathcal{IG}(\nu_0/2, \delta_0/2)$, and $\sigma^2_1 \sim \mathcal{IG}(\nu_1/2, \delta_1/2)$.

Inference centers on $p_0$, $p_1$, $\sigma^2_0$, and $\sigma^2_1$. The Metropolis-Hastings algorithm is used to sample from the posterior density.

References

Jonathan Wakefield. 2001. ``Ecological Inference for 2 x 2 Tables." Center for Statistics and the Social Sciences Working Paper no. 12. University of Washington.

Kevin M. Quinn. 2002. ``Ecological Inference in the Presence of Temporal Dependence." Paper prepared for Ecological Inference Conference, Harvard University, June 17-18, 2002.

Andrew D. Martin, Kevin M. Quinn, and Daniel Pemstein. 2003. Scythe Statistical Library 0.4. http://scythe.wustl.edu. Martyn Plummer, Nicky Best, Kate Cowles, and Karen Vines. 2002. Output Analysis and Diagnostics for MCMC (CODA). http://www-fis.iarc.fr/coda/.

See Also

MCMCbaselineEI, MCMChierEI, plot.mcmc,summary.mcmc

Examples

Run this code
r0 <- rpois(10, 300)
   r1 <- rpois(10, 200)
   c0 <- 100 + 1:10*7 + rpois(10, 30)
   c1 <- (r0+r1) - c0
   posterior <- MCMCdynamicEI(r0, r1, c0, c1,
      burnin=10000, mcmc=5000000, thin=1000)
   plot(posterior)
   summary(posterior)

Run the code above in your browser using DataLab