# Problem

I’m currently working on a problem where I have count data for $n$ items in the following form:

    Item   Count   Total
--------------------
1          1      30
2         10     340
3          0       0
4          1      10
...


I want to calculate the posterior of the success probability $p_i$ for each item. Additionally, I suspect that all the items have a similar or even the same success probability but I want the model to determine that.

My expectations for the resulting model are:

• Estimate the underlying global parameter of all the items and
• Shrink the success probabilities $p_i$ towards this global parameter in order to deal with items which have few observations.
• Items with observations indicating a different success probability should be able to “override” the global success probability.

# Stan code

This is the Stan code I implemented so far:

    data {
int<lower=0> J; // amount of items
int count[J];
int total[J];
}
parameters {
real<lower=0,upper=1> mu;
real<lower=0.001> ss;
real<lower=0,upper=1> prob[J];
}
model {
mu ~ uniform(0, 1);
ss ~ gamma(0.01, 0.01);
prob ~ beta(mu*ss, (1-mu)*ss);
count ~ binomial(total, prob);
}


Here I achieve the shrinkage by assuming that the success probabilities of the items are sampled from a beta distribution (parametrized using mean and sample size) where the mean is the global success probability.

The problem occurs when specifying the prior of the sample size. If I choose a very vague prior, Stan requires a lot of iterations until the effective sample size is large enough. Also I’m not sure if this is the right approach since it is not clear if the amount of shrinkage this model provides is always appropriate.

I’m interested if there is an elegant way to implement my expected behavior or if there is some flaw in my overall model.

• I've been struggling with a related question, regarding priors on hyperparameters for negative binomials. I don't want to confuse you, but it might be helpful to look at a related model. stats.stackexchange.com/questions/94923/… May 6, 2014 at 17:05

Try this: $p_i\stackrel{iid}{\sim} Be(\alpha,\beta)$ and $p(\alpha,\beta)\propto (\alpha+\beta)^{-5/2}$.

I believe the issue you are running into is that $ss\sim Ga(\gamma,\gamma)$ with $\gamma\to 0$ results in the improper $p(ss)\propto 1/ss$ prior. This prior, together with your uniform prior on mu, results in an improper posterior. Despite the fact that you are using a proper prior, it is close enough to this improper prior to cause issues.

You might want to take a look at page 110 of Bayesian Data Analysis (3rd ed) for a discussion of priors for this model as well as the prior suggested above.

• I came here to say that. Have an upvote.
– Sycorax
May 5, 2014 at 21:36
• Is there some clever way of constructing/approximating this prior using standard probability distributions? I'm asking because in programs like JAGS or OpenBugs it is much easier to use the built in distributions rather than using a "custom" prior. For example, instead of putting a uniform on $log(\theta)$ one can approximate this with theta ~ dgamma(0.001, 0.001) Jun 19, 2014 at 9:00
• In JAGS and BUGS, you can use the zeros-trick. See users.aims.ac.za/~mackay/BUGS/Manuals/Tricks.html. Jun 23, 2014 at 14:31

This Blog entry here http://lingpipe-blog.com/2009/09/23/bayesian-estimators-for-the-beta-binomial-model-of-batting-ability/ shows another possibility for modelling the prior on the parameters of the Beta population distribution, namely

$\frac{\alpha}{\alpha + \beta}$ ~ Uniform(0,1) = Beta(1,1)

$(a+b) \sim$ Pareto(1.5, 1)

Here is another resource regarding your question (slide 5): http://www.stat.cmu.edu/~brian/724/week06/lec15-mcmc.pdf.

• Indeed, choosing the pareto distribution as hyperprior for the sample size works and has the advantage that it’s easy to implement in Stan. I implemented both approaches and as expected the pareto prior provides slightly more shrinkage. May 6, 2014 at 14:52