Priors for hierarchial Bayesian binomial model

Question

I intend to use a hierarchical binomial model which would look like the following in BUGS language :

model{
for(i in 1:I){
   y[i] ~ dbin(theta[k[i]],n[k[i]])
}
for(j in 1:J){
   theta[j] ~ dprior(alpha, beta)
}
alpha ~ dhyperprior.alpha
beta ~ dhyperprior.beta
}

I wonder about the choice of the prior on $\theta$ and the hyperpriors. Do you know how to do this choice in order that:

1) There's no need to use MCMC (classical simulations suffice)

and/or

2) This is a good "noninformative" choice, in the sense of the frequentist-matching property (the frequentist coverage of the posterior credible intervals is close to the credibility level)

Answer 1

In answer to (1):

Set

$\theta_i \sim \text{Beta}(\alpha, \beta)$

Then the marginal data distribution, $\int_0^1\Pr(y_i|\theta_i)p(\theta_i|\alpha, \beta)d\theta_i$, is a standard distribution:

$y_i \sim \text{Beta-binomial}(\alpha, \beta)$

(Wikipedia article: Beta-binomial distribution; see especially the section on "Further Bayesian considerations".) Since this results in a bivariate posterior distribution, you can easily evaluate it on a finely-spaced grid; thus MCMC can be avoided for any choice of hyperpriors.

More specifically, you can find the region of high posterior density by numerical optimization, using the ML estimate as the initial guess. After finding the region containing all but a negligible amount of posterior probability mass and evaluating the posterior density on the fine grid, you can generate random $(\alpha, \beta)$ pairs by picking one of the points on the grid with probability proportional to the evaluated density and then adding uniform jitter centered at the selected point with width equal to the grid spacing.

If you're interested in inferring the $\theta_i$ variables, you can sample them by first sampling a set of $(\alpha, \beta)$ values from their posterior and then sampling the $\theta_i$ from their conditional posteriors,

$\theta_i|y_i, \alpha, \beta \sim \text{Beta}(\alpha + y_i, \beta + n_i - y_i)$

Answer 2

For my binomial models, I use normal priors - on the logit scale:

for (...) {
    y[i, j, ...] ~ dbin(p, M[i])
    ....
}

# priors
p <- 1/(1+exp(-logit_p))
logit_p ~ dnorm(0, tau)
tau <- 1/(4 * 4)

But note that I'm not an expert, on the contrary - I just considered this prior to be better than usual uniform (p ~ dunif(0, 1)). I plan to try better prior mentioned at http://psiexp.ss.uci.edu/research/programs_data/jags/Example3.html#1

model {
   # Prior on Rate
   theta ~ dbeta(1,1)
   # Observed Counts
   k ~ dbin(theta,n)
}

Which quite corresponds to what Cyan is proposing in his answer.