import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

rng = np.random.default_rng(6)

NAVY="#0a1230"; INK="#1a2138"; INK_SOFT="#4a5578"
CYAN="#0891b2"; PURPLE="#7c3aed"; AMBER="#d97706"
GREEN="#059669"; PINK="#db2777"; BLUE="#2563eb"; GRID="#e6e9f2"
plt.rcParams.update({
    "figure.facecolor":"white","axes.facecolor":"white","figure.dpi":110,"font.size":11,
    "axes.edgecolor":GRID,"axes.linewidth":1.2,"axes.grid":True,"grid.color":GRID,"axes.axisbelow":True,
    "axes.spines.top":False,"axes.spines.right":False,"axes.titlesize":15,"axes.titleweight":"bold","axes.titlecolor":INK,
    "axes.labelcolor":INK_SOFT,"axes.labelsize":11.5,"xtick.color":INK_SOFT,"ytick.color":INK_SOFT,"legend.frameon":False,
})
def titlecard(ax, title, subtitle=None):
    ax.set_title(title, loc="left", pad=18)
    if subtitle:
        ax.text(0, 1.02, subtitle, transform=ax.transAxes, fontsize=10.5, color=INK_SOFT, va="bottom")
print("✅ Environment ready.")

✅ Environment ready.

# Pretend we COULD measure everyone: a population of 50,000 adult heights (cm)
population = rng.normal(170, 10, 50_000)

# PARAMETERS (truth about the whole population — usually unknown in real life)
mu    = population.mean()                 # μ  population mean
sigma = population.std()                  # σ  population standard deviation
p     = np.mean(population > 180)          # p  proportion taller than 180 cm

# Take ONE random sample of 50 people and compute STATISTICS
sample = rng.choice(population, size=50, replace=False)
xbar   = sample.mean()                    # x̄  sample mean
s      = sample.std(ddof=1)               # s   sample standard deviation
phat   = np.mean(sample > 180)             # p̂  sample proportion

print(pd.DataFrame({
    "quantity":     ["mean", "std dev", "proportion > 180cm"],
    "PARAMETER (pop)":[f"μ = {mu:.1f}", f"σ = {sigma:.1f}", f"p = {p:.2f}"],
    "STATISTIC (n=50)":[f"x̄ = {xbar:.1f}", f"s = {s:.1f}", f"p̂ = {phat:.2f}"],
}).to_string(index=False))

          quantity PARAMETER (pop) STATISTIC (n=50)
              mean       μ = 170.0       x̄ = 166.9
           std dev        σ = 10.0          s = 9.9
proportion > 180cm        p = 0.16        p̂ = 0.08

sample_means = np.array([rng.choice(population, 50, replace=False).mean() for _ in range(1500)])

fig, ax = plt.subplots(figsize=(9,4.6))
ax.hist(sample_means, bins=35, color=PURPLE, alpha=0.85, edgecolor="white")
ax.axvline(mu, color=PINK, lw=2.6, ls="--")
ax.text(mu, ax.get_ylim()[1]*0.93, f"  true μ = {mu:.1f}", color=PINK, fontweight="bold")
titlecard(ax, "1,500 sample means (each from n=50)", "statistics vary — but center on the parameter")
ax.set_xlabel("sample mean  x̄ (cm)"); ax.set_ylabel("how often")
plt.tight_layout(); plt.show()

print(f"Average of all the sample means: {sample_means.mean():.2f}  (≈ μ = {mu:.2f})")
print(f"Spread of the sample means     : {sample_means.std():.2f} cm  (the standard error)")

Average of all the sample means: 170.01  (≈ μ = 170.02)
Spread of the sample means     : 1.40 cm  (the standard error)

sizes = [10, 30, 100, 300, 1000, 3000]
errors = []
for n in sizes:
    means = np.array([rng.choice(population, n, replace=False).mean() for _ in range(400)])
    errors.append(means.std())   # typical distance of x̄ from μ

fig, ax = plt.subplots(figsize=(9,4.6))
ax.plot(sizes, errors, "o-", color=AMBER, lw=2.4, markersize=8, markeredgecolor="white")
for n, e in zip(sizes, errors):
    ax.text(n, e+0.05, f"{e:.2f}", ha="center", fontsize=9, color=INK_SOFT)
titlecard(ax, "Estimation error vs sample size", "standard error falls as n grows (∝ 1/√n)")
ax.set_xscale("log"); ax.set_xlabel("sample size n (log scale)"); ax.set_ylabel("standard error of x̄ (cm)")
plt.tight_layout(); plt.show()

print("Rule of thumb: to halve the error you need ~4× the sample (error ∝ 1/√n).")

Rule of thumb: to halve the error you need ~4× the sample (error ∝ 1/√n).

# A population of 2,000,000 voters; the TRUE support (parameter) p:
true_support = 0.54
voters = rng.random(2_000_000) < true_support
p_param = voters.mean()

# A "census" would poll all 2,000,000. Instead we poll a SAMPLE of 1,000.
poll = rng.choice(voters, size=1000, replace=False)
phat = poll.mean()
moe  = 1.96 * np.sqrt(phat*(1-phat)/1000)   # 95% margin of error

print(f"CENSUS  (all 2,000,000): true support p = {p_param*100:.1f}%  (slow, expensive)")
print(f"SAMPLE  (poll of 1,000): estimate p̂  = {phat*100:.1f}%  ± {moe*100:.1f}%  (fast, cheap)")
print(f"\n95% confidence interval: {(phat-moe)*100:.1f}% to {(phat+moe)*100:.1f}%  -> contains the truth ✔" )

CENSUS  (all 2,000,000): true support p = 54.0%  (slow, expensive)
SAMPLE  (poll of 1,000): estimate p̂  = 53.4%  ± 3.1%  (fast, cheap)

95% confidence interval: 50.3% to 56.5%  -> contains the truth ✔

fig, ax = plt.subplots(figsize=(9,3.4))
ax.axvline(p_param*100, color=PINK, lw=2.6, ls="--", label=f"true support (census) = {p_param*100:.1f}%")
ax.errorbar(phat*100, 0, xerr=moe*100, fmt="o", color=BLUE, markersize=12, capsize=8,
            elinewidth=2.5, markeredgecolor="white", label=f"poll estimate = {phat*100:.1f}% ± {moe*100:.1f}%")
titlecard(ax, "A 1,000-person poll vs the full census", "the cheap estimate brackets the expensive truth")
ax.set_xlabel("% support"); ax.set_yticks([]); ax.set_ylim(-1,1); ax.legend(loc="upper left", fontsize=9)
plt.tight_layout(); plt.show()

#	Demo	Concept
1	Compute parameters, then statistics	Parameter (μ, σ, p) vs statistic (x̄, s, p̂)
2	Take many samples	Sampling variability, statistics jiggle around the parameter
3	Grow the sample	Bigger n → better estimate
4	Run an opinion poll	Census vs sampling + margin of error

🎯 What you'll build in this notebook¶

⚙️ Setup, imports & the book's plotting style¶