ML8.4 · LESSON 4 OF 6

Posterior prediction

Update a prior, compare posterior means with MAP estimates and interpret predictive uncertainty.

THE IDEA

Predict by averaging over theta

What is the chance the next flip is heads? Given theta it would be theta, but you do not know theta; you have a posterior. The honest answer averages theta over that posterior, which for a Beta(a, b) is simply its mean: P(heads next) = a / (a + b). With a flat prior and 3 heads, 0 tails the posterior is Beta(4, 1), so the answer is 0.8.

The shortcut of picking one estimate and plugging it in (here the frequency 3/3 = 1) treats theta as known exactly. It says tails are impossible, so it is certain about the next flip when three flips cannot justify that.

posterior predictive (averages over theta)plug in the point estimate

heads seen3tails seen0next flips m10

prior Beta(1, 1) -> posterior Beta(4, 1) next flip is heads: predictive = 4/5 = 0.8 plug-in = 3/3 = 1 all 10 next flips heads: predictive = 0.286 plug-in = 1 spread of the head count (sd): predictive = 2 plug-in = 0

The point estimate says tails are impossible. The posterior still allows them.

Blue: predictive distribution of heads in the next m flips, averaged over theta. Orange: the same question with the point estimate plugged in. Try 3 heads, 0 tails, then 3 and 3, then 40 and 40.

Where the extra spread comes from

A plugged-in binomial only has coin-flip noise. The predictive adds uncertainty about theta itself, so the head count is more spread out. With 3 heads and 3 tails (posterior Beta(4, 4)) and 10 future flips, the predictive standard deviation of the head count is 2.236against 1.581 for the plug-in. With lots of data the two nearly coincide: the parameter becomes known and only flip noise is left.

Predict first. Flat prior, 3 heads and 0 tails. How likely is a run of 10 more heads?

Common mistake

Reporting the plug-in probability as if it were the predictive one. Look at the case of 0 heads in 4 flips: the point estimate says the next flip is heads with probability 0, while the flat-prior predictive says 1/6. Any decision that divides by that probability, or bets everything on it, fails on the first heads.

Chapter context

Probabilistic modeling represents uncertainty explicitly instead of presenting every fitted number as certain. Priors, likelihoods and approximation choices are assumptions that should be exposed rather than hidden behind a library call.

TRY IT YOURSELF

Part 1. Do it hands-on

Flat Beta(1, 1) prior. Move the sliders until each task ticks.

posterior predictive (averages over theta)plug in the point estimate

heads seen0tails seen0next flips m5

prior Beta(1, 1) -> posterior Beta(1, 1) next flip is heads: predictive = 1/2 = 0.5 plug-in = undefined (no flips) all 5 next flips heads: predictive = 0.167 plug-in = undefined spread of the head count (sd): predictive = 1.708 plug-in = undefined

To do: Find data where the point estimate says the next flip is heads with certainty (100%) but the posterior predictive is below 0.85
To do: Make the point estimate give at least 0.5 for a run of all heads in the next m flips, while the posterior predictive gives less than 0.35 (use m of 10 or more)
To do: With m = 10, find data (at least one head and one tail) where the predictive spread of the head count is at least 20% larger than the point-estimate binomial

Part 2. Do it in code

predictive_heads(a, b, heads, tails): probability the next flip is heads under a Beta(a, b) prior: (a + heads) / (a + b + heads + tails).
plugin_heads(heads, tails): the point estimate heads / (heads + tails). Raise ValueError when there are no flips.
prob_all_heads(m, a, b, heads, tails): probability that the next m flips are all heads, averaging over the posterior. Multiply the one-step predictive probabilities, updating after each imagined head. m = 0 gives 1.
mc_predictive(a, b, heads, tails, n, seed): estimate the next-flip probability by simulation: draw n values of theta from Beta(a + heads, b + tails) with np.random.default_rng(seed), flip one coin for each, return the fraction of heads.

import numpy as np

def predictive_heads(a, b, heads, tails):
    pass

def plugin_heads(heads, tails):
    pass

def prob_all_heads(m, a, b, heads, tails):
    pass

def mc_predictive(a, b, heads, tails, n, seed):
    pass

print(predictive_heads(1, 1, 0, 4), plugin_heads(0, 4))   # 1/6 and 0

Ctrl/⌘ + Enter runs it

Stuck, or something not working? You can come back to it any time.

CHECK YOUR UNDERSTANDING

Which uncertainty disappears if a fitted parameter is treated as certain?

AThe uncertainty about future coin-flip noise, since the parameter fixes the outcome.BThe prior, because a fitted parameter no longer needs one.CThe uncertainty about that parameter, so the prediction becomes overconfident.

Chapter reference: math, worked example & code

This section belongs to the whole chapter, so it is the same on every lesson in it and may cover a different lesson's topic. The code for the experiment you just ran is inside the experiment itself.

Math and notation

For Bernoulli success probability theta and a Beta prior:

p(\theta\mid D)=\operatorname{Beta}(a+s,b+f),\quad E[\theta\mid D]=\frac{a+s}{a+b+s+f}.

For Beta(alpha,beta) with both parameters greater than one, the MAP is $(\alpha-1)/(\alpha+\beta-2)$ .

For a zero-mean Gaussian-process regression model with independent noise variance sigma-n squared, a test prediction has

\mu_*=k_*^T(K+\sigma_n^2I)^{-1}y,

v_*=k_{**}-k_*^T(K+\sigma_n^2I)^{-1}k_*.

K is the training covariance matrix, k-star connects the test input with training inputs, and k-star-star is its prior variance. This is latent-function variance; observation noise adds another term for a noisy predicted observation. Practical code solves systems instead of explicitly computing an inverse.

Worked example

Start with Beta(2,2), then observe three successes and one failure. The posterior is Beta(5,3). Its mean is 5/8=0.625. Its MAP is 4/6≈0.6667. Both summarize the same distribution but answer different questions.

Reference code

The following example uses the Python standard library and does not call a remote service.

from math import isclose
a, b = 2+3, 2+1
assert (a,b) == (5,3)
assert a/(a+b) == .625
assert isclose((a-1)/(a+b-2), 2/3)

Python examples are provided for study and local execution; this page does not execute Python.

Common mistakes

A posterior estimate is conditional on the chosen model. MAP is not the same as a posterior mean. A small uncertainty under a misspecified model need not mean a reliable real-world answer.