Skip to content

Statistics for Data Science

Turning data into decisions: describe what you have, then infer what you don't. Statistics is the discipline that separates patterns from noise.

Table of Contents

Installation / Setup

import numpy as np
import pandas as pd
from scipy import stats

Core Syntax

Two branches of statistics

Branch Question Tools
Descriptive Summarize the sample you have mean, median, IQR, std, plots
Inferential Make claims about the population estimators, CIs, hypothesis tests
Term Definition
Population Everyone/everything of interest
Sample The subset you measure
Parameter Number describing the population ($\mu$, $\sigma$) — usually unknown
Statistic Number computed from the sample ($\bar{x}$, $s$)

Descriptive statistics

Statistic Formula Robust? Use
Mean $\bar{x} = \frac{1}{n}\sum x_i$ No Symmetric data
Median middle value Yes Skewed data (income!)
Mode most frequent — Categorical data
Variance (sample) $s^2 = \frac{1}{n-1}\sum (x_i - \bar{x})^2$ No Spread
Std deviation $s = \sqrt{s^2}$ No Spread in original units
IQR $Q_3 - Q_1$ Yes Robust spread
Range max − min No Quick check
Skewness $E[(X-\mu)^3]/\sigma^3$ — Asymmetry direction
Kurtosis $E[(X-\mu)^4]/\sigma^4$ — Tails/heavy outliers

Sample variance divides by $n-1$ (Bessel's correction) so $E[s^2] = \sigma^2$ — unbiased. NumPy defaults to the population version (ddof=0); pandas/scipy ddof=1 by default.

data.mean(), data.median(), data.std(ddof=1)     # sample std
np.percentile(data, [25, 50, 75])
from scipy import stats
stats.describe(data)                             # n, minmax, mean, var, skew, kurtosis
data.skew(), data.kurt()                         # pandas

Rule of thumb: mean ± 1 std ≈ 68%, ± 2 ≈ 95%, ± 3 ≈ 99.7% for (near-)normal data. For skewed data prefer median + IQR.

Standard error (why samples wobble)

$$\mathrm{SE} = \frac{\sigma}{\sqrt{n}} \quad \text{(use } s \text{ when } \sigma \text{ unknown)}$$

  • The sample mean $\bar{x}$ is a random variable: repeat sampling → different means, spread $\approx$ SE.
  • SE shrinks like $1/\sqrt{n}$ — quadruple the sample to halve the error.
  • SE ≠ std of the data; SE describes the estimator, std describes the data.

Estimators & their quality

MLE intuition: pick the parameter that makes the observed data most likely.

$$\hat{\theta}{MLE} = \arg\max\theta \prod_i p(x_i \mid \theta) = \arg\max_\theta \sum_i \log p(x_i \mid \theta)$$

  • Sample mean is the MLE of a Normal's mean; logistic/softmax training is MLE for classification.
  • Good estimators are unbiased ($E[\hat{\theta}] = \theta$) and low-variance.

Bias–variance tradeoff (the ML core)

$$\underbrace{E[(y - \hat{f}(x))^2]}{\text{expected error}} = \underbrace{\mathrm{Bias}^2}}} + \underbrace{\mathrm{Variance}{\text{overfit}} + \underbrace{\sigma^2}$$}

  • Bias: error from a too-simple model (linear on curved truth).
  • Variance: error from a too-flexible model (memorizing noise).
  • Cross-validation exists to find the sweet spot; regularization trades variance for bias.

Confidence intervals

A 95% CI for the mean (large $n$, or normal data):

$$\bar{x} \pm z_{0.975}\,\frac{s}{\sqrt{n}} \quad \text{or} \quad \bar{x} \pm t_{n-1, 0.975}\,\frac{s}{\sqrt{n}}$$

  • Correct reading: if you repeated the experiment many times, ~95% of such intervals would contain $\mu$.
  • Wrong reading: "$\mu$ is in this interval with 95% probability" ($\mu$ is fixed, not random).
  • Width shrinks with $n$, tightens with higher confidence ($z$ grows), and widens with data spread.
