a motivating example

Suppose I have two trick coins:

• one (coin A) comes up heads 75% of the time, and
• the other (coin B) only 25% of the time.

Possible answers:

1. Er, probably coin (A)?

2. Well, \[\begin{aligned} \P\{ \text{6 H in 10 flips} \given \text{coin A} \} &= \binom{10}{6} (.75)^6 (.25)^4 \\ &= 0.146 \end{aligned}\] and \[\begin{aligned} \P\{ \text{6 H in 10 flips} \given \text{coin B} \} &= \binom{10}{6} (.25)^6 (.75)^4 \\ &= 0.016 \end{aligned}\] … so, probably coin (A)?

For a precise answer…

3. Before flipping, each coin seems equally likely. Then

   \[\begin{aligned} \P\{ \text{coin A} \given \text{6 H in 10 flips} \} &= \frac{ \frac{1}{2} \times 0.146 }{ \frac{1}{2} \times 0.146 + \frac{1}{2} \times 0.016 } \\ &= 0.9 \end{aligned}\]

Bayes’ rule

\[\begin{aligned} \P\{B \given A\} = \frac{\P\{B\} \P\{A \given B\}}{ \P\{A\} } , \end{aligned}\]

where

• \(B\): possible model
• \(A\): data
• \(\P\{B\}\): prior weight on model \(B\)
• \(\P\{A \given B\}\): likelihood of data under \(B\)
• \(\P\{B\} \P\{A \given B\}\): posterior weight on \(B\)
• \(\P\{A\}\): total sum of posterior weights

Suppose instead I had 9 coins, with probabilities 10%, 20%, …, 90%; as before I flipped one 10 times and got 6 heads. For each \(\theta\) in \(0.1, 0.2, \ldots, 0.8, 0.9,\) find \[\begin{aligned} \P\{\text{ coin has prob $\theta$ } \given \text{ 6 H in 10 flips } \} . \end{aligned}\]

Question: which coin(s) is it, and how sure are we? (And, what does it mean when we say how sure we are?)

Uniform prior

Weak prior

Strong prior

The likelihood: 6 H in 10 flips

The likelihood: 30 H in 50 flips

The likelihood: 60 H in 100 flips

The likelihood: 6,000 H in 10,000 flips

What is the right answer to the “coin question”?

theta <- (1:9)/10
prior <- rep(1/9, 9)
likelihood <- dbinom(6, size=10, prob=theta)
posterior <- prior*likelihood/sum(prior*likelihood)
names(posterior) <- theta
barplot(posterior, xlab='true prob of heads', main='posterior probability')

Motivating example

Now suppose we want to estimate the probability of heads for a coin without knowing the possible values. (or, a disease incidence, or error rate in an experiment, …)

We flip it \(n\) times and get \(z\) Heads.

The likelihood of this, given the prob-of-heads \(\theta\), is: \[p(z \given \theta) = \binom{n}{z}\theta^z (1-\theta)^{n-z} . \]

How to weight the possible \(\theta\)? Need a flexible set of weighting functions, i.e., prior distributions on \([0,1]\).

Beta-Binomial Bayesian analysis

If \[\begin{aligned} P &\sim \text{Beta}(a,b) \\ Z &\sim \text{Binom}(n,P) , \end{aligned}\]

then “miraculously”,

\[\begin{aligned} (P \given Z = z) \sim \text{Beta}(a+z, b+n-z) . \end{aligned}\]

Discuss:

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

1. What prior to use?

2. Plot the prior and the posterior.

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

4. Best guess at \(\theta\)?

5. How far off are we, probably?

Tools include: rbeta( )

IN CLASS:

# prior: theta could be anything,
# but a little less likely to be 0% or 100%
alpha <- 1.1
beta <- 1.1

prior_samples <- rbeta(1e6, shape1=alpha, shape2=beta)
hist(prior_samples, breaks=100)

3. It’s not totally improbable that theta is 0.5, but it doesn’t look likely.

4. Around 64%.

5. We’re off by around 5-10%, probably.

How do we communicate results?

If we want a point estimate:

1. posterior mean,
2. posterior median, or
3. maximum a posteriori estimate (“MAP”: highest posterior density).

These all convey “where the posterior distribution is”, more or less.

Credible region

Definition: A 95% credible region is a portion of parameter space having a total of 95% of the posterior probability.

Interpretation #1

If we construct a 95% credible interval for \(\theta\) for each of many datasets; and the coin in each dataset has \(\theta\) drawn independently from the prior, then the true \(\theta\) will fall in its credible interval for 95% of the datasets.

Interpretation #2

If we construct a 95% credible interval for \(\theta\) with a dataset, and the distribution of the coin’s true \(\theta\) across many parallel universes is given by the prior, then the true \(\theta\) will be in the credible interval in 95% of those universes.

Interpretation #3

Given my prior beliefs (prior distribution), the posterior distribution is the most rational\({}^*\) way to update my beliefs to account for the data.

But which credible interval?

Definition: The “95% highest density interval” is the 95% credible interval with the highest posterior probability density at each point.