Jacob Bernoulli’s law of large numbers: frequencies settle down

Bernoulli’s “golden theorem” (1713): in \(n\) independent trials of an event with probability \(p\), the relative frequency of the event converges to \(p\) as \(n\) grows. Here a biased coin with \(p = 0.3\) is tossed a million times, and the running fraction of heads homes in on 0.3 while its typical error shrinks like \(\sqrt{p(1-p)/n}\).

import matplotlib.pyplot as plt
import numpy as np

from mathematicskit.probability import Binomial
from mathematicskit.probability.systems.limit_theorems import law_of_large_numbers_trace

p = 0.3
coin = Binomial(n=1, p=p)  # one Bernoulli trial: 1 = heads, 0 = tails

The relative frequency of heads converges to p

n_values = np.unique(np.logspace(0, 6, 400).astype(int))
fig, ax = plt.subplots()
for seed in range(5):
    ax.semilogx(n_values, law_of_large_numbers_trace(coin, n_values, seed=seed), lw=1)
band = 2.0 * np.sqrt(p * (1 - p) / n_values)
ax.semilogx(n_values, p + band, "k--", n_values, p - band, "k--", lw=1)
ax.axhline(p, color="k", lw=1.5)
ax.set_ylim(0.0, 1.0)
ax.set_xlabel("number of tosses $n$")
ax.set_ylabel("relative frequency of heads")
ax.set_title(r"Law of large numbers: frequency $\to p = 0.3$ (dashed: $\pm 2$ s.d.)")

Text(0.5, 1.0, 'Law of large numbers: frequency $\\to p = 0.3$ (dashed: $\\pm 2$ s.d.)')

Bernoulli’s statement: \(P(|\bar{X}_n - p| > \varepsilon) \to 0\)

checkpoints = [10, 100, 1000, 10000, 1000000]
for n, freq in zip(checkpoints, law_of_large_numbers_trace(coin, checkpoints, seed=0)):
    print(f"n={n:>8}: relative frequency = {freq:.4f} (error {abs(freq - p):.4f})")

eps = 0.02
for n in (100, 1000, 10000):
    freqs = Binomial(n=n, p=p).sample(size=20000, seed=1) / n
    print(f"n={n:>6}: P(|frequency - p| > {eps}) ~ {np.mean(np.abs(freqs - p) > eps):.4f}")

n=      10: relative frequency = 0.3000 (error 0.0000)
n=     100: relative frequency = 0.2300 (error 0.0700)
n=    1000: relative frequency = 0.2910 (error 0.0090)
n=   10000: relative frequency = 0.2972 (error 0.0028)
n= 1000000: relative frequency = 0.3005 (error 0.0005)
n=   100: P(|frequency - p| > 0.02) ~ 0.6656
n=  1000: P(|frequency - p| > 0.02) ~ 0.1731
n= 10000: P(|frequency - p| > 0.02) ~ 0.0000