import scipy.stats as stats
mean = x.mean(); se = x.std(ddof=1) / np.sqrt(len(x))
ci = stats.t.interval(0.95, df=len(x) - 1, loc=mean, scale=se)

Correlation vs causation

  • Correlation measures association, not effect. Confounders (lurking variables) can create or hide it.
  • Simpson's paradox: a trend present in every subgroup can reverse when groups are pooled — always check grouping variables.

Common Operations

# Five-number summary (box plot data)
pd.Series(x).describe()                  # count, mean, std, min, quartiles, max

# Standard error of the mean
se = x.std(ddof=1) / np.sqrt(len(x))

# Bootstrap CI (no normality assumption)
rng = np.random.default_rng(0)
boot = np.array([rng.choice(x, len(x), replace=True).mean() for _ in range(10_000)])
ci = np.percentile(boot, [2.5, 97.5])    # percentile interval

# Z-score standardization
z = (x - x.mean()) / x.std(ddof=1)

Practical Examples

1. Compare two groups descriptively before testing

import pandas as pd
grp_a = df[df["variant"] == "A"]["revenue"]
grp_b = df[df["variant"] == "B"]["revenue"]
print(grp_a.describe(), grp_b.describe())
# means differ by 3.1; sd ~15 each; n=500 -> SE ~0.67 each
# difference of ~3 is ~3 SEs away from 0 — promising (formal: hypothesis-testing.md)

2. Bootstrap a confidence interval for the median

boot_med = np.array([
    np.median(rng.choice(income, len(income), replace=True))
    for _ in range(5_000)
])
print(np.percentile(boot_med, [2.5, 97.5]))

3. Spot the bias–variance tradeoff empirically

from sklearn.model_selection import cross_val_score
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures

for degree in [1, 3, 15]:
    Xp = PolynomialFeatures(degree).fit_transform(X)
    scores = cross_val_score(LinearRegression(), Xp, y, cv=5)
    print(degree, scores.mean(), scores.std())
    # degree 1: high bias (low train & test)
    # degree 15: low train error, high test error (variance)

Common Patterns

# "Is this difference real?" quick sanity check
diff = grp_b.mean() - grp_a.mean()
se_diff = np.sqrt(grp_a.var(ddof=1)/len(grp_a) + grp_b.var(ddof=1)/len(grp_b))
print(diff, "±", 1.96 * se_diff)        # interval excludes 0 -> promising

# Outlier policy before statistics: winsorize or log-transform skewed columns
x_log = np.log1p(x)

# Always pair point estimates with uncertainty (CI/SE), not just numbers

Best Practices

  • Report spread with location: mean + std or median + IQR — never mean alone on skewed data.
  • State the sampling story — SE/CI claims require (near-)random sampling.
  • Prefer bootstrap CIs when the sampling distribution is unknown or the statistic is exotic (median, ratio).
  • Plot first, test second — histograms/boxplots catch assumption violations that tests miss.
  • Treat p-values as one input: report effect sizes and CIs too (see hypothesis-testing.md).
  • For regression-style inference use statsmodels; it reports SEs and CIs properly.

Common Mistakes

Mistake Consequence Fix
Confusing SE with data std Overconfident/wrong error bars SE = std/√n; two different things
numpy std() (ddof=0) for sample Slightly biased, tiny effect at scale ddof=1 or pandas .std()
Mean on skewed data Misleading center Median + IQR, log-transform
"95% chance $\mu$ is in this CI" Misinterpreted interval Long-run coverage language
Correlation ⇒ causation Bad decisions Seek experiments or control confounders
Cherry-picking the biggest subgroup effect Simpson's paradox Pool and stratify both
Reporting p-values without effect size Statistical ≠ practical Add Cohen's d / raw difference

Troubleshooting

Symptom Cause Fix
Bootstrap CI far from normal CI Skewed data / small n Trust bootstrap; it doesn't assume normality
Tiny p-value but tiny effect Huge n amplifies differences Report effect size; n drives significance
Std > mean with negative values Data has negatives/wide spread Reconsider distribution; don't force log
SE ~ 0 with n=3 Too few samples More data; treat any CI as unreliable

Official Documentation

Last Updated

2026-09-06 — Initial version.