BrevFeed
Sign in

24. The central limit theorem

🎬 Video · 12 min
💡 Every code box below is live — edit it and hit Run.

The Law of Large Numbers says \bar X_n converges to \mu. The Central Limit Theorem says how — and the answer is the same for essentially every starting distribution.

This is the most consequential theorem in statistics. Almost every confidence interval and hypothesis test in the rest of this course is an application of it.

The statement

Let X_1, X_2, \dots be independent, identically distributed with mean \mu and finite variance \sigma^2. Then

\frac{\bar X_n - \mu}{\sigma/\sqrt n} \;\xrightarrow{d}\; N(0, 1)

Equivalently, for large n:

\bar X_n \;\approx\; N\!\left(\mu,\; \frac{\sigma^2}{n}\right), \qquad \sum_{i=1}^n X_i \;\approx\; N\!\left(n\mu,\; n\sigma^2\right)

Read the standardisation carefully. We subtract the mean and divide by the standard error \sigma/\sqrt n — the exact quantity the LLN produced last lesson. The LLN said that error shrinks; the CLT says what it looks like on the way down.

Why it's remarkable

The starting distribution barely matters. Skewed, discrete, bimodal, uniform — as long as the variance is finite, averages become Normal. The original shape is forgotten.

This explains why the Normal is everywhere in nature. Height, measurement error, and test scores are all sums of many small independent contributions, and any such sum is pushed toward the same bell curve.

It's also why the whole apparatus of statistical inference is possible. We rarely know the distribution of the data, but we don't need to — we only need the distribution of the average, and the CLT hands us that for free.

The proof sketch

We did the work in §1. Standardise so \mu = 0, \sigma = 1, and expand the characteristic function near 0:

\varphi_X(t) = 1 - \frac{t^2}{2} + o(t^2)

For S_n = \frac{1}{\sqrt n}\sum X_j, independence turns the sum into a product:

\varphi_{S_n}(t) = \left[\varphi_X\!\left(\frac{t}{\sqrt n}\right)\right]^n = \left[1 - \frac{t^2}{2n} + o\!\left(\tfrac1n\right)\right]^n \longrightarrow e^{-t^2/2}

which is the standard Normal's characteristic function. Done.

Notice where the assumptions enter. Independence gives the product. Finite variance gives the t^2 term — without it the expansion fails, and so does the theorem.

When it fails, or needs care

Infinite variance. The Cauchy has none, and its averages never become Normal (§1: they stay Cauchy forever). More subtly, distributions with finite variance but very heavy tails converge slowly.

Dependence. Strongly correlated observations break the product step. There are CLTs for weakly dependent sequences, but plain independence is what the basic theorem needs.

"n \ge 30" is a rule of thumb, not a theorem. How fast convergence happens depends on skewness:

  • Symmetric starting distribution: n = 10 is often plenty.
  • Moderately skewed: n \approx 30 is reasonable.
  • Severely skewed (e.g. Exponential, or Bernoulli with p = 0.01): you may need hundreds or thousands.

For a Bernoulli, the usual check is np \ge 10 and n(1-p) \ge 10 — which for p = 0.01 demands n \ge 1000.

It's about the average, not the data. The CLT never claims your observations are Normal. A sample of incomes is right-skewed no matter how large n is; it's the sampling distribution of the mean that becomes Normal. Confusing the two is the most common misreading of the theorem.

Worked example

A die is rolled 100 times. What's the probability the total exceeds 380?

For one die, \mu = 3.5 and \sigma^2 = 35/12 \approx 2.9167 (§1).

For the sum of 100 rolls:

E[S] = 100 \times 3.5 = 350, \qquad \operatorname{Var}(S) = 100 \times 2.9167 = 291.67

\sigma_S = \sqrt{291.67} \approx 17.078

Standardise:

z = \frac{380 - 350}{17.078} \approx 1.757

P(S > 380) \approx 1 - \Phi(1.757) \approx 0.0395

About 4%.

A refinement: the sum is an integer, so "more than 380" means "at least 381". The continuity correction uses 380.5 instead of 380:

z = \frac{380.5 - 350}{17.078} \approx 1.786 \implies P \approx 0.0371

That correction matters when approximating a discrete variable by a continuous one, especially for small n.

Doing it in Python

Watch a badly non-Normal starting distribution become Normal as n grows:

import numpy as np
from scipy.stats import skew, kurtosis

rng = np.random.default_rng(0)
trials = 100_000

print("Starting from Exponential(1) -- heavily right-skewed (skew = 2):\n")
print(f"{'n':>6} {'mean':>8} {'sd':>8} {'sd theory':>11} {'skew':>8} {'ex.kurt':>9}")
for n in (1, 2, 5, 30, 100, 1000):
    means = rng.exponential(1.0, size=(trials, n)).mean(axis=1)
    print(f"{n:>6} {means.mean():>8.4f} {means.std():>8.4f} {1/np.sqrt(n):>11.4f} "
          f"{skew(means):>8.4f} {kurtosis(means):>9.4f}")

print("\nSkewness marches toward 0 -- the shape is becoming Normal.")
print("The sd tracks sigma/sqrt(n) exactly, as the LLN said it would.")

The die example, both with and without the continuity correction:

import numpy as np
from scipy.stats import norm

mu, var = 3.5, 35 / 12
n = 100
mean_s, sd_s = n * mu, np.sqrt(n * var)

plain = 1 - norm.cdf(380, mean_s, sd_s)
corrected = 1 - norm.cdf(380.5, mean_s, sd_s)

rng = np.random.default_rng(1)
sims = rng.integers(1, 7, size=(500_000, n)).sum(axis=1)
exact = (sims > 380).mean()

print(f"CLT, no correction : {plain:.4f}")
print(f"CLT, +0.5 continuity: {corrected:.4f}")
print(f"simulated truth     : {exact:.4f}")

And the failure case — where the CLT simply doesn't apply:

import numpy as np
from scipy.stats import skew

rng = np.random.default_rng(2)
trials = 50_000

print(f"{'n':>7} {'normal-source sd':>18} {'cauchy-source sd':>18}")
for n in (10, 100, 1_000, 10_000):
    from_normal = rng.standard_normal((trials, n)).mean(axis=1)
    from_cauchy = rng.standard_cauchy((trials, n)).mean(axis=1)
    print(f"{n:>7} {from_normal.std():>18.4f} {from_cauchy.std():>18.2f}")

print("\nNormal-source averages shrink like 1/sqrt(n).")
print("Cauchy-source averages do not shrink at all -- infinite variance, no CLT.")

Your turn

1. Heights have \mu = 170, \sigma = 10. For a sample of 25, what's the distribution of \bar X? What's P(\bar X > 174)?

2. A coin is flipped 1,000 times. Approximate P(\text{heads} > 525).

3. Why can't you use the CLT for the sample mean of Cauchy variables?

Solutions

1. By the CLT (and exactly, if heights are Normal to begin with):

\bar X \sim N\!\left(170, \frac{100}{25}\right) = N(170, 4), \qquad \text{SE} = 2

z = \frac{174 - 170}{2} = 2 \implies P(\bar X > 174) = 1 - \Phi(2) \approx 0.0228

About 2.3%. Note how much tighter the mean is than a single observation: an individual above 174 cm has probability 1 - \Phi(0.4) \approx 34\%. Averaging 25 people makes 174 a 2σ event instead of a 0.4σ one.

2. Heads \sim Binomial(1000, 0.5), so \mu = 500 and \sigma = \sqrt{1000 \times 0.25} = 15.81.

With the continuity correction (the count is an integer, and "more than 525" means "at least 526", so use 525.5):

z = \frac{525.5 - 500}{15.81} \approx 1.613 \implies P \approx 0.0534

Roughly 5%. (Without the correction you'd get z = 1.581 and $P \approx 0.0569$ — the exact binomial answer is 0.0534, so the correction genuinely helps.)

3. Because the CLT requires finite variance, and the Cauchy has neither a finite variance nor even a finite mean.

The proof breaks at the Taylor expansion: $\varphi_X(t) = 1 - \frac{t^2}{2} + o(t^2)$ needs E[X^2] < \infty to produce that t^2 coefficient. For the Cauchy, \varphi_X(t) = e^{-|t|}, whose expansion near 0 has a |t| term rather than a t^2 term — not differentiable at the origin.

Running the same argument with $\varphi_{\bar X_n}(t) = \left[e^{-|t|/n}\right]^n = e^{-|t|}$ shows the average has the same distribution as one observation, for every n. There is no convergence to anything, Normal or otherwise.

Heavy tails are not a technicality here — they are the whole obstruction.

Check yourself in code

Show that averages of Exponential(1) samples become Normal: the skewness of the sample mean should fall from 2 toward 0 as n grows, while the standard deviation tracks 1/\sqrt{n}.

Print exactly this:

n=1 skew 2.01 sd 1.0
n=30 skew 0.37 sd 0.18
n=1000 skew 0.05 sd 0.03
skew shrinking: True

Use default_rng(0) and 100000 trials per n. Round the skewness to 2 decimal places and the standard deviation to 2. Report skew shrinking: True if the skewness strictly decreases across the three values of n.

import numpy as np
from scipy.stats import skew

rng = np.random.default_rng(0)
trials = 100_000
skews = []

for n in (1, 30, 1000):
    means = rng.exponential(1.0, size=(trials, n)).mean(axis=1)
    s = round(float(skew(means)), 2)
    skews.append(s)
    print(f"n={n} skew {s} sd {round(float(means.std()), 2)}")

# Report whether the skewness is strictly decreasing.

import numpy as np
from scipy.stats import skew

rng = np.random.default_rng(0)
trials = 100_000
skews = []

for n in (1, 30, 1000):
    means = rng.exponential(1.0, size=(trials, n)).mean(axis=1)
    s = round(float(skew(means)), 2)
    skews.append(s)
    print(f"n={n} skew {s} sd {round(float(means.std()), 2)}")

print("skew shrinking:", all(a > b for a, b in zip(skews, skews[1:])))

Standardised sample means converge in distribution to N(0,1), whatever the underlying distribution, provided the variance is finite. That is what makes inference possible without knowing what you're sampling from — and the finite variance requirement is exactly where it breaks.

Next: how to carry these limits through a function, so you can get the distribution of g(\bar X_n) and not just \bar X_n.

← Previous Next →