Level 5 · Statistics & Probability: Inference

Inference is the step from the data you have to the world you care about. You measured 40 apples, not the whole harvest; you tested a new teaching method on two classes, not on every school. Level 4 showed that any statistic computed from a sample wobbles from sample to sample, by an amount you can measure (the standard error). This level turns that wobble into answers: a confidence interval says how precisely a quantity has been pinned down, and a hypothesis test asks whether a pattern in the data is larger than chance alone would produce.

These methods are powerful and widely misused, so this level is as much about interpretation as about computation. You will simulate what a p-value is, see that a "95% interval" is a statement about a method rather than about one interval, count how often tests are wrong in each direction, and see why a study with too little data can miss a real effect. Then come the standard tests (t-tests, chi-square, one-way ANOVA), each computed by hand and checked by a simulation that needs no formula at all. The level ends with the other school of inference, Bayesian updating, where the answer is a probability distribution for the unknown quantity itself.

In the Machine Learning path the same questions appear under other names: is model A really better than model B on the test set, or did it get lucky? Cross-validation (ML Level 5) is a way of measuring the wobble of an accuracy score, just as the bootstrap measured the wobble of a median.

Estimators and confidence intervals

An estimator is a rule that turns a sample into a guess about the population: the sample mean estimates the population mean, the sample proportion estimates a probability, and the n - 1 variance of Level 2 estimates the population variance. A single number is a point estimate; it says nothing about how precise it is. A confidence interval adds that.

For a mean, Level 4 gave the ingredients: the sample mean is roughly normal (the central limit theorem) and its standard error is s/\sqrt{n}, with s the sample standard deviation. About 95% of a normal distribution lies within 1.96 standard deviations of its centre, so the 95% confidence interval for the mean is

Here a simulated harvest of 10,000 apples plays the role of the population, so you can see whether the interval works. In real life you would only have the sample.

import math
import random

rng = random.Random(50)
population = [rng.gauss(150, 20) for _ in range(10000)]
sample = rng.sample(population, 40)

n = len(sample)
m = sum(sample) / n
s = math.sqrt(sum((x - m) ** 2 for x in sample) / (n - 1))
se = s / math.sqrt(n)
print(f"sample mean {m:.1f} g, standard error {se:.2f}")
print(f"95% interval: {m - 1.96 * se:.1f} to {m + 1.96 * se:.1f}")
print(f"true mean:    {sum(population) / len(population):.1f}")

sample mean 146.2 g, standard error 2.73
95% interval: 140.9 to 151.6
true mean:    150.1

The interval catches the true mean. What does "95%" promise? Not that this interval has a 95% chance of containing the truth: the truth is a fixed number, and this interval either contains it or not. The promise is about the method: if you repeated the whole study many times, 95% of the intervals it produced would contain the truth. Simulation makes that concrete. Draw 10,000 samples, build an interval from each and count the hits:

import math
import random

def coverage(n, crit, trials=10000):
    rng = random.Random(40)
    hits = 0
    for _ in range(trials):
        s = [rng.gauss(50, 10) for _ in range(n)]
        m = sum(s) / n
        sd = math.sqrt(sum((x - m) ** 2 for x in s) / (n - 1))
        half = crit * sd / math.sqrt(n)
        hits += m - half <= 50 <= m + half
    return hits / trials

print("n=30, 1.96: ", coverage(30, 1.96))
print("n=5,  1.96: ", coverage(5, 1.96))
print("n=5,  2.776:", coverage(5, 2.776))

n=30, 1.96:  0.9342
n=5,  1.96:  0.8738
n=5,  2.776: 0.9467

With 30 values the method delivers about 93%, a little short of its promise. With 5 values, 1.96 covers only about 87%: the interval is too narrow, because s itself is a rough estimate from so few values and is often too small. The fix, found by William Gosset in 1908 while testing barley for a brewery, is to use a slightly wider multiplier from the t distribution with n - 1 degrees of freedom: 2.776 for 4 degrees of freedom, 2.045 for 29, and approaching 1.96 as n grows. The t-test section below computes these numbers.

