\[%% % Add your macros here; they'll be included in pdf and html output. %% \newcommand{\R}{\mathbb{R}} % reals \newcommand{\E}{\mathbb{E}} % expectation \renewcommand{\P}{\mathbb{P}} % probability \DeclareMathOperator{\sd}{sd} \DeclareMathOperator{\var}{var} \DeclareMathOperator{\Normal}{Normal} \DeclareMathOperator{\Poisson}{Poisson} \DeclareMathOperator{\Beta}{Beta} \DeclareMathOperator{\Binom}{Binomial} \DeclareMathOperator{\Gam}{Gamma} \DeclareMathOperator{\Exp}{Exponential} \newcommand{\given}{\;\vert\;} \]

Question: We flip an odd-looking coin 100 times, and get 65 heads. What is it’s true* probability of heads?

What is True?

The long-run proportion of heads (if we flipped it a lot).

What prior to use?

We’ll use \(\alpha = \beta = 1\) : a flat prior.

Analysis:

We know that the posterior distribution on \(P\), the true probability of heads, is Beta(\(\alpha + z\), \(\beta + n - z\)). Here \(z=65\) and \(n=100\). We’ll make a plot of the density function.

z <- 65
n <- 100
prior_alpha <- 1
prior_beta <- 1
post_alpha <- prior_alpha + z
post_beta <- prior_beta + n - z
x <- seq(0, 1, length.out=1000)
dens <- dbeta(x, shape1=post_alpha, shape2=post_beta)
plot(x, dens, type='l',
     xlab='prob of heads', ylab='posterior density')
plot of chunk bayes
plot of chunk bayes

Is it reasonable that \(\theta=1/2\)?

It looks like not: the posterior probability that \(P \le 1/2\) is 0.0013269.

Best guess at \(\theta\)?

We’ll take the maximum a posteriori estimate (MAP), which is 0.65.

How far off are we, probably?

To communicate our uncertainty in that estimate, we’ll compute a 95% “highest density” credible interval (HDI). We do this with the HDInterval package.

library(HDInterval)
hdi <- hdi(qbeta, credMass=.95, shape1=post_alpha, shape2=post_beta)

We can be 95% certain that the true \(\theta\) is between 0.5542689 and 0.7382436.