Confidence intervals

An estimate without an error bar is half an answer. A confidence interval says what range the true parameter plausibly lies in, and it is also how traders are asked to quote: “give me a 90% interval for the number of windows in this building.”

Definition

A confidence interval of level $1-\alpha$ for a parameter $\theta$ is a random interval $\mathcal I$, computed from the data, with

Careful with the wording. The interval is random; the parameter is fixed. “95% confidence” is a property of the procedure: if you repeated the experiment many times, about 95% of the intervals you built would contain $\theta$. A single interval either does or doesn’t.

Building one from the CLT

The estimator is usually an average, so its centre comes from the LLN and its width from the CLT. With $q_{\alpha/2}$ the $(1-\alpha/2)$ quantile of $\mathcal N(0,1)$ (1.645 for 90%, 1.96 for 95%, 2.576 for 99%):

There is a catch: $\sigma$ usually depends on the unknown parameter. For a Bernoulli, $\sigma = \sqrt{p(1-p)}$. Two standard fixes:

• Conservative bound. $p(1-p) \le \tfrac14$ always, so $\hat p \pm \dfrac{q_{\alpha/2}}{2\sqrt n}$ is valid whatever $p$ is. It is a little too wide.
• Plug-in (Slutsky’s theorem). Replace $p$ by $\hat p$ inside the variance: $\hat p \pm q_{\alpha/2}\sqrt{\hat p(1-\hat p)/n}$. Slutsky says the CLT still holds asymptotically after this swap.

Intervals on a trading floor

In interviews, an interval is your honesty about uncertainty. Ask for a 95% interval on a normal quantity with mean 100 and sd 10 and the answer is about 80 to 120 (±2σ). Ask for one on “the number of windows in this building” and there is no formula: your interval should be wider when you know less, and still realistic. “Between 0 and a billion” contains the truth and is useless, exactly like a market that is 0 @ 1 billion (lesson 17).

In code

import numpy as np
from scipy import stats
rng = np.random.default_rng(5)

p, n, level = 0.55, 50, 0.95
q = stats.norm.ppf(1 - (1 - level) / 2)          # 1.96

hits = 0
for _ in range(10_000):
    ph = rng.binomial(n, p) / n
    half = q * np.sqrt(ph * (1 - ph) / n)        # plug-in (Slutsky)
    hits += (ph - half <= p <= ph + half)
print(hits / 10_000)                              # ≈ 0.94: slightly under 95% at n = 50