The width of an interval shrinks like 1/\sqrt{n}: to halve it, collect four times as much data. A wide interval is not a failure; it is an honest report that the data does not pin the answer down.

Common mistakes: saying "there is a 95% probability the true mean is in this interval" (that is a Bayesian statement, and needs the Bayesian machinery at the end of this level). Using 1.96 for small samples. Using the population standard deviation formula (divide by n) for s. And forgetting what the interval does not cover: it measures sampling wobble only, not a biased sample, a broken scale or a question asked the wrong way.

Practice: How Long Do the Batteries Last?, Two Towns, One Question, Which Way Do the Intervals Miss?

Hypothesis tests and p-values

A friend's coin lands heads 60 times in 100 tosses. Is it biased, or is this ordinary luck? A hypothesis test answers a precise version of that question. It starts from the null hypothesis H₀, the boring explanation ("the coin is fair, p = 0.5"), and asks: if H₀ were true, how surprising would data like this be?

Measure the surprise with a test statistic, here the distance of the head count from 50. The p-value is the probability, assuming H₀, of a result at least as extreme as the one observed. Simulation computes it directly: toss a fair coin 100 times, over and over, and count how often the result is 10 or more away from 50.

import math
import random

rng = random.Random(41)
trials = 20000
extreme = 0
for _ in range(trials):
    heads = sum(rng.random() < 0.5 for _ in range(100))
    if abs(heads - 50) >= 10:
        extreme += 1
print("simulated p-value:", extreme / trials)

exact = sum(math.comb(100, k) for k in range(101) if abs(k - 50) >= 10) / 2 ** 100
print("exact p-value:    ", round(exact, 4))

simulated p-value: 0.0574
exact p-value:     0.0569

A fair coin gives a result this lopsided (60 or more heads, or 40 or fewer) about 5.7% of the time. That is unusual but not rare: roughly one friend in eighteen with a perfectly fair coin would see it. This is a two-sided test, because 40 heads would be just as suspicious as 60; a one-sided test counts only one tail and should be chosen before seeing the data.

The exact line of the code is a sum with a condition. In symbols the p-value is

where the condition under the \sum picks which head counts k take part. A condition under a sum becomes an if clause at the end of the generator expression. Press Run to see the p-value as an area: the bars are the probabilities of each head count for a fair coin, and the two tails at least as extreme as 60 are drawn in a second colour. The p-value is the total height of those tail bars.

import math

ks = range(25, 76)
prob = {k: math.comb(100, k) / 2 ** 100 for k in ks}
inner = [k for k in ks if abs(k - 50) < 10]
tails = [k for k in ks if abs(k - 50) >= 10]
plot(inner, [prob[k] for k in inner], kind="bar", label="closer to 50", xlabel="heads in 100 tosses", ylabel="probability if fair")
plot(tails, [prob[k] for k in tails], kind="bar", label="10 or more from 50: the p-value")

The convention is to fix a significance level α in advance, usually 0.05, and to reject H₀ when the p-value is below it. Here 0.057 is just above 0.05, so the data is not strong enough evidence to call the coin biased. That is not the same as showing it is fair: with 100 tosses a coin that lands heads 55% of the time would often pass this test too.

What a p-value is not:

• It is not the probability that H₀ is true. It is computed assuming H₀ is true. (Level 3's warning: P(data | H₀) is not P(H₀ | data).)
• It is not the size or importance of an effect. With a million tosses, a coin landing heads 50.2% of the time gives a tiny p-value, although the bias is too small to matter.
• A p-value of 0.049 and one of 0.051 are nearly the same evidence. The 0.05 line is a convention, not a law of nature.

Every test in this level has the same structure: a null hypothesis, a test statistic that measures departure from it, and the distribution of that statistic when the null is true. Simulation can supply that distribution for almost any test, and formulas supply it quickly for the standard ones.

Common mistakes: choosing one-sided or two-sided (or α) after looking at the data. Reading "not significant" as "no effect". Reading "significant" as "important". And testing many things and reporting only the one that came out below 0.05; the next section shows why that guarantees false discoveries.

Practice: Is the Defect Rate Really 10%?

Type I and type II errors, and power

A test can be wrong in two ways:

The significance level α is the type I error rate: with α = 0.05, a test of a true null hypothesis rejects it 5% of the time, by design. The probability of detecting a real effect is the power of the test, and it depends on how big the effect is and how much data you have. Simulate a coin test at α = 0.05 (using the normal approximation: reject when the head count is more than 1.96 standard errors from n / 2) for fair and biased coins and three sample sizes:

import math
import random

def rejects(heads, n):
    z = (heads - n / 2) / math.sqrt(n / 4)
    return abs(z) > 1.96

rng = random.Random(42)
print("p     n=50   n=100  n=400")
for p in (0.5, 0.55, 0.6):
    row = []
    for n in (50, 100, 400):
        trials = 2000
        r = sum(rejects(sum(rng.random() < p for _ in range(n)), n) for _ in range(trials))
        row.append(r / trials)
    print(p, "  ".join(f"{x:.3f}" for x in row))

p     n=50   n=100  n=400
0.5 0.060  0.058  0.053
0.55 0.140  0.185  0.500
0.6 0.356  0.568  0.984

Read the table row by row. For the fair coin (first row) the test rejects about 5% of the time at every sample size: those are all type I errors, exactly as promised. For a coin with a 55% bias, 100 tosses detect it less than one time in five; the test almost always misses it (a type II error), and even 400 tosses give only a coin flip's chance. A 60% bias is found reliably with 400 tosses.

This is why studies plan their sample size before collecting data: choose the smallest effect worth detecting, and simulate (or calculate) the n that gives, say, 80% power. An underpowered study that finds "no significant effect" has learnt very little. Worse, when an underpowered study does reach significance, it tends to overstate the effect, because only the luckiest samples cleared the bar.

The type I rate also multiplies. Test 20 fair coins at α = 0.05 and the chance that at least one of them looks "significant" is the complement from Level 2:

print(round(1 - 0.95 ** 20, 3))
print(round(1 - (1 - 0.05 / 20) ** 20, 3))

0.642
0.049

Almost two chances in three of a false discovery. The simplest remedy, the Bonferroni correction, tests each of m hypotheses at α / m, which brings the overall false-alarm rate back under 5% (second line), at the cost of power.

Common mistakes: reading a non-significant result from a small study as evidence of no effect. Running many tests (many outcomes, many subgroups, many time points) and reporting the significant ones. And choosing α after seeing the p-value.

Practice: Eight Study Tips, One Lucky Result?, Checking the Experiment Every Day, How Many Patients Does the Trial Need?

Comparing two means: the t-test and the permutation test

Two classes were taught the same topic, one with the old method and one with a new one, and then took the same test:

old = [62, 70, 58, 65, 71, 60, 66, 55, 68, 63, 59, 67]
new = [68, 74, 63, 70, 79, 66, 61, 72, 75, 69, 64, 71]
print(round(sum(old) / 12, 2), round(sum(new) / 12, 2))

63.67 69.33

The new class averaged about 5.7 points higher. Could a difference that size arise by chance if the method made no difference at all? That is H₀.

The cleanest test needs no formula. If the method makes no difference, the labels "old" and "new" are arbitrary: any student's score could just as well have appeared in the other group. So shuffle the 24 scores, split them into two groups of 12 again, and see how often a random split produces a difference at least as large as the real one. This is a permutation test:

import random

old = [62, 70, 58, 65, 71, 60, 66, 55, 68, 63, 59, 67]
new = [68, 74, 63, 70, 79, 66, 61, 72, 75, 69, 64, 71]
observed = sum(new) - sum(old)          # totals: the same test as means, and exact integers
pooled = old + new
rng = random.Random(60)
trials = 10000
count = 0
for _ in range(trials):
    rng.shuffle(pooled)
    diff = sum(pooled[12:]) - sum(pooled[:12])
    if abs(diff) >= abs(observed):
        count += 1
print("observed difference in totals:", observed)
print("permutation p-value:", count / trials)

observed difference in totals: 68
permutation p-value: 0.0146

Only about 1.5% of random relabellings give a gap this large, so the difference is significant at α = 0.05. (Comparing totals rather than means gives the same test, because both groups have 12 students, and avoids floating-point ties.)

The classical two-sample t-test reaches the same conclusion with a formula. Its statistic is the difference in means divided by its standard error, here using a pooled variance s_p^2 that assumes the two groups have similar spread:

Here \bar{x}, s^2 and n are each group's mean, sample variance and size; the test has n_{\text{old}} + n_{\text{new}} - 2 degrees of freedom. In the code below sp2 is s_p^2.

Under H₀ (and assuming roughly normal data), t follows the t distribution. The standard library has no t distribution, but its density is a short formula, and Simpson's rule integrates it:

import math

def t_pdf(t, df):
    c = math.exp(math.lgamma((df + 1) / 2) - math.lgamma(df / 2)) / math.sqrt(df * math.pi)
    return c * (1 + t * t / df) ** (-(df + 1) / 2)

def t_cdf(t, df, steps=1000):
    # P(T <= t) for t >= 0: the left half (0.5) plus the area from 0 to t
    h = t / steps
    s = t_pdf(0, df) + t_pdf(t, df)
    for i in range(1, steps):
        s += (4 if i % 2 else 2) * t_pdf(i * h, df)
    return 0.5 + s * h / 3

def t_crit(df, conf=0.95):
    lo, hi = 0.0, 50.0            # bisection: find t with t_cdf(t) = 0.975
    for _ in range(40):
        mid = (lo + hi) / 2
        if t_cdf(mid, df) < 1 - (1 - conf) / 2:
            lo = mid
        else:
            hi = mid
    return lo

def var(v):
    m = sum(v) / len(v)
    return sum((x - m) ** 2 for x in v) / (len(v) - 1)

old = [62, 70, 58, 65, 71, 60, 66, 55, 68, 63, 59, 67]
new = [68, 74, 63, 70, 79, 66, 61, 72, 75, 69, 64, 71]
na, nb = len(old), len(new)
sp2 = ((na - 1) * var(old) + (nb - 1) * var(new)) / (na + nb - 2)
t = (sum(new) / nb - sum(old) / na) / math.sqrt(sp2 * (1 / na + 1 / nb))
df = na + nb - 2
print("t =", round(t, 3), " df =", df, " p =", round(2 * (1 - t_cdf(abs(t), df)), 4))
print("t critical values:", [round(t_crit(d), 3) for d in (4, 22, 29, 1000)])

t = 2.691  df = 22  p = 0.0133
t critical values: [2.776, 2.074, 2.045, 1.962]

The t-test's p-value (0.013) is close to the permutation test's, and the last line produces the multipliers promised in the confidence-interval section: 2.776 for 4 degrees of freedom, approaching 1.96 for large samples. Two variants are worth knowing. Welch's t-test drops the equal-spread assumption (it is the safer default when the groups' standard deviations differ a lot). A paired t-test is for two measurements on the same people (before and after): take each person's difference and do a one-sample test of whether their mean is 0, which removes the person-to-person variation and is far more powerful.

