Data Science with Python

Sampling and Confidence Intervals

Statistics · Lesson 10 of 30

You almost never have data on everyone. You have a sample, and you use it to estimate something about the whole population. Statistics tells you how far to trust that estimate.

Sample statistics vary

Take several samples from the same population and each gives a slightly different mean.

import numpy as np

rng = np.random.default_rng(7)
population = rng.exponential(scale=40, size=200_000)   # skewed, like order values
print("population mean:", round(population.mean(), 2))

for i in range(5):
    sample = rng.choice(population, size=50)
    print("sample", i + 1, "mean:", round(sample.mean(), 2))

The central limit theorem

Now take thousands of samples and look at the distribution of their means.

import numpy as np

rng = np.random.default_rng(7)
population = rng.exponential(scale=40, size=200_000)

for n in (5, 50, 500):
    means = np.array([rng.choice(population, size=n).mean() for _ in range(3000)])
    print(f"n={n:3d}  mean of sample means={means.mean():6.2f}  spread={means.std():5.2f}")

Two things happen as the sample size grows:

  • The sample means centre on the population mean.
  • Their spread shrinks.

The central limit theorem says that the distribution of the sample mean approaches a normal distribution as the sample size increases, even when the population itself is not normal. This is why so many statistical methods rely on the normal distribution.

Standard error

The spread of the sample mean is called the standard error. It is estimated as the sample standard deviation divided by the square root of the sample size.

To halve the standard error you need four times as much data.

Confidence intervals

A confidence interval is a range of plausible values for the population figure.

import numpy as np
from scipy import stats

rng = np.random.default_rng(3)
sample = rng.normal(loc=52, scale=12, size=40)

mean = sample.mean()
sem = stats.sem(sample)
low, high = stats.t.interval(0.95, df=len(sample) - 1, loc=mean, scale=sem)

print(f"sample mean: {mean:.2f}")
print(f"95% confidence interval: {low:.2f} to {high:.2f}")

The interval uses the t distribution, which accounts for the extra uncertainty of estimating the standard deviation from a small sample.

What 95% means

If you repeated the whole procedure many times, about 95% of the intervals would contain the true population mean. You can check that by simulation:

import numpy as np
from scipy import stats

rng = np.random.default_rng(11)
true_mean, hits, trials = 52, 0, 2000

for _ in range(trials):
    s = rng.normal(true_mean, 12, 40)
    low, high = stats.t.interval(0.95, df=39, loc=s.mean(), scale=stats.sem(s))
    hits += low <= true_mean <= high

print("share of intervals containing the true mean:", hits / trials)

Sampling bias

All of this assumes the sample was chosen at random. A large sample that is biased is worse than a small one that is not. A survey answered only by your most engaged customers says nothing reliable about the rest, however many respond.

Practice

Change the sample size in the confidence interval example from 40 to 400. How does the width of the interval change?