Counting ways
A posterior is a list of candidate answers with a plausibility attached to each. For Northgate the candidates are redemption rates between 0 and 1, and the plausibility of each is how many ways that rate could have produced the data, weighted by how much you believed it beforehand. Statistical Rethinking's second lecture, "Garden of Forking Data", builds the whole idea from that count. Tonight starts there because once you have seen a posterior as a count on a grid, sampling from it stops being mysterious: a draw is just a rate picked in proportion to its count.
★ What you walk out with tonight. Two ways to build the same posterior, a grid and a set of 20,000 draws, and the check that they agree with the exact Beta to the decimals reported: mean 0.100, interval 0.020 to 0.222. The habit of reading every summary as a quantile of the draws. And the one piece of tooling history that will save you from misreporting a number: ArviZ 1 moved its default interval to 89 percent equal-tailed, and a notebook that upgraded without noticing changed its answers while its prose stayed the same.
LiveThe garden of forking data, on a grid4 min▶
Lay out candidate rates on a grid: 0.00, 0.05, 0.10, and so on up to 1.00. For each candidate, ask two things. How plausible was this rate before the data, which is the prior evaluated at that point. And how many ways could this rate have produced one redemption in eighteen vouchers, which is the likelihood evaluated at that point. Multiply the two. The candidate at 0.00 gets nothing, because a rate of zero cannot produce a redemption at all. The candidate at 0.05 gets a great deal, because one in eighteen is roughly what a rate of 0.05 does. The candidate at 0.50 gets almost nothing, because nine redemptions were expected and one arrived.
Divide every product by the sum of all of them so the column adds to one, and the column is the posterior. Read the mean off it as a weighted average of the grid points. With a flat prior and a fine enough grid the mean lands at 0.100, which is exactly what Beta(2, 18) says, and it should: the grid is computing the same object the algebra computes, one point at a time, without knowing the algebra exists.
Halloway's analyst had used the conjugate formula for a year without being able to say what a posterior was. The grid, done once in a notebook for Northgate, was the first time she saw the thing itself: a column of plausibilities, tallest near 0.05, still meaningfully tall at 0.15, close to nothing above 0.30. "The posterior is that column" turned out to be the sentence that made session 4 make sense.
Self-studyThe words that will keep coming back4 min▶
- Grid approximation. Evaluate prior times likelihood at a finite set of candidate rates, normalise so they sum to one, read the posterior off. Exact in the limit of a fine grid; the teaching tool, not the production tool.
- Normalise. Divide every product by their total. Nothing about the shape changes; the column now sums to one so its entries are probabilities.
- Draw. One candidate rate picked at random in proportion to its posterior plausibility. 20,000 of them is the working set for every number in these sessions.
- Quantile. The value below which a stated fraction of the draws fall. The 5.5 percent quantile of Northgate's draws is 0.020; the 94.5 percent quantile is 0.222.
- Equal-tailed interval. Cut the same fraction from each end. For 89 percent that is 5.5 percent below and 5.5 percent above. ArviZ calls it
"eti". - Highest density interval, HDI. The narrowest region holding the stated fraction. On a skewed posterior it sits in a different place from the equal-tailed one.
- stats.ci_prob, stats.ci_kind. The two ArviZ 1 rcParams that decide what a summary reports by default: 0.89 and
"eti". Both changed at the 1.0 release. - Credible interval. The Bayesian interval, the one that carries the plain sentence: there is an 89 percent probability the rate is in here. The frequentist version is a different sentence, and learn-statistics-with-phoebe owns it.
Prior times likelihood, one candidate at a time
10 min liveNothing tonight needs the Beta formula. A grid of candidate rates, the prior at each, the likelihood of one in eighteen at each, the product normalised: that is the posterior, and making the grid finer makes it smoother without changing what it is. The mean read off a fine grid is 0.100 to three decimals, the same as the exact Beta(2, 18). The picture below is the whole method.
LiveEvaluate, multiply, normalise, read4 min▶
Four moves and no formula. Evaluate the prior at every grid point; with Beta(1, 1) that is the same value everywhere, which is what "flat" means. Evaluate the likelihood at every grid point: the rate to the power of the redemptions, times one minus the rate to the power of the non-redemptions, so for Northgate rate to the one times one minus rate to the seventeen. Multiply the two columns. Divide by the total so the column sums to one.
Every summary is now a calculation on that column. The mean is each grid rate times its posterior plausibility, summed. A probability that the rate exceeds some value is the plausibility summed over the grid points above it. An interval is found by walking up the column until 5.5 percent of the mass has been passed, then again until 94.5 percent has. On a fine grid the mean is 0.100. The exact Beta says 0.100. The grid did not know it was computing a Beta.
Self-studyWhy the grid is a teaching tool and not the tool3 min▶
The grid works beautifully for one parameter and stops working almost immediately after. Add a second parameter and every candidate for the first has to be paired with every candidate for the second, so the number of cells is the grid size multiplied by itself. Add a third and it is multiplied again. The hierarchical model in session 6 has a parameter per outlet plus two more for how the outlets relate, and no grid fine enough to be useful fits in memory.
That is the reason session 4 exists. Markov chain Monte Carlo produces draws from a posterior without ever evaluating it everywhere, and the draws are read exactly the way tonight's are: mean, quantiles, probabilities above a threshold. Learn to read draws tonight, on a posterior whose exact answer you can check, and session 4 only has to teach you how the draws are made and how to tell when the maker has failed.
Every summary is a quantile
11 min liveTake 20,000 draws from Beta(2, 18). Their mean is 0.100. Sort them and the 5.5 percent quantile is 0.020 and the 94.5 percent quantile is 0.222. That is the 89 percent equal-tailed interval, and it matches the exact Beta to the decimals reported. The interval is not the only way to hold 89 percent of a posterior, though, and on a skewed curve like Northgate's the two standard ways land in visibly different places. Both are shown below; one is taught.
LiveFrom draws to numbers: the mean and two quantiles4 min▶
A draw from the posterior is a candidate rate picked in proportion to its plausibility, so a large pile of draws is the posterior in a form that any summary can be read from by sorting and counting. The mean is the average of the pile. The 5.5 percent quantile is the value that 5.5 percent of draws fall below; the 94.5 percent quantile is the value that 94.5 percent fall below. Between them sit 89 percent of the draws, and that is the interval reported on every page here: 89 percent, equal-tailed, from 20,000 draws.
For Northgate under the flat prior the pile gives mean 0.100 and quantiles 0.020 and 0.222. The exact Beta(2, 18) gives the same numbers to the decimals shown, which is the check that matters tonight: sampling a posterior you can write down agrees with the algebra, so sampling a posterior you cannot write down, which is every posterior from session 6 onward, is not a leap of faith. Westfield under the flat prior, 127 redemptions in 900, has a posterior mean of 0.142 by the same route.
Halloway's analyst reported Northgate to the regional manager as "mean 0.100, 89 percent interval 0.020 to 0.222, equal-tailed, from 20,000 draws". The manager only used the first number. Six weeks later, when Northgate's true rate turned out to be near the top of that interval, the sentence was still in the deck, and it had said all along that 0.192 was inside what the data allowed.
LiveWhy 89 and why equal-tailed: the default moved, and it matters4 min▶
ArviZ 1 ships two settings that decide what every summary reports when you do not say otherwise. stats.ci_prob is 0.89, where it used to be 0.94. stats.ci_kind is "eti", equal-tailed, where it used to be the highest density interval. The migration guide gives the reason for the level in one sentence: "Using 0.89 produces intervals with lower variability, leading to more stable summaries." Nothing about 89 is sacred; it is a default chosen for stability across runs, and the argument it makes is a quiet one that intervals are summaries of draws and should not wobble between reruns more than they have to.
The change of kind is the one that bites. Equal-tailed cuts 5.5 percent from each end. HDI finds the narrowest region holding 89 percent, which on Northgate's right-skewed posterior starts closer to zero and ends lower, as the figure shows. A third party, the CoursIA repository in issue #15592, documented two notebooks whose intervals silently became equal-tailed on upgrade while the prose beside them still said HDI. The numbers changed; the words did not; nobody was told. So the rule here is simple. Both kinds exist and you should be able to name both. The equal-tailed one is taught and reported on every page, and every page says so, because the kind has to travel with the number or the number is not fully stated.
Self-studyThe posterior is the answer; the mean is a summary of it3 min▶
Nothing you report is the posterior. The mean is one number pulled from it, chosen because it is easy to say. The interval is two numbers pulled from it, chosen because they carry the width the mean throws away. When somebody asks for "the number", the honest reply is that the answer is a distribution and here is the summary that fits the decision they are making; for a budget the mean, for a risk the lower quantile, for a promise the interval.
The sentence that carries the interval into a meeting is: "there is an 89 percent probability the rate is between 0.020 and 0.222". That is a statement about this posterior under this prior and this model, and it is exactly what it sounds like, which is the reason it is worth having. Whether 89 percent statements like it turn out right 89 percent of the time across many forecasts is a different question, a measurement rather than a derivation, and learn-decision-intelligence-with-phoebe owns it. Tonight's job is to produce the sentence correctly, kind and level attached. Session 9 is about saying it to the people who act on it.
Grid it, then draw it
16 min build-alongPython, numpy and scipy, nothing else. Build Northgate's posterior on a grid and read its mean, then draw 20,000 samples from the exact Beta(2, 18) and read the same mean and the two quantiles. The two routes agree, and from session 4 onward only the second one will be available.
Grid approximation for Northgate. One redemption in eighteen, flat prior. The likelihood is the binomial kernel with its constant dropped; the constant would vanish in the normalising step anyway.
import numpy as np # Northgate: 1 redemption in 18 vouchers, flat Beta(1, 1) prior y, n = 1, 18 grid = np.linspace(0, 1, 1001) prior = np.ones_like(grid) # Beta(1, 1): the same everywhere likelihood = grid**y * (1 - grid)**(n - y) # ways each rate produces 1 in 18 posterior = prior * likelihood posterior = posterior / posterior.sum() # normalise: the column sums to 1 mean = (grid * posterior).sum() # weighted average of the grid print(round(mean, 3)) # 0.100
Draws from the exact posterior. Beta(1 + 1, 1 + 17) is Beta(2, 18). Take 20,000 draws and read the mean and the 5.5 and 94.5 percent quantiles. The bench's figures are 0.020 and 0.222; a different seed can move the last digit by one, which is the "lower variability" the ArviZ default is about.
import numpy as np from scipy.stats import beta posterior = beta(2, 18) # flat prior + 1 in 18 draws = posterior.rvs(20000) lo, hi = np.quantile(draws, [0.055, 0.945]) # 89 percent, equal-tailed print(round(draws.mean(), 3)) # 0.100 print(round(lo, 3), round(hi, 3)) # 0.020 0.222
Check the draws against the algebra. The exact quantiles come straight from the distribution object. They should match the sampled ones to the decimals reported; if they do not, the draw count is too small, not the method.
exact_lo, exact_hi = posterior.ppf([0.055, 0.945]) print(round(posterior.mean(), 3)) # 0.100 print(round(exact_lo, 3), round(exact_hi, 3)) # 0.020 0.222
Westfield by the same route. 127 redemptions in 900 under the flat prior is Beta(128, 774). The mean is 0.142. Note how much less the draws spread: the interval that was 0.2 wide at Northgate is a sliver here, and that difference is the whole of session 2's argument.
westfield = beta(1 + 127, 1 + 900 - 127) # Beta(128, 774) w = westfield.rvs(20000) print(round(w.mean(), 3)) # 0.142
Write the sentence, kind and level attached. "Northgate, flat prior: mean 0.100, 89 percent equal-tailed interval 0.020 to 0.222, from 20,000 draws." Every later session reports in exactly this form, and the two words after the percentage are not optional.
Set the kind explicitly even when it matches the default. A notebook that relies on stats.ci_kind being "eti" is a notebook whose intervals will change the next time the default does, and the prose will not know. Quantiles at 0.055 and 0.945, written out, cannot be upgraded out from under you.
Three questions
A colleague computes Northgate's posterior on a grid of 1001 rates and gets a mean of 0.100, then says grid approximation must be a different method from the Beta(2, 18) they learned in session 1. What is true?
Prior times likelihood, normalised, is the posterior by definition. The Beta is what that product happens to be for a Beta prior and a binomial likelihood; the grid computes the same product point by point without knowing the name. The mean agrees at 0.100 because they are one object.
You upgrade a notebook to ArviZ 1 and its summary intervals change, though the code did not. The prose beside them still says "94 percent HDI". What happened, and what should the prose now say?
This is the documented CoursIA case: intervals silently became equal-tailed on upgrade while the text still said HDI. The fix is to state kind and level with every number, and to write the quantiles out rather than lean on a default that can move.
On Northgate's right-skewed posterior the 89 percent equal-tailed interval is 0.020 to 0.222. Where does the 89 percent HDI sit relative to it, and which is reported here?
Equal-tailed cuts 5.5 percent from each end regardless of shape. HDI chases the highest part of the curve, which on a skewed posterior shifts the whole interval toward the peak and makes it shorter. Same 89 percent, different interval, so the kind has to travel with the number.
Before session 4
30 minRun tonight's grid for Northgate under each of session 2's other three priors, Beta(2, 14), Beta(12, 88) and Beta(40, 60), and read the mean and the 5.5 and 94.5 percent quantiles off each grid. Check them against the Northgate table from session 2. Then look at the shape of the grid column under Beta(12, 88) and ask yourself where the truth of 0.192 sits on it. Session 4 replaces the grid with a sampler; bring the habit of reading quantiles, because that is the part that does not change.