The permutation test assumed only that the labels were exchangeable under H₀; the t-test also assumed roughly normal data. When they disagree, the permutation test is usually the one to trust, and simulation made it possible.

Common mistakes: using a two-sample test on paired data (it throws away the pairing and loses power). Using the pooled test when one group is far more spread out than the other. And forgetting that a significant difference between two classes is only evidence about the method if students were assigned to the classes at random; otherwise the classes may differ in other ways (a confounder again).

Practice: Does Coffee Speed Up Reactions?, Two Couriers and One Very Late Parcel, When Pooling the Variances Goes Wrong

The chi-square test for counts

The t-test compares means. For counts in categories the standard tool is the chi-square (χ²) test. It compares the observed count O in each category with the count E expected under H₀:

In the code below that is sum((o - expected) ** 2 / expected for o in observed): the loop runs over the categories, and the summand is the same fraction.

Each term is a squared gap scaled by how big the count should be, so a gap of 6 matters more when 20 were expected than when 2000 were.

Goodness of fit. A die rolled 120 times gives the counts below. A fair die should give about 20 of each. Is this die loaded?

import random

observed = [14, 22, 17, 26, 13, 28]
n = sum(observed)
expected = n / 6
x2 = sum((o - expected) ** 2 / expected for o in observed)
print("X^2 =", round(x2, 2))

