April 17, 2013

The simplest possible Bayesian model

In chapter 5, Kruschke goes through the process of using Bayes' rule for updating belief in the 'lopsidedness' of a coin. Instead of using R I decided to try implementing the model in Excel.

The Model

Coin tosses are IID with a $\mathrm{Bernoulli}(\theta)$ distribution, thus the probability of observing $N$ of which $z$ are heads is
$$p(z|\theta, N) = \theta^z (1-\theta)^{(N-z)}$$
Using Bayes' rule we get
$$p(\theta|z,N) = \frac{\theta^z (1-\theta)^{N-z} p(\theta)}{p(z|N)}$$
Note that if $p(\theta)$ is already of the form
$$f(\theta;a,b) = \frac{\theta^{a-1} (1-\theta)^{b-1}}{B(a,b)}$$
where the denominator is a normalizing constant then $p(\theta|N) = f(\theta;a+z,b+N-z)$.

Excel Implementation

This is very straight-forward using the BETADIST function, which is cumulative, so we divide the unit interval into equally spaced subintervals, calculate the probability of $\theta$ being in the interval and update on the basis of the new observations.

Here's the result:


Discrete model

Using a discrete distribution for $\theta$ rather than a continuous parametric form we can allow for much more freedom in the choice of initial prior distribution. The update rule becomes
$$p(\theta|z_N,N)=\frac{p(z_N|z_{N-1},N-1,\theta)p(\theta|z_{N-1},N-1)}{\sum_{i=1}^n p(z_N|z_{N-1},N-1,\theta_i) p(\theta_i|z_{N-1},N-1)}$$
in other words
$$p_{N,j} = \frac{\tilde p_{N,j}}{\sum_{i=1}^n \tilde p_{N,i}}$$
where
$$\tilde p_{N,k} = (\theta_k y_N + (1-\theta_k)(1-y_N)) p_{N-1,j}$$
and $y_N = z_N - z_{N-1}$ is the result of the $N$th coin toss.

As expected the result is very similar to the continuous model.

In the discrete model we can assign any initial prior. Here is an example with an initial prior that assigns non-zero probability to only 3 values of $\theta$ (22.5%, 52.5% and 77.5%) and random tosses generated with a real theta of 20%. The distribution quickly converges to probability 1 that $\theta=22.5\%$.


No comments:

Post a Comment