rng = random.Random(70)
trials = 10000
count = 0
for _ in range(trials):
    counts = [0] * 6
    for _ in range(n):
        counts[rng.randrange(6)] += 1
    if sum((c - expected) ** 2 / expected for c in counts) >= x2:
        count += 1
print("simulated p-value:", count / trials)

X^2 = 9.9
simulated p-value: 0.0778

The counts look uneven (13 against 28), but a fair die produces an X² this large about 8% of the time, so the evidence is not strong enough at α = 0.05. Fair dice are more uneven over 120 rolls than intuition expects. Without a simulation you would compare X² with the chi-square distribution with (categories − 1) = 5 degrees of freedom, whose 5% critical value, from a table, is 11.07; 9.9 falls short of it, in agreement.

Independence in a two-way table. Level 3's homework table asked whether passing depended on doing the homework. Under H₀ (independence), each cell's expected count is row total × column total / grand total; the test then measures how far the table is from that:

import math

table = [[90, 30],     # homework:    passed, failed
         [40, 40]]     # no homework: passed, failed
rows = [sum(r) for r in table]
cols = [sum(c) for c in zip(*table)]
total = sum(rows)

x2 = 0
for i in range(2):
    for j in range(2):
        e = rows[i] * cols[j] / total
        print(f"cell {i},{j}: observed {table[i][j]}, expected {e:.0f}")
        x2 += (table[i][j] - e) ** 2 / e
print("X^2 =", round(x2, 2))
print("p-value =", f"{math.erfc(math.sqrt(x2 / 2)):.5f}")

cell 0,0: observed 90, expected 78
cell 0,1: observed 30, expected 42
cell 1,0: observed 40, expected 52
cell 1,1: observed 40, expected 28
X^2 = 13.19
p-value = 0.00028

A 2 × 2 table has (2 − 1)(2 − 1) = 1 degree of freedom, and for one degree of freedom the chi-square tail probability has a closed form, erfc(sqrt(x / 2)), because X² is then the square of a single standard normal. The p-value is about 0.0003: passing is clearly associated with doing the homework. Associated, not caused: students who do their homework may differ in many other ways.

Common mistakes: running the test on percentages instead of counts (X² depends on the sample size, and percentages hide it). Using it when expected counts are small (a common rule: every expected count at least 5; otherwise simulate, or merge categories). Counting the same person in several cells. And forgetting that a large table with thousands of people will find "significant" associations that are far too small to matter; report the proportions as well as the p-value.

Practice: Does the Commute Depend on the District?

One-way ANOVA: more than two groups

Three fertilisers were each given to five plants, and the heights (cm) after a month were:

groups = {"A": [20, 22, 19, 24, 21], "B": [25, 27, 23, 26, 24], "C": [21, 23, 22, 20, 25]}
for name, g in groups.items():
    print(name, sum(g) / len(g))

A 21.2
B 25.0
C 22.2

Why not run three t-tests (A–B, A–C, B–C)? Because each has a 5% false-alarm rate, and together the chance of at least one false alarm is well above 5% (the multiple-testing problem from earlier). One-way analysis of variance (ANOVA) asks a single question instead: are all the group means equal?

The idea is in the name. Split the variation in the data into two parts: between groups (how far the group means are from the overall mean) and within groups (how far plants are from their own group's mean). If the fertilisers make no difference, the group means should differ only as much as the within-group noise implies. The F statistic compares the two, each divided by its degrees of freedom. With k groups and N values in all, n_g values in group g, group means \bar{x}_g and overall mean \bar{x}:

The double sum in SS_{\text{within}} (over groups, then over the values in each group) becomes a sum inside a sum in f_stat below, one for for each \sum.

Compute it, then get a p-value with a permutation test: shuffle the 15 heights among the three groups and see how often F comes out at least as large.

import random

def f_stat(groups):
    values = [x for g in groups for x in g]
    N, k = len(values), len(groups)
    grand = sum(values) / N
    ss_between = sum(len(g) * (sum(g) / len(g) - grand) ** 2 for g in groups)
    ss_within = sum(sum((x - sum(g) / len(g)) ** 2 for x in g) for g in groups)
    return (ss_between / (k - 1)) / (ss_within / (N - k))

groups = [[20, 22, 19, 24, 21], [25, 27, 23, 26, 24], [21, 23, 22, 20, 25]]
F = f_stat(groups)
print("F =", round(F, 3))

rng = random.Random(80)
pooled = [x for g in groups for x in g]
trials = 10000
count = 0
for _ in range(trials):
    rng.shuffle(pooled)
    if f_stat([pooled[0:5], pooled[5:10], pooled[10:15]]) >= F - 1e-9:
        count += 1
print("permutation p-value:", count / trials)

F = 5.879
permutation p-value: 0.021

(The - 1e-9 makes a shuffle that reproduces exactly the same F count as "at least as large" despite floating-point rounding.) With p ≈ 0.02 the means are unlikely to be all equal. The classical test compares F with the F distribution with (2, 12) degrees of freedom, whose 5% critical value is 3.89, and reaches the same verdict.

ANOVA says that the means differ, not which. Here B looks like the odd one out; confirming that needs follow-up comparisons with a correction for multiple testing. ANOVA's classical form also assumes roughly normal data with similar spread in each group; the permutation version needs only that, under H₀, the group labels are exchangeable.

Common mistakes: running all the pairwise t-tests instead. Concluding from a significant F that every group differs from every other. Using ANOVA on paired or repeated measurements of the same plants (that needs a repeated-measures design). And, as always, attributing the difference to the fertiliser when plants were not assigned to fertilisers at random.

Practice: Three Feeds for the Tomato Plants

Bayesian updating

Everything so far treated the unknown quantity (a mean, a coin's bias) as fixed and asked how the data would vary. Bayesian inference turns this around: it describes what you believe about the unknown quantity with a probability distribution, and updates that distribution with Bayes' rule as data arrives, exactly as Level 3 updated the probability of a disease after a test result:

The sign \propto means "proportional to": the two sides are equal once the right-hand side is rescaled so that it adds up to 1.

The simplest way to compute it is on a grid: list 101 possible values of a basketball player's free-throw success rate θ from 0 to 1, give each a prior weight, multiply by the probability of each observed shot, and renormalise so the weights add to 1. Start with a flat prior (every θ equally plausible) and update one shot at a time:

grid = [i / 100 for i in range(101)]

def normalise(w):
    total = sum(w)
    return [x / total for x in w]

def summary(post):
    mean = sum(p * w for p, w in zip(grid, post))
    cum, lo, hi = 0, None, None
    for p, w in zip(grid, post):
        cum += w
        if lo is None and cum >= 0.025:
            lo = p
        if hi is None and cum >= 0.975:
            hi = p
    return f"mean {mean:.3f}, most likely {grid[post.index(max(post))]}, 95% interval {lo} to {hi}"

shots = [1, 1, 0, 1, 1, 1, 0, 1, 0, 1]      # 7 made out of 10
post = normalise([1] * 101)
for s in shots:
    post = normalise([w * (p if s else 1 - p) for p, w in zip(grid, post)])
print("flat prior:    ", summary(post))

league = normalise([p ** 10 * (1 - p) ** 10 for p in grid])   # belief centred on 0.5
print("league prior:  ", summary(league))
post2 = normalise([w * p ** 7 * (1 - p) ** 3 for p, w in zip(grid, league)])
print("after 10 shots:", summary(post2))
post3 = normalise([w * p ** 70 * (1 - p) ** 30 for p, w in zip(grid, league)])
print("after 100 shots:", summary(post3))

flat prior:     mean 0.667, most likely 0.7, 95% interval 0.39 to 0.89
league prior:   mean 0.500, most likely 0.5, 95% interval 0.3 to 0.7
after 10 shots: mean 0.562, most likely 0.57, 95% interval 0.39 to 0.73
after 100 shots: mean 0.664, most likely 0.67, 95% interval 0.58 to 0.74

With a flat prior, the most likely value is 0.7 (the observed rate) and the 95% credible interval runs from 0.39 to 0.89. Unlike a confidence interval, this one does mean "given the data and the prior, there is a 95% probability that θ is in here", which is what people usually wanted a confidence interval to mean.

The price is the prior. A scout who knows that most players shoot around 50% might use the "league" prior, as confident as having already seen 10 made out of 20. After 10 shots the posterior (mean about 0.56) sits between the prior and the data; after 70 made out of 100 it has moved to about 0.66, most of the way to the observed 0.7, and the interval has narrowed to 0.58–0.74. With little data the prior matters a lot; with plenty of data, reasonable priors end up in nearly the same place. Updating shot by shot and updating once with all ten shots give identical posteriors, because multiplication does not care about order.

Press Run to see the prior's influence fade as the data grows:

grid = [i / 100 for i in range(101)]

def normalise(w):
    total = sum(w)
    return [x / total for x in w]

league = normalise([p ** 10 * (1 - p) ** 10 for p in grid])
plot(grid, league, label="league prior", xlabel="success rate θ", ylabel="weight")
plot(grid, normalise([w * p ** 7 * (1 - p) ** 3 for p, w in zip(grid, league)]), label="after 7 of 10")
plot(grid, normalise([w * p ** 70 * (1 - p) ** 30 for p, w in zip(grid, league)]), label="after 70 of 100")

Level 6 shows that for this kind of problem the grid is unnecessary: a Beta prior gives a Beta posterior in closed form (a conjugate prior). The grid, though, works for any prior and any likelihood you can write down. The Naive Bayes classifier of the Machine Learning path is a Bayesian update too, with class labels in place of θ.

Common mistakes: a prior that gives zero weight to a value (no amount of data can revive it; posteriors multiply). Forgetting to renormalise. Choosing the prior after seeing the data to get the answer you wanted. And a grid too coarse to resolve a narrow posterior: with thousands of observations, 101 points are too few.

Practice: The New Goalkeeper's Save Rate, How Likely Is the Seed Better Than Half?

Practice

These problems are required on this path at this level; the path adds one problem from each of the level's topics, picked for you and kept.

• How Long Do the Batteries Last?
• Two Towns, One Question
• How Many Wristbands Were Sold?
• Is the Defect Rate Really 10%?
• Does Coffee Speed Up Reactions?
• When Pooling the Variances Goes Wrong

Checklist

You can now:

• Explain what an estimator is and build a confidence interval for a mean, with 1.96 for large samples and a t multiplier for small ones.
• Describe what "95% confidence" promises, and show it by simulating the coverage of many intervals.
• Compute a p-value by simulation and exactly, and say what a p-value is and is not.
• Name type I and type II errors, simulate the power of a test for different effect sizes and sample sizes, and explain why many tests need a correction.
• Compare two groups with a permutation test and with a two-sample t-test computed by hand, and choose a paired test for paired data.
• Run a chi-square test for goodness of fit and for independence in a two-way table.
• Compute the F statistic of a one-way ANOVA and get its p-value by permutation.
• Update a prior to a posterior on a grid, read a credible interval, and explain how the prior's influence fades with more data.