Hands-on lab — Day 138: Data Ethics, Bias, and Provenance

"""Shared fixtures for the Day 138 lab.

All four fixtures are built from a fixed seed, so every test in the suite
sees exactly the same synthetic table and every assertion below is a claim
about a specific, reproducible set of numbers rather than about a
distribution "on average".
"""

from __future__ import annotations

import pytest

import ethics as et

SEED = 138


@pytest.fixture(scope="session")
def register():
    """The synthetic register: quasi-identifiers plus one sensitive column."""
    return et.synthetic_register(n=5_000, seed=SEED)


@pytest.fixture(scope="session")
def generalised(register):
    """The same register with birth year coarsened to a decade band."""
    return et.generalise_quasi_ids(register)


@pytest.fixture(scope="session")
def fairness_pop():
    """The integer-exact, perfectly calibrated two-group population."""
    return et.fairness_population()


@pytest.fixture(scope="session")
def versions():
    """Two releases of one dataset, differing only in group composition."""
    return et.build_versions(n=2_000, seed=SEED)

"""Bias you can measure — the measurement and documentation functions.

Every function here is deterministic. The ones that draw data take an
explicit seed and use ``numpy.random.default_rng``, so the same seed gives
the same numbers on every machine and every run. Nothing in this module
reads a real dataset about real people: every table is constructed here,
in code, with properties chosen so that the thing being demonstrated is a
fact about the construction rather than a claim about anyone.

The module has four groups of functions:

* measurement of bias  -- ``biased_frame_experiment``, ``coverage_report``,
  ``proxy_gap_experiment``, ``aggregation_bias_experiment``
* fairness accounting  -- ``fairness_population``, ``evaluate_policy``,
  ``calibration_by_score_bin``, ``fairness_incompatibility``
* disclosure risk      -- ``unique_row_count``, ``generalise_quasi_ids``,
  ``k_anonymity_level``, ``suppress_small_classes``, ``homogeneous_classes``
* documentation        -- ``check_datasheet``, ``diff_versions``
"""

from __future__ import annotations

from dataclasses import dataclass, field
from typing import Any, Iterable, Mapping, Sequence

import numpy as np
import pandas as pd

# ---------------------------------------------------------------------------
# 1. More data does not fix a biased frame
# ---------------------------------------------------------------------------

#: The world the sample is drawn from. Group "B" is one tenth of the
#: population and its outcome sits INTERCEPT_GAP units above group "A" at
#: every value of x. The slope is identical for both groups, so the only
#: thing a pooled straight-line fit can get wrong about group B is the
#: intercept -- which makes the bias analytically predictable and therefore
#: checkable rather than merely plausible.
POPULATION_SHARE_B = 0.10
FRAME_SHARE_B = 0.01
INTERCEPT_A = 2.0
INTERCEPT_B = 8.0
INTERCEPT_GAP = INTERCEPT_B - INTERCEPT_A
TRUE_SLOPE = 1.0
NOISE_SD = 1.0


def draw_from_frame(n: int, rng: np.random.Generator) -> pd.DataFrame:
    """Draw ``n`` rows the way the biased sampling frame draws them.

    The frame reaches group B with probability ``FRAME_SHARE_B`` rather than
    its true population share ``POPULATION_SHARE_B``. Everything else about
    the two groups is identical: same x distribution, same slope, same noise.
    """
    group = np.where(rng.random(n) < FRAME_SHARE_B, "B", "A")
    x = rng.uniform(0.0, 10.0, size=n)
    intercept = np.where(group == "B", INTERCEPT_B, INTERCEPT_A)
    y = intercept + TRUE_SLOPE * x + rng.normal(0.0, NOISE_SD, size=n)
    return pd.DataFrame({"group": group, "x": x, "y": y})


def fit_pooled_line(frame: pd.DataFrame) -> tuple[float, float]:
    """Ordinary least squares of y on x, ignoring the group column entirely.

    Returns ``(slope, intercept)``. This is the model almost everyone fits
    first: one line for everybody.
    """
    slope, intercept = np.polyfit(frame["x"].to_numpy(), frame["y"].to_numpy(), 1)
    return float(slope), float(intercept)


def signed_group_error(slope: float, intercept: float, group: str) -> float:
    """Signed prediction error of a pooled line for one group's true line.

    Because both groups share ``TRUE_SLOPE``, the error is the same at every
    x, so a single number describes it exactly: predicted intercept minus
    true intercept, plus the slope error times a representative x. Evaluated
    at x = 5, the midpoint of the x range.
    """
    true_intercept = INTERCEPT_B if group == "B" else INTERCEPT_A
    x_mid = 5.0
    predicted = intercept + slope * x_mid
    truth = true_intercept + TRUE_SLOPE * x_mid
    return float(predicted - truth)


def biased_frame_experiment(
    n: int, replicates: int = 40, seed: int = 138
) -> dict[str, float]:
    """Fit the pooled line ``replicates`` times at sample size ``n``.

    Returns the mean signed error for each group across replicates (the
    BIAS -- what more data cannot fix) and the standard deviation of that
    error across replicates (the VARIANCE -- what more data does fix), plus
    the in-sample RMSE, which is the number that makes the model look fine.
    """
    root = np.random.default_rng(seed)
    seeds = root.integers(0, 2**32 - 1, size=replicates)

    errors_a: list[float] = []
    errors_b: list[float] = []
    rmses: list[float] = []
    b_shares: list[float] = []

    for s in seeds:
        rng = np.random.default_rng(int(s))
        sample = draw_from_frame(n, rng)
        slope, intercept = fit_pooled_line(sample)
        errors_a.append(signed_group_error(slope, intercept, "A"))
        errors_b.append(signed_group_error(slope, intercept, "B"))
        predicted = intercept + slope * sample["x"].to_numpy()
        residual = sample["y"].to_numpy() - predicted
        rmses.append(float(np.sqrt(np.mean(residual**2))))
        b_shares.append(float((sample["group"] == "B").mean()))

    return {
        "n": float(n),
        "replicates": float(replicates),
        "bias_a": float(np.mean(errors_a)),
        "bias_b": float(np.mean(errors_b)),
        "sd_b": float(np.std(errors_b, ddof=1)),
        "sd_a": float(np.std(errors_a, ddof=1)),
        "in_sample_rmse": float(np.mean(rmses)),
        "sample_share_b": float(np.mean(b_shares)),
    }


def bias_variance_ladder(
    sizes: Sequence[int] = (500, 5_000, 50_000),
    replicates: int = 40,
    seed: int = 138,
) -> pd.DataFrame:
    """Run ``biased_frame_experiment`` at each sample size and tabulate it."""
    rows = [biased_frame_experiment(n, replicates=replicates, seed=seed) for n in sizes]
    table = pd.DataFrame(rows)
    table["n"] = table["n"].astype("int64")
    return table


# ---------------------------------------------------------------------------
# 2. Coverage mismatch, computed
# ---------------------------------------------------------------------------


def coverage_report(
    sample_counts: Mapping[str, int],
    reference_shares: Mapping[str, float],
    flag_below: float = 0.5,
    flag_above: float = 2.0,
) -> dict[str, Any]:
    """Compare a sample's group composition against a reference population.

    ``total_variation_distance`` is one half the sum of absolute differences
    between the two share vectors: 0.0 when the sample mirrors the reference,
    1.0 when they share no mass at all. It is a single summary, so it is
    reported alongside the thing that actually names the problem: a
    per-group ``representation_ratio`` of sample share divided by reference
    share. A ratio of 0.1 means the group appears at one tenth of the rate
    the reference population says it should.

    A group is flagged when its ratio falls outside
    ``[flag_below, flag_above]``.
    """
    groups = sorted(set(sample_counts) | set(reference_shares))
    total = sum(sample_counts.values())
    if total <= 0:
        raise ValueError("sample_counts must contain at least one row")

    sample_shares = {g: sample_counts.get(g, 0) / total for g in groups}
    ratios: dict[str, float] = {}
    flagged: list[str] = []
    for g in groups:
        reference = reference_shares.get(g, 0.0)
        if reference <= 0.0:
            ratios[g] = float("inf") if sample_shares[g] > 0 else 1.0
        else:
            ratios[g] = sample_shares[g] / reference
        if ratios[g] < flag_below or ratios[g] > flag_above:
            flagged.append(g)

    tvd = 0.5 * sum(
        abs(sample_shares[g] - reference_shares.get(g, 0.0)) for g in groups
    )

    return {
        "groups": groups,
        "sample_shares": sample_shares,
        "reference_shares": dict(reference_shares),
        "representation_ratio": ratios,
        "total_variation_distance": float(tvd),
        "flagged": flagged,
        "worst_group": min(ratios, key=lambda g: ratios[g]),
    }


# ---------------------------------------------------------------------------
# 3. The proxy gap
# ---------------------------------------------------------------------------

#: Group B's need is recorded at ACCESS_FACTOR_B of its true level, because
#: the proxy measures a recorded interaction rather than the need itself and
#: group B interacts with the recording system less often. Group A's need is
#: recorded faithfully.
ACCESS_FACTOR_A = 1.0
ACCESS_FACTOR_B = 0.85


def proxy_gap_experiment(
    n_a: int = 3_000, n_b: int = 1_000, budget: int = 500, seed: int = 138
) -> dict[str, Any]:
    """Rank by a proxy, then check what that costs against the real target.

    Both groups are drawn from the SAME true-need distribution, so a fair
    procedure would select them in proportion to their size. The proxy
    understates group B's need systematically. Selecting the top ``budget``
    rows by proxy is compared against selecting the top ``budget`` by true
    need, and both the composition and the total true need served are
    reported for each.
    """
    rng = np.random.default_rng(seed)
    n = n_a + n_b
    group = np.array(["A"] * n_a + ["B"] * n_b)
    true_need = rng.normal(50.0, 10.0, size=n)
    access = np.where(group == "B", ACCESS_FACTOR_B, ACCESS_FACTOR_A)
    proxy = true_need * access + rng.normal(0.0, 2.0, size=n)

    frame = pd.DataFrame({"group": group, "true_need": true_need, "proxy": proxy})

    by_proxy = frame.nlargest(budget, "proxy")
    by_target = frame.nlargest(budget, "true_need")

    population_share_b = n_b / n
    return {
        "budget": budget,
        "population_share_b": float(population_share_b),
        "proxy_correlation": float(frame["proxy"].corr(frame["true_need"])),
        "selected_share_b_by_proxy": float((by_proxy["group"] == "B").mean()),
        "selected_share_b_by_target": float((by_target["group"] == "B").mean()),
        "need_served_by_proxy": float(by_proxy["true_need"].sum()),
        "need_served_by_target": float(by_target["true_need"].sum()),
        "b_need_served_by_proxy": float(
            by_proxy.loc[by_proxy["group"] == "B", "true_need"].sum()
        ),
        "b_need_served_by_target": float(
            by_target.loc[by_target["group"] == "B", "true_need"].sum()
        ),
    }


# ---------------------------------------------------------------------------
# 4. Aggregation bias
# ---------------------------------------------------------------------------


def aggregation_bias_experiment(
    n_per_group: int = 2_000, seed: int = 138
) -> dict[str, Any]:
    """One pooled line against two per-group lines on the same rows.

    Group A lives at low x with a high intercept; group B lives at high x
    with a much higher intercept. Both groups have a NEGATIVE true slope.
    Pooled, the group offset dominates and the fitted slope comes out
    positive: a single model that is wrong in direction for every subgroup
    in the data, while looking like a confident fit overall.
    """
    rng = np.random.default_rng(seed)
    x_a = rng.uniform(0.0, 5.0, size=n_per_group)
    y_a = 10.0 - 1.0 * x_a + rng.normal(0.0, 1.0, size=n_per_group)
    x_b = rng.uniform(5.0, 10.0, size=n_per_group)
    y_b = 30.0 - 1.0 * x_b + rng.normal(0.0, 1.0, size=n_per_group)

    frame = pd.DataFrame(
        {
            "group": np.array(["A"] * n_per_group + ["B"] * n_per_group),
            "x": np.concatenate([x_a, x_b]),
            "y": np.concatenate([y_a, y_b]),
        }
    )

    pooled_slope, pooled_intercept = fit_pooled_line(frame)

    per_group: dict[str, dict[str, float]] = {}
    pooled_rmse: dict[str, float] = {}
    grouped_rmse: dict[str, float] = {}
    for name, part in frame.groupby("group", observed=True):
        slope, intercept = np.polyfit(part["x"].to_numpy(), part["y"].to_numpy(), 1)
        per_group[str(name)] = {"slope": float(slope), "intercept": float(intercept)}
        own = intercept + slope * part["x"].to_numpy()
        pool = pooled_intercept + pooled_slope * part["x"].to_numpy()
        grouped_rmse[str(name)] = float(
            np.sqrt(np.mean((part["y"].to_numpy() - own) ** 2))
        )
        pooled_rmse[str(name)] = float(
            np.sqrt(np.mean((part["y"].to_numpy() - pool) ** 2))
        )

    return {
        "pooled_slope": pooled_slope,
        "pooled_intercept": pooled_intercept,
        "per_group_fit": per_group,
        "pooled_rmse_by_group": pooled_rmse,
        "per_group_rmse_by_group": grouped_rmse,
        "pooled_worse_for_every_group": all(
            pooled_rmse[g] > grouped_rmse[g] for g in pooled_rmse
        ),
        "sign_flip": pooled_slope > 0
        and all(fit["slope"] < 0 for fit in per_group.values()),
    }


# ---------------------------------------------------------------------------
# 5. The fairness tension
# ---------------------------------------------------------------------------

#: A fully specified, integer-exact population. Each entry is
#: (group, score, count). Within every (group, score) cell exactly
#: ``score`` of the rows are positive, so the score is PERFECTLY calibrated
#: for both groups by construction -- no fitting, no randomness, no
#: approximation. The two groups differ only in where their mass sits, which
#: gives them different base rates: that difference is the whole engine of
#: the incompatibility.
FAIRNESS_CELLS: tuple[tuple[str, float, int], ...] = (
    ("A", 0.1, 50),
    ("A", 0.3, 100),
    ("A", 0.5, 200),
    ("A", 0.7, 300),
    ("A", 0.9, 350),
    ("B", 0.1, 350),
    ("B", 0.3, 300),
    ("B", 0.5, 200),
    ("B", 0.7, 100),
    ("B", 0.9, 50),
)


def fairness_population() -> pd.DataFrame:
    """The constructed population as a tidy frame of score cells.

    Columns: group, score, count, positives. ``positives`` is exactly
    ``score * count`` and every one of those products is a whole number, so
    nothing here is rounded.
    """
    rows = []
    for group, score, count in FAIRNESS_CELLS:
        positives = score * count
        assert abs(positives - round(positives)) < 1e-9, "cell must be integer-exact"
        rows.append(
            {
                "group": group,
                "score": score,
                "count": count,
                "positives": int(round(positives)),
            }
        )
    return pd.DataFrame(rows)


def base_rates(population: pd.DataFrame | None = None) -> dict[str, float]:
    """Fraction of each group that is actually positive."""
    pop = fairness_population() if population is None else population
    grouped = pop.groupby("group", observed=True)[["count", "positives"]].sum()
    return {
        str(g): float(row["positives"] / row["count"])
        for g, row in grouped.iterrows()
    }


def calibration_by_score_bin(population: pd.DataFrame | None = None) -> pd.DataFrame:
    """Observed positive rate in every (group, score) cell.

    Calibration means the score can be read as a probability: among rows
    scored 0.7, seven in ten really are positive, in EVERY group. This
    returns the observed rate so the claim can be checked rather than
    asserted.
    """
    pop = fairness_population() if population is None else population
    out = pop.copy()
    out["observed_rate"] = out["positives"] / out["count"]
    out["deviation"] = (out["observed_rate"] - out["score"]).abs()
    return out[["group", "score", "count", "observed_rate", "deviation"]]


def _select_from_top(part: pd.DataFrame, target: float, key: str) -> float:
    """Select ``target`` units of ``key`` from the highest scores downward.

    Returns the number of ROWS selected. When the target lands inside a score
    band, that band is taken proportionally -- the standard way to hit an
    exact rate when a whole band would overshoot. Returns a float because a
    proportional slice of a band generally is not a whole number of rows.
    """
    ordered = part.sort_values("score", ascending=False)
    remaining = target
    rows_selected = 0.0
    for _, band in ordered.iterrows():
        available = float(band[key])
        if available <= 0:
            continue
        if remaining >= available - 1e-9:
            rows_selected += float(band["count"])
            remaining -= available
        else:
            fraction = remaining / available
            rows_selected += fraction * float(band["count"])
            remaining = 0.0
        if remaining <= 1e-9:
            break
    return rows_selected


def _metrics_for_selection(part: pd.DataFrame, rows_selected: float) -> dict[str, float]:
    """Selection rate, true-positive rate and precision for one group.

    ``rows_selected`` is a count of rows taken from the top score band
    downward; the positives inside that selection are accumulated band by
    band, taking a proportional slice of whichever band the cut falls in.
    """
    ordered = part.sort_values("score", ascending=False)
    remaining = rows_selected
    positives_selected = 0.0
    for _, band in ordered.iterrows():
        count = float(band["count"])
        if remaining >= count - 1e-9:
            positives_selected += float(band["positives"])
            remaining -= count
        else:
            fraction = remaining / count if count else 0.0
            positives_selected += fraction * float(band["positives"])
            remaining = 0.0
        if remaining <= 1e-9:
            break

    total = float(part["count"].sum())
    total_positives = float(part["positives"].sum())
    return {
        "selection_rate": rows_selected / total,
        "true_positive_rate": positives_selected / total_positives,
        "precision": positives_selected / rows_selected if rows_selected else 0.0,
        "rows_selected": rows_selected,
        "positives_selected": positives_selected,
    }


def evaluate_policy(
    policy: str, population: pd.DataFrame | None = None
) -> dict[str, Any]:
    """Score one decision policy on all three fairness criteria at once.

    Policies:

    * ``"single_threshold"`` -- one cut at score >= 0.5 for everybody. This
      is the policy that leaves the calibrated score untouched.
    * ``"demographic_parity"`` -- give both groups the same selection rate,
      set to the rate the single threshold produces overall.
    * ``"equal_opportunity"`` -- give both groups the same true-positive
      rate, set to the rate that policy targets.

    Reported gaps are absolute differences between the two groups:
    ``parity_gap`` (selection rate), ``equal_opportunity_gap`` (true-positive
    rate) and ``precision_gap`` (precision among the selected, the
    calibration-of-the-decision criterion). None of the three is privileged
    here; all three are printed so the trade can be seen rather than
    argued about.
    """
    pop = fairness_population() if population is None else population
    groups = sorted(pop["group"].unique())
    parts = {g: pop[pop["group"] == g] for g in groups}

    if policy == "single_threshold":
        rows = {
            g: float(parts[g].loc[parts[g]["score"] >= 0.5, "count"].sum())
            for g in groups
        }
    elif policy == "demographic_parity":
        selected_at_threshold = sum(
            float(parts[g].loc[parts[g]["score"] >= 0.5, "count"].sum())
            for g in groups
        )
        total = float(pop["count"].sum())
        rate = selected_at_threshold / total
        rows = {g: rate * float(parts[g]["count"].sum()) for g in groups}
    elif policy == "equal_opportunity":
        target_tpr = 0.80
        rows = {
            g: _select_from_top(
                parts[g], target_tpr * float(parts[g]["positives"].sum()), "positives"
            )
            for g in groups
        }
    else:  # pragma: no cover - guarded by the caller
        raise ValueError(f"unknown policy: {policy!r}")

    per_group = {g: _metrics_for_selection(parts[g], rows[g]) for g in groups}
    a, b = groups[0], groups[1]
    return {
        "policy": policy,
        "per_group": per_group,
        "parity_gap": abs(
            per_group[a]["selection_rate"] - per_group[b]["selection_rate"]
        ),
        "equal_opportunity_gap": abs(
            per_group[a]["true_positive_rate"] - per_group[b]["true_positive_rate"]
        ),
        "precision_gap": abs(per_group[a]["precision"] - per_group[b]["precision"]),
    }


FAIRNESS_POLICIES = ("single_threshold", "demographic_parity", "equal_opportunity")


def fairness_incompatibility(tolerance: float = 1e-9) -> dict[str, Any]:
    """Evaluate every policy and report which criteria each one satisfies.

    Returns the three policies' gap triples plus ``any_policy_satisfies_all``,
    which is the claim this exercise exists to test. It is False here, and
    that is a consequence of the two groups having different base rates --
    not of the policies being badly chosen.
    """
    results = {p: evaluate_policy(p) for p in FAIRNESS_POLICIES}
    satisfied = {
        p: {
            "demographic_parity": r["parity_gap"] <= tolerance,
            "equal_opportunity": r["equal_opportunity_gap"] <= tolerance,
            "equal_precision": r["precision_gap"] <= tolerance,
        }
        for p, r in results.items()
    }
    return {
        "results": results,
        "satisfied": satisfied,
        "base_rates": base_rates(),
        "any_policy_satisfies_all": any(
            all(flags.values()) for flags in satisfied.values()
        ),
    }


# ---------------------------------------------------------------------------
# 6-7. Disclosure risk: uniqueness, generalisation, k-anonymity and its limit
# ---------------------------------------------------------------------------

QUASI_IDS = ("birth_year", "postcode", "sex")
GENERALISED_QUASI_IDS = ("birth_decade", "postcode", "sex")


def synthetic_register(n: int = 5_000, seed: int = 138) -> pd.DataFrame:
    """A synthetic register: no real person is represented here.

    Three quasi-identifiers -- birth year, postcode, sex -- plus one
    sensitive attribute, ``diagnosis``. Nothing in this table was collected;
    every value is drawn from a seeded generator, which is the point: the
    uniqueness counts below are properties of the shape of the table, not
    facts about anybody.
    """
    rng = np.random.default_rng(seed)
    return pd.DataFrame(
        {
            "birth_year": rng.integers(1940, 2006, size=n),
            "postcode": rng.integers(1000, 1060, size=n).astype(str),
            "sex": rng.choice(["F", "M"], size=n),
            "diagnosis": rng.choice(
                ["none", "asthma", "diabetes", "hypertension"],
                size=n,
                p=[0.55, 0.15, 0.15, 0.15],
            ),
        }
    )


def unique_row_count(frame: pd.DataFrame, quasi_ids: Iterable[str]) -> int:
    """Rows that are the ONLY row with their combination of quasi-identifiers.

    A unique row can be singled out by anyone who knows those few fields
    about a person, whether or not the table carries a name.
    """
    keys = list(quasi_ids)
    sizes = frame.groupby(keys, observed=True)[keys[0]].transform("size")
    return int((sizes == 1).sum())


def generalise_quasi_ids(frame: pd.DataFrame) -> pd.DataFrame:
    """Replace exact birth year with the decade it falls in.

    One field, coarsened one step. Nothing else changes -- so any drop in
    uniqueness is attributable to this single edit.
    """
    out = frame.copy()
    out["birth_decade"] = (out["birth_year"] // 10 * 10).astype("int64")
    return out


def k_anonymity_level(frame: pd.DataFrame, quasi_ids: Iterable[str]) -> int:
    """The size of the smallest equivalence class.

    A table is k-anonymous when every combination of quasi-identifiers that
    appears at all appears at least k times, so this minimum IS the k the
    table achieves.
    """
    keys = list(quasi_ids)
    return int(frame.groupby(keys, observed=True)[keys[0]].size().min())


def suppress_small_classes(
    frame: pd.DataFrame, quasi_ids: Iterable[str], k: int
) -> tuple[pd.DataFrame, int]:
    """Drop every row whose equivalence class is smaller than ``k``.

    Returns the surviving frame and the number of rows suppressed. This is
    what k-anonymity costs: it is bought with rows, and the rows it costs
    are exactly the unusual ones -- which is worth saying out loud, because
    the people in rare categories are frequently the people an analysis was
    supposed to be about.
    """
    keys = list(quasi_ids)
    sizes = frame.groupby(keys, observed=True)[keys[0]].transform("size")
    kept = frame[sizes >= k]
    return kept.reset_index(drop=True), int(len(frame) - len(kept))


#: A small, hand-built table where k-anonymity holds and still discloses.
#: Every equivalence class has four rows, so the table is 4-anonymous. The
#: class (1970, "1001", "F") has four rows that all carry the same
#: diagnosis, so knowing only that someone is a 1970s-born woman in postcode
#: 1001 and appears in this table reveals her diagnosis exactly. No
#: re-identification was needed; the anonymity held and the secret still
#: came out.
HOMOGENEOUS_TABLE_ROWS: tuple[dict[str, Any], ...] = tuple(
    [
        {"birth_decade": 1970, "postcode": "1001", "sex": "F", "diagnosis": "diabetes"}
    ]
    * 4
    + [
        {"birth_decade": 1970, "postcode": "1001", "sex": "M", "diagnosis": "none"},
        {"birth_decade": 1970, "postcode": "1001", "sex": "M", "diagnosis": "asthma"},
        {"birth_decade": 1970, "postcode": "1001", "sex": "M", "diagnosis": "none"},
        {
            "birth_decade": 1970,
            "postcode": "1001",
            "sex": "M",
            "diagnosis": "hypertension",
        },
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "none"},
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "asthma"},
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "none"},
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "diabetes"},
    ]
)


def homogeneous_table() -> pd.DataFrame:
    """The 4-anonymous table that still leaks, as a frame."""
    return pd.DataFrame(list(HOMOGENEOUS_TABLE_ROWS))


def homogeneous_classes(
    frame: pd.DataFrame, quasi_ids: Iterable[str], sensitive: str, k: int
) -> list[dict[str, Any]]:
    """Equivalence classes of size >= k whose sensitive values are all identical.

    These are the classes where k-anonymity is satisfied and the sensitive
    attribute is disclosed anyway. In the l-diversity vocabulary these are
    the classes with l = 1.
    """
    keys = list(quasi_ids)
    leaks: list[dict[str, Any]] = []
    for key, part in frame.groupby(keys, observed=True):
        values = set(part[sensitive])
        if len(part) >= k and len(values) == 1:
            leaks.append(
                {
                    "class": tuple(key) if isinstance(key, tuple) else (key,),
                    "size": int(len(part)),
                    "distinct_sensitive_values": len(values),
                    "disclosed_value": next(iter(values)),
                }
            )
    return sorted(leaks, key=lambda item: str(item["class"]))


# ---------------------------------------------------------------------------
# 8. A datasheet contract
# ---------------------------------------------------------------------------

#: The provenance fields a dataset must carry before anything is built on
#: it. Each one answers a question that cannot be recovered from the data
#: once it is lost: who made these decisions, when, for what, and about whom.
DATASHEET_REQUIRED_FIELDS: tuple[str, ...] = (
    "collector",
    "collection_period",
    "purpose",
    "population_definition",
    "sampling_frame",
    "inclusion_criteria",
    "exclusion_criteria",
    "known_gaps",
    "licence",
    "version",
    "changelog",
)


def check_datasheet(
    record: Mapping[str, Any], required: Sequence[str] = DATASHEET_REQUIRED_FIELDS
) -> dict[str, Any]:
    """Check a dataset's documentation against the required provenance fields.

    A field counts as missing when it is absent, ``None``, or an empty
    string or collection -- an empty ``known_gaps`` list is a real answer
    only if it is deliberate, so this contract requires it to be written
    down rather than left off. Returns the missing fields BY NAME, because
    "documentation incomplete" is not actionable and
    "no exclusion_criteria recorded" is.
    """
    missing = []
    for name in required:
        if name not in record:
            missing.append(name)
            continue
        value = record[name]
        if value is None:
            missing.append(name)
        elif isinstance(value, str) and not value.strip():
            missing.append(name)
        elif isinstance(value, (list, tuple, dict, set)) and len(value) == 0:
            missing.append(name)
    return {
        "complete": not missing,
        "missing": missing,
        "checked": list(required),
    }


# ---------------------------------------------------------------------------
# 9. Version drift
# ---------------------------------------------------------------------------


@dataclass(frozen=True)
class DatasetVersion:
    """One released version of a dataset: the rows plus the datasheet."""

    version: str
    frame: pd.DataFrame
    datasheet: dict[str, Any] = field(default_factory=dict)


def build_versions(n: int = 2_000, seed: int = 138) -> tuple[DatasetVersion, DatasetVersion]:
    """Two releases of the same dataset, differing only in who is in it.

    The measured column is byte-for-byte identical between the two versions.
    What changed is the sampling frame: version 2 cut group B's quota from
    half the sample to a quarter. Because both groups' values are drawn from
    the same distribution, every summary statistic of ``value`` -- count,
    mean, standard deviation, quartiles, extremes -- is EXACTLY unchanged. A
    reader comparing ``describe()`` output between the two releases sees
    nothing at all.
    """
    rng = np.random.default_rng(seed)
    value = rng.normal(50.0, 10.0, size=n)

    labels_v1 = np.where(np.arange(n) % 2 == 0, "A", "B")
    labels_v2 = np.where(np.arange(n) % 4 == 3, "B", "A")

    common = {
        "collector": "Course fixture generator (synthetic)",
        "collection_period": "2026-01-01/2026-06-30",
        "purpose": "Teaching example for measuring composition drift",
        "population_definition": "Synthetic adults in groups A and B",
        "inclusion_criteria": ["synthetic record generated by this module"],
        "exclusion_criteria": ["none"],
        "licence": "CC0-1.0",
        "known_gaps": ["contains no real people by design"],
    }

    v1 = DatasetVersion(
        version="1.0.0",
        frame=pd.DataFrame({"group": labels_v1, "value": value}),
        datasheet={
            **common,
            "sampling_frame": "equal quota: 50% group A, 50% group B",
            "version": "1.0.0",
            "changelog": ["1.0.0 initial release"],
        },
    )
    v2 = DatasetVersion(
        version="2.0.0",
        frame=pd.DataFrame({"group": labels_v2, "value": value}),
        datasheet={
            **common,
            "sampling_frame": "revised quota: 75% group A, 25% group B",
            "version": "2.0.0",
            "changelog": [
                "1.0.0 initial release",
                "2.0.0 sampling frame re-scoped: group B quota reduced from 50% to 25%",
            ],
        },
    )
    return v1, v2


def summary_stats(frame: pd.DataFrame, column: str = "value") -> dict[str, float]:
    """The summary a reader would actually look at when comparing releases."""
    series = frame[column]
    return {
        "count": float(series.count()),
        "mean": float(series.mean()),
        "std": float(series.std()),
        "min": float(series.min()),
        "median": float(series.median()),
        "max": float(series.max()),
    }


def diff_versions(
    old: DatasetVersion,
    new: DatasetVersion,
    group_column: str = "group",
    value_column: str = "value",
    shift_threshold: float = 0.05,
) -> dict[str, Any]:
    """Compare two releases on summary statistics AND on composition.

    ``summary_identical`` is True when every summary statistic matches to
    within floating-point tolerance -- the case where a release-note-free
    upgrade looks like no change at all. ``composition_shift`` is the largest
    absolute change in any group's share. ``changelog_explains`` reports
    whether the new version's changelog gained an entry the old one did not
    have, which is the only part of this comparison that says WHY.
    """
    old_stats = summary_stats(old.frame, value_column)
    new_stats = summary_stats(new.frame, value_column)
    summary_identical = all(
        abs(old_stats[k] - new_stats[k]) < 1e-9 for k in old_stats
    )

    old_shares = old.frame[group_column].value_counts(normalize=True).to_dict()
    new_shares = new.frame[group_column].value_counts(normalize=True).to_dict()
    groups = sorted(set(old_shares) | set(new_shares))
    shifts = {
        g: float(new_shares.get(g, 0.0) - old_shares.get(g, 0.0)) for g in groups
    }
    composition_shift = max(abs(v) for v in shifts.values()) if shifts else 0.0

    old_log = list(old.datasheet.get("changelog", []))
    new_log = list(new.datasheet.get("changelog", []))
    new_entries = [entry for entry in new_log if entry not in old_log]

    return {
        "old_version": old.version,
        "new_version": new.version,
        "summary_old": old_stats,
        "summary_new": new_stats,
        "summary_identical": summary_identical,
        "group_shares_old": {g: float(old_shares.get(g, 0.0)) for g in groups},
        "group_shares_new": {g: float(new_shares.get(g, 0.0)) for g in groups},
        "group_share_shift": shifts,
        "composition_shift": float(composition_shift),
        "material_shift": composition_shift >= shift_threshold,
        "summary_hides_the_change": summary_identical
        and composition_shift >= shift_threshold,
        "changelog_new_entries": new_entries,
        "changelog_explains": bool(new_entries),
        "frame_changed": old.datasheet.get("sampling_frame")
        != new.datasheet.get("sampling_frame"),
    }

"""Fixture data for the Day 138 lab. All of it is invented for teaching.

There is no real dataset anywhere in this lab, and that is deliberate. A
lesson about who is missing from a dataset should not be taught by taking
a dataset about real people whose consent nobody in this room obtained,
and a synthetic population has a property no real one has: its true
composition is known exactly, so "the sample under-represents group B by a
factor of ten" is a checkable fact rather than an estimate.

Every name here is invented. Every number is chosen to make a mechanism
visible.
"""

from __future__ import annotations

from typing import Any

# ---------------------------------------------------------------------------
# Exercise 2 -- the reference population the sample claims to describe
# ---------------------------------------------------------------------------

#: What the published statistical abstract says the population looks like.
#: This is the "reference distribution" the coverage check needs: without
#: one, a sample's composition cannot be judged at all, only described.
REFERENCE_SHARES: dict[str, float] = {
    "north": 0.40,
    "south": 0.35,
    "east": 0.15,
    "west": 0.10,
}

#: What the survey actually collected. The west is reached at a tenth of
#: the rate the reference population says it should be -- the survey ran by
#: telephone during working hours, and the west's shift patterns are
#: different. Every other region is close to its reference share.
SURVEY_COUNTS: dict[str, int] = {
    "north": 4_200,
    "south": 3_700,
    "east": 1_600,
    "west": 100,
}

# ---------------------------------------------------------------------------
# Exercise 8 -- datasheets, complete and incomplete
# ---------------------------------------------------------------------------

#: A dataset shipped the way most datasets are shipped: a file, a title,
#: and a vague sense that somebody official produced it. Seven of the eleven
#: required provenance fields are absent, one is present but empty, and one
#: is present but null. Each of those is a different way of not answering.
INCOMPLETE_DATASHEET: dict[str, Any] = {
    "collector": "Regional Statistics Office",
    "collection_period": "2025",
    "purpose": "",
    "population_definition": None,
    "licence": "CC-BY-4.0",
}

#: The same dataset, documented. Nothing here is harder to produce than the
#: file itself; it is simply written down at the time rather than
#: reconstructed from memory two years later by somebody who was not there.
COMPLETE_DATASHEET: dict[str, Any] = {
    "collector": "Regional Statistics Office, Household Surveys Unit",
    "collection_period": "2025-01-06/2025-12-19",
    "purpose": (
        "Quarterly regional labour-force estimates for internal planning; "
        "not designed for small-area or sub-annual estimates"
    ),
    "population_definition": (
        "Residents aged 16 and over in private households in the four "
        "administrative regions"
    ),
    "sampling_frame": (
        "Landline and mobile telephone numbers on the national directory, "
        "called between 09:00 and 17:00 on weekdays"
    ),
    "inclusion_criteria": [
        "aged 16 or over on the first day of the reference week",
        "resident in a private household",
        "answered and consented to the full interview",
    ],
    "exclusion_criteria": [
        "communal establishments (halls of residence, care homes, barracks)",
        "households with no telephone number on the directory",
        "partial interviews terminated before the employment block",
    ],
    "known_gaps": [
        "the west region is reached at roughly a tenth of its population "
        "share because daytime calling misses its dominant shift pattern",
        "no coverage of communal establishments at all",
    ],
    "licence": "CC-BY-4.0",
    "version": "3.1.0",
    "changelog": [
        "3.0.0 first public release",
        "3.1.0 corrected regional weights for the east; no change to raw responses",
    ],
}

#: The eleven required fields, restated here so a learner can see the target
#: without opening ``ethics.py``. Kept in the same order as
#: ``ethics.DATASHEET_REQUIRED_FIELDS``.
EXPECTED_MISSING_FROM_INCOMPLETE: list[str] = [
    "purpose",
    "population_definition",
    "sampling_frame",
    "inclusion_criteria",
    "exclusion_criteria",
    "known_gaps",
    "version",
    "changelog",
]

# ---------------------------------------------------------------------------
# Named proxies, for the exercise that asks you to say what you are actually
# measuring. Each pair is (the proxy in the data, the thing it stands in for).
# ---------------------------------------------------------------------------

PROXY_PAIRS: tuple[tuple[str, str], ...] = (
    ("arrests recorded", "offences committed"),
    ("clicks", "interest"),
    ("diagnosis codes billed", "disease present"),
    ("prior spending on care", "need for care"),
    ("hours logged in the tracker", "work done"),
)

"""Reference answers for Day 138 — "Bias You Can Measure".

Nine exercises, all solved. Read `starter/test_ethics.py` and try them
first; this file is the answer key, and reading it early costs you the
only part of the day that is actually difficult.

Every assertion here is a claim about a synthetic table constructed in
`ethics.py` with known properties. No real dataset about real people is
read anywhere in this lab.
"""

from __future__ import annotations

import pytest

import ethics as et
import fixtures as fx


# ---------------------------------------------------------------------------
# 1. More data does not fix a biased frame
# ---------------------------------------------------------------------------


def test_01_more_data_does_not_fix_a_biased_frame():
    ladder = et.bias_variance_ladder(sizes=(500, 5_000, 50_000), replicates=40, seed=138)

    assert list(ladder["n"]) == [500, 5_000, 50_000]

    # The model looks fine. In-sample RMSE is about the size of the noise the
    # data was generated with, at every sample size -- nothing here says
    # "something is wrong".
    for rmse in ladder["in_sample_rmse"]:
        assert 1.0 < rmse < 1.4

    # The BIAS for the under-represented group does not move. It sits at the
    # intercept gap scaled by how little of group B the frame lets through:
    # 6.0 * (1 - 0.01) = 5.94, and it is still 5.94 after a hundredfold more
    # data.
    biases = [abs(b) for b in ladder["bias_b"]]
    for bias in biases:
        assert 5.85 < bias < 6.00
    assert max(biases) - min(biases) < 0.05, (
        "the bias must be flat across a hundredfold increase in n"
    )

    # The VARIANCE does move, and it moves the way sampling theory says it
    # should: roughly a factor of sqrt(10) per tenfold increase in n.
    sds = list(ladder["sd_b"])
    for smaller_n_sd, larger_n_sd in zip(sds, sds[1:]):
        assert smaller_n_sd / larger_n_sd > 2.2
    assert sds[0] / sds[-1] > 5.0

    # Meanwhile group A -- the group the frame does reach -- is fine. The
    # pooled model is wrong for exactly the group nobody measured.
    for bias_a in ladder["bias_a"]:
        assert abs(bias_a) < 0.15

    # And the sample's own composition never converges on the population's:
    # the frame is 1% group B whatever n is, against a true 10%.
    for share in ladder["sample_share_b"]:
        assert abs(share - et.FRAME_SHARE_B) < 0.005
        assert share < et.POPULATION_SHARE_B / 5


# ---------------------------------------------------------------------------
# 2. Coverage mismatch, computed
# ---------------------------------------------------------------------------


def test_02_coverage_mismatch_is_computable():
    report = et.coverage_report(fx.SURVEY_COUNTS, fx.REFERENCE_SHARES)

    # The single summary number is real but small, which is exactly why it is
    # not enough on its own: one region can be almost absent while the
    # distance between the two distributions still reads as "0.09".
    assert report["total_variation_distance"] == pytest.approx(0.0895833333, abs=1e-6)

    # The per-group ratio names the problem and quantifies it. The west
    # appears at roughly a tenth of its population share.
    ratios = report["representation_ratio"]
    assert ratios["west"] == pytest.approx(0.10416666, abs=1e-6)
    assert ratios["west"] < 0.15

    # Every other region is close to where it should be.
    for region in ("north", "south", "east"):
        assert 0.9 < ratios[region] < 1.2

    assert report["flagged"] == ["west"]
    assert report["worst_group"] == "west"

    # A sample that mirrors the reference produces a distance of exactly zero
    # and flags nobody -- the check is capable of saying "fine".
    matching = {"north": 400, "south": 350, "east": 150, "west": 100}
    clean = et.coverage_report(matching, fx.REFERENCE_SHARES)
    assert clean["total_variation_distance"] == pytest.approx(0.0, abs=1e-12)
    assert clean["flagged"] == []


# ---------------------------------------------------------------------------
# 3. The proxy gap
# ---------------------------------------------------------------------------


def test_03_optimising_a_proxy_hurts_the_group_it_mismeasures():
    result = et.proxy_gap_experiment(seed=138)

    # The proxy looks excellent by the check most people run.
    assert result["proxy_correlation"] > 0.90

    # Both groups are drawn from the same true-need distribution, so a fair
    # selection would take group B at roughly its population share of 0.25.
    # Ranking by the target does exactly that.
    assert result["selected_share_b_by_target"] == pytest.approx(0.25, abs=0.05)

    # Ranking by the proxy does not. Group B is selected at a small fraction
    # of its share, because the proxy records its need at 85% of the true
    # level and that is enough to push it below the cut.
    assert result["selected_share_b_by_proxy"] < 0.5 * result["population_share_b"]
    assert result["selected_share_b_by_proxy"] < 0.10

    # The true outcome is worse for group B, and by a lot.
    assert result["b_need_served_by_proxy"] < 0.35 * result["b_need_served_by_target"]

    # And the aggregate barely moves, which is why nobody notices: total need
    # served falls by under 3% while group B's share of it collapses.
    total_ratio = result["need_served_by_proxy"] / result["need_served_by_target"]
    assert 0.97 < total_ratio < 1.0

    # A named proxy and the thing it stands in for are two different things,
    # and the fixture list says so out loud.
    assert ("arrests recorded", "offences committed") in fx.PROXY_PAIRS
    for proxy, target in fx.PROXY_PAIRS:
        assert proxy != target


# ---------------------------------------------------------------------------
# 4. Aggregation bias
# ---------------------------------------------------------------------------


def test_04_one_pooled_model_is_worse_for_every_subgroup():
    result = et.aggregation_bias_experiment(n_per_group=2_000, seed=138)

    # Both subgroups genuinely trend downward.
    for fit in result["per_group_fit"].values():
        assert fit["slope"] < -0.9
        assert fit["slope"] > -1.1

    # The pooled fit trends upward. It is not merely imprecise; it has the
    # wrong sign for every subgroup in the data.
    assert result["pooled_slope"] > 1.5
    assert result["sign_flip"] is True

    pooled = result["pooled_rmse_by_group"]
    own = result["per_group_rmse_by_group"]
    assert set(pooled) == {"A", "B"}

    # Worse for EVERY subgroup, not worse on average.
    for group in ("A", "B"):
        assert pooled[group] > own[group]
        assert pooled[group] > 4.0
        assert own[group] < 1.2
    assert result["pooled_worse_for_every_group"] is True


# ---------------------------------------------------------------------------
# 5. The fairness tension, demonstrated
# ---------------------------------------------------------------------------


def test_05_fairness_criteria_are_in_genuine_tension(fairness_pop):
    # The score is perfectly calibrated for both groups: in every score bin,
    # the observed positive rate IS the score.
    calibration = et.calibration_by_score_bin(fairness_pop)
    assert float(calibration["deviation"].max()) == pytest.approx(0.0, abs=1e-12)

    # The two groups have different base rates. That single fact is the whole
    # engine of what follows.
    rates = et.base_rates(fairness_pop)
    assert rates["A"] == pytest.approx(0.66, abs=1e-9)
    assert rates["B"] == pytest.approx(0.34, abs=1e-9)
    assert rates["A"] != rates["B"]

    single = et.evaluate_policy("single_threshold", fairness_pop)
    parity = et.evaluate_policy("demographic_parity", fairness_pop)
    opportunity = et.evaluate_policy("equal_opportunity", fairness_pop)

    # One threshold for everyone leaves the calibrated score untouched and
    # violates both of the other two criteria.
    assert single["per_group"]["A"]["selection_rate"] == pytest.approx(0.85, abs=1e-9)
    assert single["per_group"]["B"]["selection_rate"] == pytest.approx(0.35, abs=1e-9)
    assert single["parity_gap"] == pytest.approx(0.50, abs=1e-9)
    assert single["equal_opportunity_gap"] > 0.30
    assert single["precision_gap"] > 0.10

    # Enforcing demographic parity closes the selection-rate gap EXACTLY and
    # opens the other two.
    assert parity["parity_gap"] == pytest.approx(0.0, abs=1e-9)
    assert parity["equal_opportunity_gap"] > 0.10
    assert parity["precision_gap"] > 0.30

    # Enforcing equal opportunity closes the true-positive-rate gap EXACTLY
    # and reopens the selection-rate gap.
    assert opportunity["equal_opportunity_gap"] == pytest.approx(0.0, abs=1e-9)
    assert opportunity["per_group"]["A"]["true_positive_rate"] == pytest.approx(
        0.80, abs=1e-9
    )
    assert opportunity["per_group"]["B"]["true_positive_rate"] == pytest.approx(
        0.80, abs=1e-9
    )
    assert opportunity["parity_gap"] > 0.10
    assert opportunity["precision_gap"] > 0.25

    # The assertion this exercise exists for: NO policy satisfies all three.
    # This is the incompatibility, not a preference between the policies.
    summary = et.fairness_incompatibility()
    assert summary["any_policy_satisfies_all"] is False
    for policy, flags in summary["satisfied"].items():
        assert not all(flags.values()), f"{policy} cannot satisfy all three"


# ---------------------------------------------------------------------------
# 6. Re-identification risk
# ---------------------------------------------------------------------------


def test_06_uniqueness_on_a_few_quasi_identifiers(register, generalised):
    assert len(register) == 5_000
    assert "diagnosis" in register.columns

    exact = et.unique_row_count(register, et.QUASI_IDS)
    coarse = et.unique_row_count(generalised, et.GENERALISED_QUASI_IDS)

    # Over half the rows are the only row with their birth year, postcode and
    # sex -- in a table that carries no name, no address and no identifier.
    assert exact == 2_723
    assert exact / len(register) > 0.50

    # Coarsening ONE field -- exact birth year to a ten-year band -- collapses
    # that by two orders of magnitude.
    assert coarse == 13
    assert coarse < exact / 100

    # The register itself is unchanged; generalisation added a column rather
    # than editing anyone's record in place.
    assert "birth_year" in generalised.columns
    assert "birth_decade" in generalised.columns
    assert len(generalised) == len(register)


# ---------------------------------------------------------------------------
# 7. k-anonymity as a check with limits
# ---------------------------------------------------------------------------


def test_07_k_anonymity_holds_and_still_discloses(generalised):
    # Generalisation alone does not buy k-anonymity: some bands still hold a
    # single person.
    assert et.k_anonymity_level(generalised, et.GENERALISED_QUASI_IDS) == 1

    kept, suppressed = et.suppress_small_classes(generalised, et.GENERALISED_QUASI_IDS, 5)

    # After suppression the table really is 5-anonymous.
    assert et.k_anonymity_level(kept, et.GENERALISED_QUASI_IDS) >= 5
    assert suppressed == 794
    assert len(kept) == 4_206
    # And it cost 794 rows -- roughly one in six people, all of them the
    # people in unusual combinations.
    assert 0.10 < suppressed / 5_000 < 0.25

    # Now the limit. This table is 4-anonymous by construction.
    table = et.homogeneous_table()
    assert et.k_anonymity_level(table, et.GENERALISED_QUASI_IDS) == 4

    # And it still discloses, because one of its classes is homogeneous in the
    # sensitive attribute: every member of it shares the same diagnosis, so
    # membership alone reveals it. No re-identification was required.
    leaks = et.homogeneous_classes(table, et.GENERALISED_QUASI_IDS, "diagnosis", 4)
    assert len(leaks) == 1
    leak = leaks[0]
    assert leak["size"] == 4
    assert leak["distinct_sensitive_values"] == 1
    assert leak["disclosed_value"] == "diabetes"
    assert leak["class"] == (1970, "1001", "F")


# ---------------------------------------------------------------------------
# 8. A datasheet contract
# ---------------------------------------------------------------------------


def test_08_datasheet_contract_names_every_missing_field():
    bad = et.check_datasheet(fx.INCOMPLETE_DATASHEET)

    assert bad["complete"] is False
    # Named, not counted: "documentation incomplete" is not actionable and
    # "no exclusion_criteria recorded" is.
    assert bad["missing"] == fx.EXPECTED_MISSING_FROM_INCOMPLETE
    assert len(bad["missing"]) == 8

    # An empty string and an explicit None both count as unanswered, because
    # a field that is present but blank tells a reader nothing.
    assert "purpose" in bad["missing"]
    assert "population_definition" in bad["missing"]
    # A field that IS answered is not reported.
    assert "collector" not in bad["missing"]
    assert "licence" not in bad["missing"]

    good = et.check_datasheet(fx.COMPLETE_DATASHEET)
    assert good["complete"] is True
    assert good["missing"] == []
    assert set(good["checked"]) == set(et.DATASHEET_REQUIRED_FIELDS)
    assert len(et.DATASHEET_REQUIRED_FIELDS) == 11

    # The complete datasheet documents the coverage gap from exercise 2 in
    # prose, which is the whole point: the gap is a known property of the
    # sampling frame, written down where the next reader will find it.
    assert any("west" in gap for gap in fx.COMPLETE_DATASHEET["known_gaps"])


# ---------------------------------------------------------------------------
# 9. Version drift
# ---------------------------------------------------------------------------


def test_09_version_drift_that_summary_statistics_hide(versions):
    v1, v2 = versions
    diff = et.diff_versions(v1, v2)

    # Every summary statistic of the measured column is IDENTICAL between the
    # two releases -- not close, identical. A reader diffing describe() output
    # sees nothing at all.
    assert diff["summary_identical"] is True
    for key in ("count", "mean", "std", "min", "median", "max"):
        assert diff["summary_old"][key] == diff["summary_new"][key]

    # And yet a quarter of the sample changed group.
    assert diff["group_shares_old"] == {"A": 0.5, "B": 0.5}
    assert diff["group_shares_new"] == {"A": 0.75, "B": 0.25}
    assert diff["composition_shift"] == pytest.approx(0.25, abs=1e-9)
    assert diff["material_shift"] is True
    assert diff["summary_hides_the_change"] is True

    # The provenance record is what makes the change visible: the sampling
    # frame field differs, and the changelog gained an entry naming why.
    assert diff["frame_changed"] is True
    assert diff["changelog_explains"] is True
    assert len(diff["changelog_new_entries"]) == 1
    assert "quota" in diff["changelog_new_entries"][0]
    assert "50% to 25%" in diff["changelog_new_entries"][0]

    # Both releases carry a complete datasheet, so the comparison above was
    # possible at all.
    for version in (v1, v2):
        assert et.check_datasheet(version.datasheet)["complete"] is True

lesson_id: D138
day: 138
kind: guided-build
languages:
  - python
  - bash
setup_commands:
  - cd labs/sections/math-statistics-and-data/day-138-data-ethics-bias-and-provenance
  - python3 -m venv .venv
  - .venv/bin/pip install -r requirements/requirements.txt
  - '.venv/bin/python3 -c "import numpy, pandas, pytest; print(numpy.__version__, pandas.__version__, pytest.__version__)"'
run_commands:
  - .venv/bin/pytest examples
  - .venv/bin/pytest starter
test_commands:
  - bash tests/run_tests.sh
cleanup_commands:
  - >-
    find . -path ./.venv -prune -o -type d -name '__pycache__' -print -exec rm -rf -- {} +
  - rm -rf .pytest_cache
  - 'rm -rf .venv  # optional: removes the lab virtual environment'
  - 'git checkout -- starter/  # optional: reset your work'
requires_network: false
requires_api_key: false
estimated_minutes: 50
last_executed: '2026-08-20'
executed_on: >-
  macOS 26.5.2 (Apple Silicon, arm64), Python 3.14.0, NumPy 2.5.2, pandas 3.0.5,
  pytest 9.1.1, bash 3.2.57 -- bash tests/run_tests.sh -> 51 checks, 0 failure(s),
  exit=0. pytest examples -> 9 passed in 0.18s. pytest starter -> 9 skipped in 0.01s
  (untouched checkout). Everything was run through a real lab-local .venv created by
  the documented setup commands. The lab is fully offline after that one pip install:
  it opens no socket, starts no server and reads no file, and section 7 of the harness
  greps the lab code for urllib, requests, socket and http URLs and asserts it finds
  none. Section 6 solves every exercise in a scratch copy (9 passed), deliberately
  breaks exercise 6's exact uniqueness assertion (exact == 2_723 -> == 9_999), confirms
  a non-zero exit and a printed failure, restores the file, and confirms 9 passed
  again. Separately, a real assertion inside examples/ethics.py was broken directly --
  the 0.5 factor in coverage_report's total-variation-distance was changed to 1.0 --
  and the whole harness was re-run: it reported 51 checks, 7 failure(s) and exited 1
  across the behaviour script, examples/, and both halves of the fail-then-restore
  proof (the coverage check itself printed the wrong value 0.179167); the file was
  restored and the harness returned to 51 checks, 0 failure(s), exit=0. Section 5
  confirms directly that `pytest examples starter` in one invocation aborts collection
  with `import file mismatch` (both directories define a module named test_ethics.py)
  rather than silently letting one shadow the other. MEASURED RESULTS. Exercise 1, the
  centrepiece: with the frame reaching group B at 1% against a true population share of
  10%, the pooled model's bias for group B was 5.9459 at n=500, 5.9384 at n=5000 and
  5.9397 at n=50000 -- a spread of 0.0075 across a hundredfold increase in n -- while
  the standard deviation of that same estimate fell 0.06025 -> 0.01796 -> 0.00524,
  ratios of 3.354 and 3.425 per tenfold, close to the sqrt(10)=3.162 sampling theory
  predicts, and 11.489x overall. In-sample RMSE was 1.1671, 1.1593 and 1.1649 at the
  three sizes, so the model looked equally fine at every one, and group A's bias stayed
  under 0.062 throughout. Exercise 2: total variation distance 0.089583 with the west
  region at a representation ratio of 0.104167, flagged by name. Exercise 3: proxy
  correlated with the target at 0.9253 yet selected group B at 0.0400 against 0.2500
  when ranking by the target, serving 0.1785 of group B's true need while total need
  served only fell to 0.9831. Exercise 4: pooled slope +1.9785 against true subgroup
  slopes of -0.9870 and -0.9991, pooled RMSE 5.0275 and 5.0545 against per-subgroup
  1.0192 and 1.0040. Exercise 5, all three metrics on an integer-exact, perfectly
  calibrated population with base rates 0.66 and 0.34 -- one threshold at 0.5:
  selection rates 0.8500/0.3500, true-positive rates 0.9470/0.6324, precision
  0.7353/0.6143 (gaps 0.5000, 0.3146, 0.1210); equalised selection rate: 0.6000/0.6000,
  0.7424/0.8529, 0.8167/0.4833 (gaps 0.0000, 0.1105, 0.3333); equalised true-positive
  rate: 0.6560/0.5400, 0.8000/0.8000, 0.8049/0.5037 (gaps 0.1160, 0.0000, 0.3012). No
  policy satisfied all three. Exercise 6: 2723 of 5000 rows unique on birth year,
  postcode and sex, falling to 13 after coarsening birth year to a decade, a reduction
  factor of 209.5. Exercise 7: k=5 achieved by suppressing 794 rows of 5000, and a
  4-anonymous table still disclosed 'diabetes' through one homogeneous class of size 4.
  Exercise 8: 8 of 11 required provenance fields named as missing on the undocumented
  record. Exercise 9: two releases with bit-identical summary statistics (mean
  49.9967280173, std 9.6894130750 in both) differing by 0.2500 in group composition,
  visible only through one new changelog entry. THREE HONESTY CALLS. FIRST: no fairness
  toolkit was installed or run. Fairlearn and AIF360 are described in the lesson from
  their published documentation only, and no output from either is reproduced anywhere
  in this lab or lesson. SECOND: every dataset in this lab is synthetic and constructed
  in ethics.py from numpy.random.default_rng with fixed seeds -- no real dataset about
  real people is read anywhere, which is stated in the lab README, in 00_brief.md and
  in security.md as a deliberate design decision rather than a limitation. THIRD: the
  numeric values above are pinned to NumPy 2.5.2 and seed=138; exercises 2, 5, 8 and 9
  contain no randomness and are identical on any machine, while exercises 1, 3, 4, 6
  and 7 depend on the NumPy random stream, and expected-output/FIELDS.md records
  exactly which is which.

# Requirements

`requirements.txt` pins the exact versions this lab was written and run
against on 2026-08-20. Everything else it uses — `dataclasses`, `typing` —
is in the Python standard library.

Install into a lab-local virtual environment so the pins cannot collide
with anything else on your machine:

```bash
cd labs/sections/math-statistics-and-data/day-138-data-ethics-bias-and-provenance
python3 -m venv .venv
.venv/bin/pip install -r requirements/requirements.txt
```

**Only that install step needs the network.** Nothing else in this lab
opens a socket, reads a file from disk, or downloads anything. Every table
it measures is constructed in `ethics.py` at import time.

## Why each pin is here

- **NumPy 2.5.2** — every random draw goes through
  `numpy.random.default_rng(seed)`, whose stream is stable across NumPy 2.x
  by NumPy's own documented policy. That stability is what turns "the bias
  stays flat" into a number you can assert rather than a trend you have to
  eyeball. `numpy.polyfit` supplies the least-squares fits in exercises 1
  and 4.
- **pandas 3.0.5** — the group-composition, equivalence-class and
  version-diff work is all `groupby` and `value_counts`. Doing it in plain
  dictionaries would be possible and would teach the wrong habit: these are
  exactly the operations you would run on a real frame.
- **pytest 9.1.1** — the exercise runner.

## What is deliberately NOT pinned here

**No fairness toolkit is installed.** Fairlearn and AIF360 both implement
the metrics exercise 5 computes by hand, and the lesson describes both from
their published documentation. Neither was installed, neither was run, and
no output from either is reproduced anywhere in this lab. Computing the
three criteria in twenty lines of pandas is also the point of exercise 5:
you can see that the incompatibility is arithmetic rather than a library's
opinion.

numpy==2.5.2
pandas==3.0.5
pytest==9.1.1

# Bias You Can Measure — the nine exercises

Today's lab is an ethics lab, and ethics labs usually fail in one of two
ways: they hand you a list of principles you cannot act on, or they preach.
Neither is on offer here. Every exercise below is either **something you
compute** or **something you write down and check**. Where a question is
genuinely a value judgement rather than a calculation — and exercise 5 is
exactly that — the lab says so explicitly and asks you to assert the
*disagreement*, not a winner.

## Everything here is synthetic, deliberately

No real dataset about real people appears anywhere in this lab. Every
table is constructed in `ethics.py` from a seeded `numpy` generator. Two
reasons, and both matter:

1. A lesson about who is missing from a dataset should not be taught on
   people who never agreed to be in one.
2. A synthetic population has a property no real one has: **its true
   composition is known exactly.** That turns "the west looks
   under-represented" into "the west appears at 0.104 of its population
   share", which is a fact you can assert rather than an impression you
   can argue about.

## The files

| File | What it is |
| --- | --- |
| `ethics.py` | Every measurement function, grouped into bias measurement, fairness accounting, disclosure risk and documentation. Read this first |
| `fixtures.py` | The survey counts, the reference population, the two datasheets, and a list of named proxies |
| `conftest.py` | Four session fixtures: `register`, `generalised`, `fairness_pop`, `versions` |
| `test_ethics.py` | Your nine exercises. Each currently calls `pytest.skip` |

Replace each `pytest.skip(...)` with real assertions and delete the skip
line. Run `pytest starter -v` as often as you like. **Never run
`pytest examples starter` in one command** — both directories hold a module
named `test_ethics.py` and pytest aborts collection on the clash. Run them
as two separate commands.

---

## Exercise 1 — more data does not fix a biased frame

The centrepiece, and the strongest argument in the whole day, because it is
arithmetic rather than opinion.

Group B is **10% of the population** and the sampling frame reaches it with
probability **1%** — a telephone survey run during working hours, say. Both
groups have the same slope and the same noise; group B's outcome simply
sits six units higher at every x.

Fit one pooled straight line and measure the signed prediction error for
each group, over 40 replicate samples, at n = 500, 5,000 and 50,000.

You are proving three things at once:

- the **in-sample RMSE looks fine** at every n, so nothing in the usual
  model report says anything is wrong;
- the **bias for group B does not move**. It sits at
  `6.0 × (1 − 0.01) = 5.94` and is still 5.94 after a hundredfold more
  data;
- the **variance does move**, falling by roughly `sqrt(10)` per tenfold
  increase — the confidence interval tightens around an answer that stays
  wrong.

This is Day 117's bias-versus-error distinction with consequences. Sampling
error shrinks with n. Sampling bias does not.

## Exercise 2 — coverage mismatch, computed

The sampling frame is a claim about who is in the world. When you have a
**reference distribution** — a census, a published statistical abstract —
the mismatch between frame and population stops being a worry and becomes
a subtraction.

`coverage_report` returns two things on purpose. The **total variation
distance** is a single summary and reads as 0.0896 here, which is small
enough to be reassuring and is exactly why it is not the number that
matters. The **representation ratio** per group is: the west appears at
0.104 of its population share, a tenfold shortfall, named by group.

Assert both, and then assert the check can say "fine" on a sample that
mirrors the reference — a check that never passes teaches nothing.

## Exercise 3 — the proxy gap

**The proxy is not the thing.** Arrests are not offences; clicks are not
interest; a billed diagnosis code is not a disease; prior spending on care
is not need for care.

Here, both groups are drawn from the *same* true-need distribution, but the
proxy records group B's need at 85% of its true level because group B
interacts with the recording system less. The proxy still correlates with
the target at above 0.9 — it passes the check most people run.

Select the top 500 by proxy, and separately the top 500 by true need.
Assert that group B's share of the selection collapses from about 0.25 to
under 0.10, that the true need served for group B falls to under 35% of
what a target-ranked selection would have served — and then assert the part
that explains why nobody notices: **total** need served falls by less than
3%. The aggregate barely moves while one group is gutted.

## Exercise 4 — aggregation bias

Day 116's Simpson's paradox, now with stakes. Group A lives at low x with a
low intercept; group B at high x with a much higher one. **Both groups
trend downward.** Pooled, the between-group offset dominates and the fitted
slope comes out positive.

Assert the sign flip, and then assert the thing that matters more: the
pooled model's RMSE is worse **for every subgroup**, not worse on average.
One model for a heterogeneous population can be wrong for everybody in it
while looking like a confident fit.

## Exercise 5 — the fairness tension, demonstrated

Read this one slowly, because it is the exercise where measurement
genuinely runs out.

The population is integer-exact and hand-built: in every score bin, the
observed positive rate **is** the score, in both groups. The score is
perfectly calibrated. The two groups have different base rates — 0.66 and
0.34 — and that single fact is the whole engine of what follows.

Three policies, three criteria:

| Policy | Selection rate | True-positive rate | Precision among selected |
| --- | --- | --- | --- |
| One threshold at 0.5 | differs by 0.50 | differs by 0.31 | differs by 0.12 |
| Equalise selection rate | **equal** | differs by 0.11 | differs by 0.33 |
| Equalise true-positive rate | differs by 0.12 | **equal** | differs by 0.30 |

Assert that each policy closes its own gap *exactly* and opens at least one
other, and that `any_policy_satisfies_all` is `False`.

**Assert the incompatibility, not a preferred answer.** When base rates
differ between groups, demographic parity, equal opportunity and equal
precision cannot generally all hold at once. That is a theorem about
arithmetic. *Which* of them your system should satisfy is a value judgement
about what the system is for and who bears the cost of being wrong — and
reasonable, well-informed people disagree about it. Your job as an analyst
is to make the choice **visible and declared**, not to let a library
default make it silently.

## Exercise 6 — re-identification risk

"Anonymised" is a claim about a threat model, not a property of a file.

The synthetic register carries no name, no address and no identifier —
only birth year, postcode, sex, and a diagnosis. Count how many rows are
the **only** row with their combination of the first three. It is 2,723 out
of 5,000: over half the table can be singled out by anyone who knows three
ordinary facts about a person.

Then coarsen **one** field — exact birth year to a ten-year band — and count
again. It falls to 13. Assert both counts and the two-orders-of-magnitude
drop. Generalisation is cheap and it works; the reason it is not applied by
default is that nobody counted first.

## Exercise 7 — k-anonymity as a check with limits

A table is *k-anonymous* when every combination of quasi-identifiers that
appears at all appears at least k times. Generalisation alone does not get
you there: after coarsening, the smallest class here still holds one
person. Suppressing every class smaller than 5 does — at a cost of **794
rows**, roughly one in six people, and every one of them somebody in an
unusual combination. Say that cost out loud; it is frequently the people
the analysis was supposed to be about.

Then the limit. `homogeneous_table()` is 4-anonymous by construction, and
it still discloses: one of its classes has four members who all share the
same diagnosis. Knowing that someone is a 1970s-born woman in postcode 1001
and appears in the table reveals her diagnosis exactly. **No
re-identification was required.** k-anonymity held; the secret came out
anyway.

Assert the technique *and* its limit. A check you trust past its stated
guarantee is worse than no check, because it buys confidence it did not
earn.

## Exercise 8 — a datasheet contract

A dataset without provenance is not a dataset; it is a rumour with columns.

`check_datasheet` tests a record against eleven required fields — collector,
collection period, purpose, population definition, sampling frame, inclusion
and exclusion criteria, known gaps, licence, version, changelog. Assert that
the incomplete record fails and that the failure **names each missing field**,
because "documentation incomplete" is not actionable and "no
`exclusion_criteria` recorded" is. Note that an empty string and an explicit
`None` both count as unanswered: a field that is present but blank tells a
reader nothing.

Then assert the complete record passes — and that its `known_gaps` field
documents exercise 2's coverage shortfall in prose, where the next reader
will actually find it.

## Exercise 9 — version drift

Two releases of one dataset. The measured column is byte-for-byte
identical, so **every** summary statistic — count, mean, standard
deviation, minimum, median, maximum — matches exactly between them. A
reader diffing `describe()` output sees nothing at all.

What changed is who is in the file: group B's quota went from half the
sample to a quarter. Assert the composition shift, assert that the summary
hides it, and then assert the only thing that makes it visible — the
sampling-frame field differs and the changelog gained one entry saying why.

That is the day's deliverable in one exercise. The documentation is not
paperwork attached to the analysis. On this comparison it is the *only*
place the change exists.

"""Shared fixtures for the Day 138 lab.

All four fixtures are built from a fixed seed, so every test in the suite
sees exactly the same synthetic table and every assertion below is a claim
about a specific, reproducible set of numbers rather than about a
distribution "on average".
"""

from __future__ import annotations

import pytest

import ethics as et

SEED = 138


@pytest.fixture(scope="session")
def register():
    """The synthetic register: quasi-identifiers plus one sensitive column."""
    return et.synthetic_register(n=5_000, seed=SEED)


@pytest.fixture(scope="session")
def generalised(register):
    """The same register with birth year coarsened to a decade band."""
    return et.generalise_quasi_ids(register)


@pytest.fixture(scope="session")
def fairness_pop():
    """The integer-exact, perfectly calibrated two-group population."""
    return et.fairness_population()


@pytest.fixture(scope="session")
def versions():
    """Two releases of one dataset, differing only in group composition."""
    return et.build_versions(n=2_000, seed=SEED)

"""Bias you can measure — the measurement and documentation functions.

Every function here is deterministic. The ones that draw data take an
explicit seed and use ``numpy.random.default_rng``, so the same seed gives
the same numbers on every machine and every run. Nothing in this module
reads a real dataset about real people: every table is constructed here,
in code, with properties chosen so that the thing being demonstrated is a
fact about the construction rather than a claim about anyone.

The module has four groups of functions:

* measurement of bias  -- ``biased_frame_experiment``, ``coverage_report``,
  ``proxy_gap_experiment``, ``aggregation_bias_experiment``
* fairness accounting  -- ``fairness_population``, ``evaluate_policy``,
  ``calibration_by_score_bin``, ``fairness_incompatibility``
* disclosure risk      -- ``unique_row_count``, ``generalise_quasi_ids``,
  ``k_anonymity_level``, ``suppress_small_classes``, ``homogeneous_classes``
* documentation        -- ``check_datasheet``, ``diff_versions``
"""

from __future__ import annotations

from dataclasses import dataclass, field
from typing import Any, Iterable, Mapping, Sequence

import numpy as np
import pandas as pd

# ---------------------------------------------------------------------------
# 1. More data does not fix a biased frame
# ---------------------------------------------------------------------------

#: The world the sample is drawn from. Group "B" is one tenth of the
#: population and its outcome sits INTERCEPT_GAP units above group "A" at
#: every value of x. The slope is identical for both groups, so the only
#: thing a pooled straight-line fit can get wrong about group B is the
#: intercept -- which makes the bias analytically predictable and therefore
#: checkable rather than merely plausible.
POPULATION_SHARE_B = 0.10
FRAME_SHARE_B = 0.01
INTERCEPT_A = 2.0
INTERCEPT_B = 8.0
INTERCEPT_GAP = INTERCEPT_B - INTERCEPT_A
TRUE_SLOPE = 1.0
NOISE_SD = 1.0


def draw_from_frame(n: int, rng: np.random.Generator) -> pd.DataFrame:
    """Draw ``n`` rows the way the biased sampling frame draws them.

    The frame reaches group B with probability ``FRAME_SHARE_B`` rather than
    its true population share ``POPULATION_SHARE_B``. Everything else about
    the two groups is identical: same x distribution, same slope, same noise.
    """
    group = np.where(rng.random(n) < FRAME_SHARE_B, "B", "A")
    x = rng.uniform(0.0, 10.0, size=n)
    intercept = np.where(group == "B", INTERCEPT_B, INTERCEPT_A)
    y = intercept + TRUE_SLOPE * x + rng.normal(0.0, NOISE_SD, size=n)
    return pd.DataFrame({"group": group, "x": x, "y": y})


def fit_pooled_line(frame: pd.DataFrame) -> tuple[float, float]:
    """Ordinary least squares of y on x, ignoring the group column entirely.

    Returns ``(slope, intercept)``. This is the model almost everyone fits
    first: one line for everybody.
    """
    slope, intercept = np.polyfit(frame["x"].to_numpy(), frame["y"].to_numpy(), 1)
    return float(slope), float(intercept)


def signed_group_error(slope: float, intercept: float, group: str) -> float:
    """Signed prediction error of a pooled line for one group's true line.

    Because both groups share ``TRUE_SLOPE``, the error is the same at every
    x, so a single number describes it exactly: predicted intercept minus
    true intercept, plus the slope error times a representative x. Evaluated
    at x = 5, the midpoint of the x range.
    """
    true_intercept = INTERCEPT_B if group == "B" else INTERCEPT_A
    x_mid = 5.0
    predicted = intercept + slope * x_mid
    truth = true_intercept + TRUE_SLOPE * x_mid
    return float(predicted - truth)


def biased_frame_experiment(
    n: int, replicates: int = 40, seed: int = 138
) -> dict[str, float]:
    """Fit the pooled line ``replicates`` times at sample size ``n``.

    Returns the mean signed error for each group across replicates (the
    BIAS -- what more data cannot fix) and the standard deviation of that
    error across replicates (the VARIANCE -- what more data does fix), plus
    the in-sample RMSE, which is the number that makes the model look fine.
    """
    root = np.random.default_rng(seed)
    seeds = root.integers(0, 2**32 - 1, size=replicates)

    errors_a: list[float] = []
    errors_b: list[float] = []
    rmses: list[float] = []
    b_shares: list[float] = []

    for s in seeds:
        rng = np.random.default_rng(int(s))
        sample = draw_from_frame(n, rng)
        slope, intercept = fit_pooled_line(sample)
        errors_a.append(signed_group_error(slope, intercept, "A"))
        errors_b.append(signed_group_error(slope, intercept, "B"))
        predicted = intercept + slope * sample["x"].to_numpy()
        residual = sample["y"].to_numpy() - predicted
        rmses.append(float(np.sqrt(np.mean(residual**2))))
        b_shares.append(float((sample["group"] == "B").mean()))

    return {
        "n": float(n),
        "replicates": float(replicates),
        "bias_a": float(np.mean(errors_a)),
        "bias_b": float(np.mean(errors_b)),
        "sd_b": float(np.std(errors_b, ddof=1)),
        "sd_a": float(np.std(errors_a, ddof=1)),
        "in_sample_rmse": float(np.mean(rmses)),
        "sample_share_b": float(np.mean(b_shares)),
    }


def bias_variance_ladder(
    sizes: Sequence[int] = (500, 5_000, 50_000),
    replicates: int = 40,
    seed: int = 138,
) -> pd.DataFrame:
    """Run ``biased_frame_experiment`` at each sample size and tabulate it."""
    rows = [biased_frame_experiment(n, replicates=replicates, seed=seed) for n in sizes]
    table = pd.DataFrame(rows)
    table["n"] = table["n"].astype("int64")
    return table


# ---------------------------------------------------------------------------
# 2. Coverage mismatch, computed
# ---------------------------------------------------------------------------


def coverage_report(
    sample_counts: Mapping[str, int],
    reference_shares: Mapping[str, float],
    flag_below: float = 0.5,
    flag_above: float = 2.0,
) -> dict[str, Any]:
    """Compare a sample's group composition against a reference population.

    ``total_variation_distance`` is one half the sum of absolute differences
    between the two share vectors: 0.0 when the sample mirrors the reference,
    1.0 when they share no mass at all. It is a single summary, so it is
    reported alongside the thing that actually names the problem: a
    per-group ``representation_ratio`` of sample share divided by reference
    share. A ratio of 0.1 means the group appears at one tenth of the rate
    the reference population says it should.

    A group is flagged when its ratio falls outside
    ``[flag_below, flag_above]``.
    """
    groups = sorted(set(sample_counts) | set(reference_shares))
    total = sum(sample_counts.values())
    if total <= 0:
        raise ValueError("sample_counts must contain at least one row")

    sample_shares = {g: sample_counts.get(g, 0) / total for g in groups}
    ratios: dict[str, float] = {}
    flagged: list[str] = []
    for g in groups:
        reference = reference_shares.get(g, 0.0)
        if reference <= 0.0:
            ratios[g] = float("inf") if sample_shares[g] > 0 else 1.0
        else:
            ratios[g] = sample_shares[g] / reference
        if ratios[g] < flag_below or ratios[g] > flag_above:
            flagged.append(g)

    tvd = 0.5 * sum(
        abs(sample_shares[g] - reference_shares.get(g, 0.0)) for g in groups
    )

    return {
        "groups": groups,
        "sample_shares": sample_shares,
        "reference_shares": dict(reference_shares),
        "representation_ratio": ratios,
        "total_variation_distance": float(tvd),
        "flagged": flagged,
        "worst_group": min(ratios, key=lambda g: ratios[g]),
    }


# ---------------------------------------------------------------------------
# 3. The proxy gap
# ---------------------------------------------------------------------------

#: Group B's need is recorded at ACCESS_FACTOR_B of its true level, because
#: the proxy measures a recorded interaction rather than the need itself and
#: group B interacts with the recording system less often. Group A's need is
#: recorded faithfully.
ACCESS_FACTOR_A = 1.0
ACCESS_FACTOR_B = 0.85


def proxy_gap_experiment(
    n_a: int = 3_000, n_b: int = 1_000, budget: int = 500, seed: int = 138
) -> dict[str, Any]:
    """Rank by a proxy, then check what that costs against the real target.

    Both groups are drawn from the SAME true-need distribution, so a fair
    procedure would select them in proportion to their size. The proxy
    understates group B's need systematically. Selecting the top ``budget``
    rows by proxy is compared against selecting the top ``budget`` by true
    need, and both the composition and the total true need served are
    reported for each.
    """
    rng = np.random.default_rng(seed)
    n = n_a + n_b
    group = np.array(["A"] * n_a + ["B"] * n_b)
    true_need = rng.normal(50.0, 10.0, size=n)
    access = np.where(group == "B", ACCESS_FACTOR_B, ACCESS_FACTOR_A)
    proxy = true_need * access + rng.normal(0.0, 2.0, size=n)

    frame = pd.DataFrame({"group": group, "true_need": true_need, "proxy": proxy})

    by_proxy = frame.nlargest(budget, "proxy")
    by_target = frame.nlargest(budget, "true_need")

    population_share_b = n_b / n
    return {
        "budget": budget,
        "population_share_b": float(population_share_b),
        "proxy_correlation": float(frame["proxy"].corr(frame["true_need"])),
        "selected_share_b_by_proxy": float((by_proxy["group"] == "B").mean()),
        "selected_share_b_by_target": float((by_target["group"] == "B").mean()),
        "need_served_by_proxy": float(by_proxy["true_need"].sum()),
        "need_served_by_target": float(by_target["true_need"].sum()),
        "b_need_served_by_proxy": float(
            by_proxy.loc[by_proxy["group"] == "B", "true_need"].sum()
        ),
        "b_need_served_by_target": float(
            by_target.loc[by_target["group"] == "B", "true_need"].sum()
        ),
    }


# ---------------------------------------------------------------------------
# 4. Aggregation bias
# ---------------------------------------------------------------------------


def aggregation_bias_experiment(
    n_per_group: int = 2_000, seed: int = 138
) -> dict[str, Any]:
    """One pooled line against two per-group lines on the same rows.

    Group A lives at low x with a high intercept; group B lives at high x
    with a much higher intercept. Both groups have a NEGATIVE true slope.
    Pooled, the group offset dominates and the fitted slope comes out
    positive: a single model that is wrong in direction for every subgroup
    in the data, while looking like a confident fit overall.
    """
    rng = np.random.default_rng(seed)
    x_a = rng.uniform(0.0, 5.0, size=n_per_group)
    y_a = 10.0 - 1.0 * x_a + rng.normal(0.0, 1.0, size=n_per_group)
    x_b = rng.uniform(5.0, 10.0, size=n_per_group)
    y_b = 30.0 - 1.0 * x_b + rng.normal(0.0, 1.0, size=n_per_group)

    frame = pd.DataFrame(
        {
            "group": np.array(["A"] * n_per_group + ["B"] * n_per_group),
            "x": np.concatenate([x_a, x_b]),
            "y": np.concatenate([y_a, y_b]),
        }
    )

    pooled_slope, pooled_intercept = fit_pooled_line(frame)

    per_group: dict[str, dict[str, float]] = {}
    pooled_rmse: dict[str, float] = {}
    grouped_rmse: dict[str, float] = {}
    for name, part in frame.groupby("group", observed=True):
        slope, intercept = np.polyfit(part["x"].to_numpy(), part["y"].to_numpy(), 1)
        per_group[str(name)] = {"slope": float(slope), "intercept": float(intercept)}
        own = intercept + slope * part["x"].to_numpy()
        pool = pooled_intercept + pooled_slope * part["x"].to_numpy()
        grouped_rmse[str(name)] = float(
            np.sqrt(np.mean((part["y"].to_numpy() - own) ** 2))
        )
        pooled_rmse[str(name)] = float(
            np.sqrt(np.mean((part["y"].to_numpy() - pool) ** 2))
        )

    return {
        "pooled_slope": pooled_slope,
        "pooled_intercept": pooled_intercept,
        "per_group_fit": per_group,
        "pooled_rmse_by_group": pooled_rmse,
        "per_group_rmse_by_group": grouped_rmse,
        "pooled_worse_for_every_group": all(
            pooled_rmse[g] > grouped_rmse[g] for g in pooled_rmse
        ),
        "sign_flip": pooled_slope > 0
        and all(fit["slope"] < 0 for fit in per_group.values()),
    }


# ---------------------------------------------------------------------------
# 5. The fairness tension
# ---------------------------------------------------------------------------

#: A fully specified, integer-exact population. Each entry is
#: (group, score, count). Within every (group, score) cell exactly
#: ``score`` of the rows are positive, so the score is PERFECTLY calibrated
#: for both groups by construction -- no fitting, no randomness, no
#: approximation. The two groups differ only in where their mass sits, which
#: gives them different base rates: that difference is the whole engine of
#: the incompatibility.
FAIRNESS_CELLS: tuple[tuple[str, float, int], ...] = (
    ("A", 0.1, 50),
    ("A", 0.3, 100),
    ("A", 0.5, 200),
    ("A", 0.7, 300),
    ("A", 0.9, 350),
    ("B", 0.1, 350),
    ("B", 0.3, 300),
    ("B", 0.5, 200),
    ("B", 0.7, 100),
    ("B", 0.9, 50),
)


def fairness_population() -> pd.DataFrame:
    """The constructed population as a tidy frame of score cells.

    Columns: group, score, count, positives. ``positives`` is exactly
    ``score * count`` and every one of those products is a whole number, so
    nothing here is rounded.
    """
    rows = []
    for group, score, count in FAIRNESS_CELLS:
        positives = score * count
        assert abs(positives - round(positives)) < 1e-9, "cell must be integer-exact"
        rows.append(
            {
                "group": group,
                "score": score,
                "count": count,
                "positives": int(round(positives)),
            }
        )
    return pd.DataFrame(rows)


def base_rates(population: pd.DataFrame | None = None) -> dict[str, float]:
    """Fraction of each group that is actually positive."""
    pop = fairness_population() if population is None else population
    grouped = pop.groupby("group", observed=True)[["count", "positives"]].sum()
    return {
        str(g): float(row["positives"] / row["count"])
        for g, row in grouped.iterrows()
    }


def calibration_by_score_bin(population: pd.DataFrame | None = None) -> pd.DataFrame:
    """Observed positive rate in every (group, score) cell.

    Calibration means the score can be read as a probability: among rows
    scored 0.7, seven in ten really are positive, in EVERY group. This
    returns the observed rate so the claim can be checked rather than
    asserted.
    """
    pop = fairness_population() if population is None else population
    out = pop.copy()
    out["observed_rate"] = out["positives"] / out["count"]
    out["deviation"] = (out["observed_rate"] - out["score"]).abs()
    return out[["group", "score", "count", "observed_rate", "deviation"]]


def _select_from_top(part: pd.DataFrame, target: float, key: str) -> float:
    """Select ``target`` units of ``key`` from the highest scores downward.

    Returns the number of ROWS selected. When the target lands inside a score
    band, that band is taken proportionally -- the standard way to hit an
    exact rate when a whole band would overshoot. Returns a float because a
    proportional slice of a band generally is not a whole number of rows.
    """
    ordered = part.sort_values("score", ascending=False)
    remaining = target
    rows_selected = 0.0
    for _, band in ordered.iterrows():
        available = float(band[key])
        if available <= 0:
            continue
        if remaining >= available - 1e-9:
            rows_selected += float(band["count"])
            remaining -= available
        else:
            fraction = remaining / available
            rows_selected += fraction * float(band["count"])
            remaining = 0.0
        if remaining <= 1e-9:
            break
    return rows_selected


def _metrics_for_selection(part: pd.DataFrame, rows_selected: float) -> dict[str, float]:
    """Selection rate, true-positive rate and precision for one group.

    ``rows_selected`` is a count of rows taken from the top score band
    downward; the positives inside that selection are accumulated band by
    band, taking a proportional slice of whichever band the cut falls in.
    """
    ordered = part.sort_values("score", ascending=False)
    remaining = rows_selected
    positives_selected = 0.0
    for _, band in ordered.iterrows():
        count = float(band["count"])
        if remaining >= count - 1e-9:
            positives_selected += float(band["positives"])
            remaining -= count
        else:
            fraction = remaining / count if count else 0.0
            positives_selected += fraction * float(band["positives"])
            remaining = 0.0
        if remaining <= 1e-9:
            break

    total = float(part["count"].sum())
    total_positives = float(part["positives"].sum())
    return {
        "selection_rate": rows_selected / total,
        "true_positive_rate": positives_selected / total_positives,
        "precision": positives_selected / rows_selected if rows_selected else 0.0,
        "rows_selected": rows_selected,
        "positives_selected": positives_selected,
    }


def evaluate_policy(
    policy: str, population: pd.DataFrame | None = None
) -> dict[str, Any]:
    """Score one decision policy on all three fairness criteria at once.

    Policies:

    * ``"single_threshold"`` -- one cut at score >= 0.5 for everybody. This
      is the policy that leaves the calibrated score untouched.
    * ``"demographic_parity"`` -- give both groups the same selection rate,
      set to the rate the single threshold produces overall.
    * ``"equal_opportunity"`` -- give both groups the same true-positive
      rate, set to the rate that policy targets.

    Reported gaps are absolute differences between the two groups:
    ``parity_gap`` (selection rate), ``equal_opportunity_gap`` (true-positive
    rate) and ``precision_gap`` (precision among the selected, the
    calibration-of-the-decision criterion). None of the three is privileged
    here; all three are printed so the trade can be seen rather than
    argued about.
    """
    pop = fairness_population() if population is None else population
    groups = sorted(pop["group"].unique())
    parts = {g: pop[pop["group"] == g] for g in groups}

    if policy == "single_threshold":
        rows = {
            g: float(parts[g].loc[parts[g]["score"] >= 0.5, "count"].sum())
            for g in groups
        }
    elif policy == "demographic_parity":
        selected_at_threshold = sum(
            float(parts[g].loc[parts[g]["score"] >= 0.5, "count"].sum())
            for g in groups
        )
        total = float(pop["count"].sum())
        rate = selected_at_threshold / total
        rows = {g: rate * float(parts[g]["count"].sum()) for g in groups}
    elif policy == "equal_opportunity":
        target_tpr = 0.80
        rows = {
            g: _select_from_top(
                parts[g], target_tpr * float(parts[g]["positives"].sum()), "positives"
            )
            for g in groups
        }
    else:  # pragma: no cover - guarded by the caller
        raise ValueError(f"unknown policy: {policy!r}")

    per_group = {g: _metrics_for_selection(parts[g], rows[g]) for g in groups}
    a, b = groups[0], groups[1]
    return {
        "policy": policy,
        "per_group": per_group,
        "parity_gap": abs(
            per_group[a]["selection_rate"] - per_group[b]["selection_rate"]
        ),
        "equal_opportunity_gap": abs(
            per_group[a]["true_positive_rate"] - per_group[b]["true_positive_rate"]
        ),
        "precision_gap": abs(per_group[a]["precision"] - per_group[b]["precision"]),
    }


FAIRNESS_POLICIES = ("single_threshold", "demographic_parity", "equal_opportunity")


def fairness_incompatibility(tolerance: float = 1e-9) -> dict[str, Any]:
    """Evaluate every policy and report which criteria each one satisfies.

    Returns the three policies' gap triples plus ``any_policy_satisfies_all``,
    which is the claim this exercise exists to test. It is False here, and
    that is a consequence of the two groups having different base rates --
    not of the policies being badly chosen.
    """
    results = {p: evaluate_policy(p) for p in FAIRNESS_POLICIES}
    satisfied = {
        p: {
            "demographic_parity": r["parity_gap"] <= tolerance,
            "equal_opportunity": r["equal_opportunity_gap"] <= tolerance,
            "equal_precision": r["precision_gap"] <= tolerance,
        }
        for p, r in results.items()
    }
    return {
        "results": results,
        "satisfied": satisfied,
        "base_rates": base_rates(),
        "any_policy_satisfies_all": any(
            all(flags.values()) for flags in satisfied.values()
        ),
    }


# ---------------------------------------------------------------------------
# 6-7. Disclosure risk: uniqueness, generalisation, k-anonymity and its limit
# ---------------------------------------------------------------------------

QUASI_IDS = ("birth_year", "postcode", "sex")
GENERALISED_QUASI_IDS = ("birth_decade", "postcode", "sex")


def synthetic_register(n: int = 5_000, seed: int = 138) -> pd.DataFrame:
    """A synthetic register: no real person is represented here.

    Three quasi-identifiers -- birth year, postcode, sex -- plus one
    sensitive attribute, ``diagnosis``. Nothing in this table was collected;
    every value is drawn from a seeded generator, which is the point: the
    uniqueness counts below are properties of the shape of the table, not
    facts about anybody.
    """
    rng = np.random.default_rng(seed)
    return pd.DataFrame(
        {
            "birth_year": rng.integers(1940, 2006, size=n),
            "postcode": rng.integers(1000, 1060, size=n).astype(str),
            "sex": rng.choice(["F", "M"], size=n),
            "diagnosis": rng.choice(
                ["none", "asthma", "diabetes", "hypertension"],
                size=n,
                p=[0.55, 0.15, 0.15, 0.15],
            ),
        }
    )


def unique_row_count(frame: pd.DataFrame, quasi_ids: Iterable[str]) -> int:
    """Rows that are the ONLY row with their combination of quasi-identifiers.

    A unique row can be singled out by anyone who knows those few fields
    about a person, whether or not the table carries a name.
    """
    keys = list(quasi_ids)
    sizes = frame.groupby(keys, observed=True)[keys[0]].transform("size")
    return int((sizes == 1).sum())


def generalise_quasi_ids(frame: pd.DataFrame) -> pd.DataFrame:
    """Replace exact birth year with the decade it falls in.

    One field, coarsened one step. Nothing else changes -- so any drop in
    uniqueness is attributable to this single edit.
    """
    out = frame.copy()
    out["birth_decade"] = (out["birth_year"] // 10 * 10).astype("int64")
    return out


def k_anonymity_level(frame: pd.DataFrame, quasi_ids: Iterable[str]) -> int:
    """The size of the smallest equivalence class.

    A table is k-anonymous when every combination of quasi-identifiers that
    appears at all appears at least k times, so this minimum IS the k the
    table achieves.
    """
    keys = list(quasi_ids)
    return int(frame.groupby(keys, observed=True)[keys[0]].size().min())


def suppress_small_classes(
    frame: pd.DataFrame, quasi_ids: Iterable[str], k: int
) -> tuple[pd.DataFrame, int]:
    """Drop every row whose equivalence class is smaller than ``k``.

    Returns the surviving frame and the number of rows suppressed. This is
    what k-anonymity costs: it is bought with rows, and the rows it costs
    are exactly the unusual ones -- which is worth saying out loud, because
    the people in rare categories are frequently the people an analysis was
    supposed to be about.
    """
    keys = list(quasi_ids)
    sizes = frame.groupby(keys, observed=True)[keys[0]].transform("size")
    kept = frame[sizes >= k]
    return kept.reset_index(drop=True), int(len(frame) - len(kept))


#: A small, hand-built table where k-anonymity holds and still discloses.
#: Every equivalence class has four rows, so the table is 4-anonymous. The
#: class (1970, "1001", "F") has four rows that all carry the same
#: diagnosis, so knowing only that someone is a 1970s-born woman in postcode
#: 1001 and appears in this table reveals her diagnosis exactly. No
#: re-identification was needed; the anonymity held and the secret still
#: came out.
HOMOGENEOUS_TABLE_ROWS: tuple[dict[str, Any], ...] = tuple(
    [
        {"birth_decade": 1970, "postcode": "1001", "sex": "F", "diagnosis": "diabetes"}
    ]
    * 4
    + [
        {"birth_decade": 1970, "postcode": "1001", "sex": "M", "diagnosis": "none"},
        {"birth_decade": 1970, "postcode": "1001", "sex": "M", "diagnosis": "asthma"},
        {"birth_decade": 1970, "postcode": "1001", "sex": "M", "diagnosis": "none"},
        {
            "birth_decade": 1970,
            "postcode": "1001",
            "sex": "M",
            "diagnosis": "hypertension",
        },
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "none"},
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "asthma"},
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "none"},
        {"birth_decade": 1980, "postcode": "1002", "sex": "F", "diagnosis": "diabetes"},
    ]
)


def homogeneous_table() -> pd.DataFrame:
    """The 4-anonymous table that still leaks, as a frame."""
    return pd.DataFrame(list(HOMOGENEOUS_TABLE_ROWS))


def homogeneous_classes(
    frame: pd.DataFrame, quasi_ids: Iterable[str], sensitive: str, k: int
) -> list[dict[str, Any]]:
    """Equivalence classes of size >= k whose sensitive values are all identical.

    These are the classes where k-anonymity is satisfied and the sensitive
    attribute is disclosed anyway. In the l-diversity vocabulary these are
    the classes with l = 1.
    """
    keys = list(quasi_ids)
    leaks: list[dict[str, Any]] = []
    for key, part in frame.groupby(keys, observed=True):
        values = set(part[sensitive])
        if len(part) >= k and len(values) == 1:
            leaks.append(
                {
                    "class": tuple(key) if isinstance(key, tuple) else (key,),
                    "size": int(len(part)),
                    "distinct_sensitive_values": len(values),
                    "disclosed_value": next(iter(values)),
                }
            )
    return sorted(leaks, key=lambda item: str(item["class"]))


# ---------------------------------------------------------------------------
# 8. A datasheet contract
# ---------------------------------------------------------------------------

#: The provenance fields a dataset must carry before anything is built on
#: it. Each one answers a question that cannot be recovered from the data
#: once it is lost: who made these decisions, when, for what, and about whom.
DATASHEET_REQUIRED_FIELDS: tuple[str, ...] = (
    "collector",
    "collection_period",
    "purpose",
    "population_definition",
    "sampling_frame",
    "inclusion_criteria",
    "exclusion_criteria",
    "known_gaps",
    "licence",
    "version",
    "changelog",
)


def check_datasheet(
    record: Mapping[str, Any], required: Sequence[str] = DATASHEET_REQUIRED_FIELDS
) -> dict[str, Any]:
    """Check a dataset's documentation against the required provenance fields.

    A field counts as missing when it is absent, ``None``, or an empty
    string or collection -- an empty ``known_gaps`` list is a real answer
    only if it is deliberate, so this contract requires it to be written
    down rather than left off. Returns the missing fields BY NAME, because
    "documentation incomplete" is not actionable and
    "no exclusion_criteria recorded" is.
    """
    missing = []
    for name in required:
        if name not in record:
            missing.append(name)
            continue
        value = record[name]
        if value is None:
            missing.append(name)
        elif isinstance(value, str) and not value.strip():
            missing.append(name)
        elif isinstance(value, (list, tuple, dict, set)) and len(value) == 0:
            missing.append(name)
    return {
        "complete": not missing,
        "missing": missing,
        "checked": list(required),
    }


# ---------------------------------------------------------------------------
# 9. Version drift
# ---------------------------------------------------------------------------


@dataclass(frozen=True)
class DatasetVersion:
    """One released version of a dataset: the rows plus the datasheet."""

    version: str
    frame: pd.DataFrame
    datasheet: dict[str, Any] = field(default_factory=dict)


def build_versions(n: int = 2_000, seed: int = 138) -> tuple[DatasetVersion, DatasetVersion]:
    """Two releases of the same dataset, differing only in who is in it.

    The measured column is byte-for-byte identical between the two versions.
    What changed is the sampling frame: version 2 cut group B's quota from
    half the sample to a quarter. Because both groups' values are drawn from
    the same distribution, every summary statistic of ``value`` -- count,
    mean, standard deviation, quartiles, extremes -- is EXACTLY unchanged. A
    reader comparing ``describe()`` output between the two releases sees
    nothing at all.
    """
    rng = np.random.default_rng(seed)
    value = rng.normal(50.0, 10.0, size=n)

    labels_v1 = np.where(np.arange(n) % 2 == 0, "A", "B")
    labels_v2 = np.where(np.arange(n) % 4 == 3, "B", "A")

    common = {
        "collector": "Course fixture generator (synthetic)",
        "collection_period": "2026-01-01/2026-06-30",
        "purpose": "Teaching example for measuring composition drift",
        "population_definition": "Synthetic adults in groups A and B",
        "inclusion_criteria": ["synthetic record generated by this module"],
        "exclusion_criteria": ["none"],
        "licence": "CC0-1.0",
        "known_gaps": ["contains no real people by design"],
    }

    v1 = DatasetVersion(
        version="1.0.0",
        frame=pd.DataFrame({"group": labels_v1, "value": value}),
        datasheet={
            **common,
            "sampling_frame": "equal quota: 50% group A, 50% group B",
            "version": "1.0.0",
            "changelog": ["1.0.0 initial release"],
        },
    )
    v2 = DatasetVersion(
        version="2.0.0",
        frame=pd.DataFrame({"group": labels_v2, "value": value}),
        datasheet={
            **common,
            "sampling_frame": "revised quota: 75% group A, 25% group B",
            "version": "2.0.0",
            "changelog": [
                "1.0.0 initial release",
                "2.0.0 sampling frame re-scoped: group B quota reduced from 50% to 25%",
            ],
        },
    )
    return v1, v2


def summary_stats(frame: pd.DataFrame, column: str = "value") -> dict[str, float]:
    """The summary a reader would actually look at when comparing releases."""
    series = frame[column]
    return {
        "count": float(series.count()),
        "mean": float(series.mean()),
        "std": float(series.std()),
        "min": float(series.min()),
        "median": float(series.median()),
        "max": float(series.max()),
    }


def diff_versions(
    old: DatasetVersion,
    new: DatasetVersion,
    group_column: str = "group",
    value_column: str = "value",
    shift_threshold: float = 0.05,
) -> dict[str, Any]:
    """Compare two releases on summary statistics AND on composition.

    ``summary_identical`` is True when every summary statistic matches to
    within floating-point tolerance -- the case where a release-note-free
    upgrade looks like no change at all. ``composition_shift`` is the largest
    absolute change in any group's share. ``changelog_explains`` reports
    whether the new version's changelog gained an entry the old one did not
    have, which is the only part of this comparison that says WHY.
    """
    old_stats = summary_stats(old.frame, value_column)
    new_stats = summary_stats(new.frame, value_column)
    summary_identical = all(
        abs(old_stats[k] - new_stats[k]) < 1e-9 for k in old_stats
    )

    old_shares = old.frame[group_column].value_counts(normalize=True).to_dict()
    new_shares = new.frame[group_column].value_counts(normalize=True).to_dict()
    groups = sorted(set(old_shares) | set(new_shares))
    shifts = {
        g: float(new_shares.get(g, 0.0) - old_shares.get(g, 0.0)) for g in groups
    }
    composition_shift = max(abs(v) for v in shifts.values()) if shifts else 0.0

    old_log = list(old.datasheet.get("changelog", []))
    new_log = list(new.datasheet.get("changelog", []))
    new_entries = [entry for entry in new_log if entry not in old_log]

    return {
        "old_version": old.version,
        "new_version": new.version,
        "summary_old": old_stats,
        "summary_new": new_stats,
        "summary_identical": summary_identical,
        "group_shares_old": {g: float(old_shares.get(g, 0.0)) for g in groups},
        "group_shares_new": {g: float(new_shares.get(g, 0.0)) for g in groups},
        "group_share_shift": shifts,
        "composition_shift": float(composition_shift),
        "material_shift": composition_shift >= shift_threshold,
        "summary_hides_the_change": summary_identical
        and composition_shift >= shift_threshold,
        "changelog_new_entries": new_entries,
        "changelog_explains": bool(new_entries),
        "frame_changed": old.datasheet.get("sampling_frame")
        != new.datasheet.get("sampling_frame"),
    }

"""Fixture data for the Day 138 lab. All of it is invented for teaching.

There is no real dataset anywhere in this lab, and that is deliberate. A
lesson about who is missing from a dataset should not be taught by taking
a dataset about real people whose consent nobody in this room obtained,
and a synthetic population has a property no real one has: its true
composition is known exactly, so "the sample under-represents group B by a
factor of ten" is a checkable fact rather than an estimate.

Every name here is invented. Every number is chosen to make a mechanism
visible.
"""

from __future__ import annotations

from typing import Any

# ---------------------------------------------------------------------------
# Exercise 2 -- the reference population the sample claims to describe
# ---------------------------------------------------------------------------

#: What the published statistical abstract says the population looks like.
#: This is the "reference distribution" the coverage check needs: without
#: one, a sample's composition cannot be judged at all, only described.
REFERENCE_SHARES: dict[str, float] = {
    "north": 0.40,
    "south": 0.35,
    "east": 0.15,
    "west": 0.10,
}

#: What the survey actually collected. The west is reached at a tenth of
#: the rate the reference population says it should be -- the survey ran by
#: telephone during working hours, and the west's shift patterns are
#: different. Every other region is close to its reference share.
SURVEY_COUNTS: dict[str, int] = {
    "north": 4_200,
    "south": 3_700,
    "east": 1_600,
    "west": 100,
}

# ---------------------------------------------------------------------------
# Exercise 8 -- datasheets, complete and incomplete
# ---------------------------------------------------------------------------

#: A dataset shipped the way most datasets are shipped: a file, a title,
#: and a vague sense that somebody official produced it. Seven of the eleven
#: required provenance fields are absent, one is present but empty, and one
#: is present but null. Each of those is a different way of not answering.
INCOMPLETE_DATASHEET: dict[str, Any] = {
    "collector": "Regional Statistics Office",
    "collection_period": "2025",
    "purpose": "",
    "population_definition": None,
    "licence": "CC-BY-4.0",
}

#: The same dataset, documented. Nothing here is harder to produce than the
#: file itself; it is simply written down at the time rather than
#: reconstructed from memory two years later by somebody who was not there.
COMPLETE_DATASHEET: dict[str, Any] = {
    "collector": "Regional Statistics Office, Household Surveys Unit",
    "collection_period": "2025-01-06/2025-12-19",
    "purpose": (
        "Quarterly regional labour-force estimates for internal planning; "
        "not designed for small-area or sub-annual estimates"
    ),
    "population_definition": (
        "Residents aged 16 and over in private households in the four "
        "administrative regions"
    ),
    "sampling_frame": (
        "Landline and mobile telephone numbers on the national directory, "
        "called between 09:00 and 17:00 on weekdays"
    ),
    "inclusion_criteria": [
        "aged 16 or over on the first day of the reference week",
        "resident in a private household",
        "answered and consented to the full interview",
    ],
    "exclusion_criteria": [
        "communal establishments (halls of residence, care homes, barracks)",
        "households with no telephone number on the directory",
        "partial interviews terminated before the employment block",
    ],
    "known_gaps": [
        "the west region is reached at roughly a tenth of its population "
        "share because daytime calling misses its dominant shift pattern",
        "no coverage of communal establishments at all",
    ],
    "licence": "CC-BY-4.0",
    "version": "3.1.0",
    "changelog": [
        "3.0.0 first public release",
        "3.1.0 corrected regional weights for the east; no change to raw responses",
    ],
}

#: The eleven required fields, restated here so a learner can see the target
#: without opening ``ethics.py``. Kept in the same order as
#: ``ethics.DATASHEET_REQUIRED_FIELDS``.
EXPECTED_MISSING_FROM_INCOMPLETE: list[str] = [
    "purpose",
    "population_definition",
    "sampling_frame",
    "inclusion_criteria",
    "exclusion_criteria",
    "known_gaps",
    "version",
    "changelog",
]

# ---------------------------------------------------------------------------
# Named proxies, for the exercise that asks you to say what you are actually
# measuring. Each pair is (the proxy in the data, the thing it stands in for).
# ---------------------------------------------------------------------------

PROXY_PAIRS: tuple[tuple[str, str], ...] = (
    ("arrests recorded", "offences committed"),
    ("clicks", "interest"),
    ("diagnosis codes billed", "disease present"),
    ("prior spending on care", "need for care"),
    ("hours logged in the tracker", "work done"),
)

"""Your exercises for Day 138 — "Bias You Can Measure".

Nine exercises. Every test below currently calls `pytest.skip(...)` —
replace the skip with real assertions and delete the skip line. Read
`00_brief.md` for the exercise-by-exercise explanation, and `ethics.py`
for the measurement functions you are testing. `fixtures.py` holds the
survey counts, the reference population and the two datasheets.

Check yourself at any point:

    pytest starter -v

Never run `pytest examples starter` in one command — both directories
define a module named `test_ethics.py` and pytest aborts collection on
the clash. Run them as two separate commands.

The reference answer key lives in `examples/test_ethics.py` — read it
AFTER you have tried, never before.

Nothing in this lab reads a real dataset about real people. Every table is
constructed in `ethics.py` from a seeded generator, which is deliberate:
a lesson about who is missing from a dataset should not be taught on
people who never agreed to be in it, and a synthetic population has one
property no real one has — its true composition is known exactly, so
"under-represented by a factor of ten" is checkable rather than estimated.
"""

from __future__ import annotations

import pytest

import ethics as et
import fixtures as fx


def test_01_more_data_does_not_fix_a_biased_frame():
    pytest.skip(
        "The centrepiece. Call et.bias_variance_ladder(sizes=(500, 5000, "
        "50000), replicates=40, seed=138). Assert that in_sample_rmse stays "
        "between 1.0 and 1.4 at every n -- the model always looks fine. Then "
        "assert that abs(bias_b) is between 5.85 and 6.00 at every n AND that "
        "the spread (max minus min) across the three sample sizes is under "
        "0.05: a hundredfold more data does not move it. Then assert the "
        "opposite for variance -- each sd_b is more than 2.2 times the next "
        "one down the ladder, and the first is more than 5 times the last. "
        "Finally assert abs(bias_a) < 0.15 at every n, and that sample_share_b "
        "stays near et.FRAME_SHARE_B rather than converging on "
        "et.POPULATION_SHARE_B"
    )


def test_02_coverage_mismatch_is_computable():
    pytest.skip(
        "Call et.coverage_report(fx.SURVEY_COUNTS, fx.REFERENCE_SHARES). "
        "Assert total_variation_distance is approximately 0.0895833333 -- real, "
        "but small enough on its own to be reassuring, which is why it is not "
        "the number that matters. Assert representation_ratio['west'] is "
        "approximately 0.10416666 and below 0.15, that every other region's "
        "ratio sits between 0.9 and 1.2, that flagged == ['west'] and that "
        "worst_group == 'west'. Then prove the check can say 'fine': call it "
        "again with counts {'north': 400, 'south': 350, 'east': 150, 'west': "
        "100} and assert the distance is 0.0 and flagged is empty"
    )


def test_03_optimising_a_proxy_hurts_the_group_it_mismeasures():
    pytest.skip(
        "Call et.proxy_gap_experiment(seed=138). Assert proxy_correlation is "
        "above 0.90 -- the proxy passes the check most people run. Assert "
        "selected_share_b_by_target is about 0.25 (group B's real share, "
        "since both groups are drawn from the same need distribution). Assert "
        "selected_share_b_by_proxy is below half of population_share_b and "
        "below 0.10. Assert b_need_served_by_proxy is under 35% of "
        "b_need_served_by_target. Then assert the part that explains why "
        "nobody notices: need_served_by_proxy divided by need_served_by_target "
        "is between 0.97 and 1.0. Finally, check fx.PROXY_PAIRS contains "
        "('arrests recorded', 'offences committed') and that no pair has the "
        "proxy equal to the target"
    )


def test_04_one_pooled_model_is_worse_for_every_subgroup():
    pytest.skip(
        "Call et.aggregation_bias_experiment(n_per_group=2000, seed=138). "
        "Assert every entry of per_group_fit has a slope between -1.1 and -0.9 "
        "-- both subgroups genuinely trend downward. Then assert pooled_slope "
        "is above 1.5 and sign_flip is True: the single model has the wrong "
        "SIGN for every subgroup. Then compare pooled_rmse_by_group against "
        "per_group_rmse_by_group and assert that for BOTH 'A' and 'B' the "
        "pooled error is larger, that each pooled error is above 4.0 and each "
        "per-group error below 1.2, and that pooled_worse_for_every_group is "
        "True"
    )


def test_05_fairness_criteria_are_in_genuine_tension(fairness_pop):
    pytest.skip(
        "Call et.calibration_by_score_bin(fairness_pop) and assert the maximum "
        "deviation is 0.0 -- the score is perfectly calibrated for both groups "
        "by construction. Call et.base_rates(fairness_pop) and assert A is 0.66 "
        "and B is 0.34. Then evaluate all three policies with "
        "et.evaluate_policy(name, fairness_pop). For 'single_threshold' assert "
        "the selection rates are 0.85 and 0.35, parity_gap is 0.50, "
        "equal_opportunity_gap exceeds 0.30 and precision_gap exceeds 0.10. For "
        "'demographic_parity' assert parity_gap is exactly 0.0 while "
        "equal_opportunity_gap exceeds 0.10 and precision_gap exceeds 0.30. For "
        "'equal_opportunity' assert equal_opportunity_gap is exactly 0.0 with "
        "both true-positive rates at 0.80, while parity_gap exceeds 0.10 and "
        "precision_gap exceeds 0.25. Finally call et.fairness_incompatibility() "
        "and assert any_policy_satisfies_all is False, and that no policy in "
        "its 'satisfied' map has all three flags True. Assert the "
        "INCOMPATIBILITY -- do not assert that one of the three policies is the "
        "right one, because that is a value judgement and not something this "
        "arithmetic can settle"
    )


def test_06_uniqueness_on_a_few_quasi_identifiers(register, generalised):
    pytest.skip(
        "Assert the register has 5000 rows and a 'diagnosis' column. Call "
        "et.unique_row_count(register, et.QUASI_IDS) and assert it equals 2723 "
        "and is more than half the table -- these are rows that are the ONLY "
        "row with their birth year, postcode and sex, in a table carrying no "
        "name and no identifier. Then call et.unique_row_count(generalised, "
        "et.GENERALISED_QUASI_IDS) and assert it equals 13 and is under one "
        "hundredth of the exact count: coarsening ONE field collapses "
        "uniqueness by two orders of magnitude. Assert the generalised frame "
        "still has both birth_year and birth_decade and the same row count"
    )


def test_07_k_anonymity_holds_and_still_discloses(generalised):
    pytest.skip(
        "Assert et.k_anonymity_level(generalised, et.GENERALISED_QUASI_IDS) is "
        "1 -- generalisation alone does not buy k-anonymity. Call "
        "et.suppress_small_classes(generalised, et.GENERALISED_QUASI_IDS, 5), "
        "assert the surviving table's k level is at least 5, that exactly 794 "
        "rows were suppressed, that 4206 remain, and that the suppressed "
        "fraction is between 0.10 and 0.25 -- the cost is paid in rows, and the "
        "rows it costs are the unusual ones. Then the limit: call "
        "et.homogeneous_table(), assert its k level is 4, and call "
        "et.homogeneous_classes(table, et.GENERALISED_QUASI_IDS, 'diagnosis', "
        "4). Assert there is exactly one leaking class, of size 4, with "
        "distinct_sensitive_values == 1, disclosed_value == 'diabetes' and "
        "class == (1970, '1001', 'F'). k-anonymity held and the sensitive "
        "attribute came out anyway"
    )


def test_08_datasheet_contract_names_every_missing_field():
    pytest.skip(
        "Call et.check_datasheet(fx.INCOMPLETE_DATASHEET). Assert complete is "
        "False and that missing equals fx.EXPECTED_MISSING_FROM_INCOMPLETE "
        "exactly, with 8 entries -- named, not counted. Assert 'purpose' (an "
        "empty string) and 'population_definition' (an explicit None) are both "
        "in missing, while 'collector' and 'licence' are not. Then call it on "
        "fx.COMPLETE_DATASHEET and assert complete is True, missing is empty, "
        "the checked list matches et.DATASHEET_REQUIRED_FIELDS and there are 11 "
        "required fields. Finally assert that at least one entry of "
        "fx.COMPLETE_DATASHEET['known_gaps'] mentions 'west' -- exercise 2's "
        "coverage gap, written down where the next reader will find it"
    )


def test_09_version_drift_that_summary_statistics_hide(versions):
    pytest.skip(
        "Unpack v1, v2 = versions and call et.diff_versions(v1, v2). Assert "
        "summary_identical is True and that count, mean, std, min, median and "
        "max are each EQUAL between the two summaries -- not close, equal. Then "
        "assert group_shares_old is {'A': 0.5, 'B': 0.5}, group_shares_new is "
        "{'A': 0.75, 'B': 0.25}, composition_shift is 0.25, material_shift is "
        "True and summary_hides_the_change is True. Then assert the provenance "
        "record is what makes it visible: frame_changed is True, "
        "changelog_explains is True, there is exactly one new changelog entry, "
        "and that entry contains both 'quota' and '50% to 25%'. Finally assert "
        "et.check_datasheet(version.datasheet)['complete'] is True for both "
        "versions -- the comparison was only possible because both releases "
        "were documented"
    )

#!/usr/bin/env bash
# Tests for the Day 138 lab. Run from the lab directory:
#   bash tests/run_tests.sh
#
# The harness proves the day's claims by driving the real measurement
# functions over the real synthetic populations -- never by reading source
# and never by asserting on a timing:
#
#   * a pooled model fitted on a biased sampling frame keeps a flat bias for
#     the under-represented group across a hundredfold increase in n, while
#     the variance of the estimate falls roughly as sqrt(n);
#   * coverage mismatch against a reference distribution is computed, and
#     the under-represented group is named with a representation ratio;
#   * optimising a proxy that under-records one group's need collapses that
#     group's share of the selection while the aggregate barely moves;
#   * one pooled linear model is worse for EVERY subgroup than per-subgroup
#     models on the same rows, and has the wrong slope sign for both;
#   * demographic parity, equal opportunity and equal precision cannot all
#     hold at once when base rates differ -- each policy closes its own gap
#     exactly and opens another, and no policy closes all three;
#   * over half the rows of a name-free register are unique on three quasi-
#     identifiers, and coarsening one field cuts that by two orders of
#     magnitude;
#   * k-anonymity is achieved, at a counted cost in suppressed rows, and a
#     4-anonymous table is shown still disclosing a sensitive attribute
#     through a homogeneous class;
#   * a datasheet check names every missing provenance field by name, and
#     passes a complete record;
#   * two releases with byte-identical summary statistics differ by a
#     quarter of the sample in composition, and only the provenance record
#     makes that visible;
#   * the reference suite (examples/) passes in full;
#   * the exercise suite (starter/) is all-skip on an untouched checkout,
#     and the harness proves it can genuinely FAIL by solving every
#     exercise in a scratch copy, breaking one assertion on purpose,
#     confirming a non-zero exit and a printed failure, then restoring it;
#   * no __pycache__ or .pytest_cache survives the run.
#
# Everything runs offline. No network, no server, no sudo. Deterministic,
# non-interactive, exits 0 only if every check passes.
set -u

export PYTHONDONTWRITEBYTECODE=1

lab_dir="$(cd "$(dirname "${BASH_SOURCE[0]}")/.." && pwd)"

find "${lab_dir}" -name '.venv' -prune -o -type d -name '__pycache__' -exec rm -rf {} + 2>/dev/null || true
find "${lab_dir}" -name '.venv' -prune -o -type d -name '.pytest_cache' -exec rm -rf {} + 2>/dev/null || true

failures=0
checks=0

check() {
  local label="$1" ok="$2"
  checks=$((checks + 1))
  if [ "${ok}" = "yes" ]; then
    echo "  ok: ${label}"
  else
    echo "  FAIL: ${label}"
    failures=$((failures + 1))
  fi
}

resolve_tool() {
  local tool="$1" override="$2"
  if [ -n "${override}" ] && [ -x "${override}" ]; then echo "${override}"; return 0; fi
  if [ -x "${lab_dir}/.venv/bin/${tool}" ]; then echo "${lab_dir}/.venv/bin/${tool}"; return 0; fi
  if command -v "${tool}" >/dev/null 2>&1; then command -v "${tool}"; return 0; fi
  return 1
}

pytest_bin="$(resolve_tool pytest "${PYTEST:-}")" || {
  echo "FAIL: pytest not found." >&2
  echo "  Install the lab's dependencies with:" >&2
  echo "    python3 -m venv .venv" >&2
  echo "    .venv/bin/pip install -r requirements/requirements.txt" >&2
  echo "  Or point this suite at an existing pytest:" >&2
  echo "    PYTEST=/path/to/pytest bash tests/run_tests.sh" >&2
  exit 1
}

python_bin="$(dirname "${pytest_bin}")/python3"
if [ ! -x "${python_bin}" ]; then
  python_bin="$(command -v python3 || true)"
fi
if [ -z "${python_bin}" ]; then
  echo "FAIL: python3 not found on PATH." >&2
  exit 1
fi

for module in numpy pandas; do
  if ! "${python_bin}" -c "import ${module}" >/dev/null 2>&1; then
    echo "FAIL: ${module} is not importable from ${python_bin}." >&2
    echo "  Install the lab's dependencies with:" >&2
    echo "    python3 -m venv .venv" >&2
    echo "    .venv/bin/pip install -r requirements/requirements.txt" >&2
    exit 1
  fi
done

echo "Day 138 — Bias You Can Measure"
echo

# --------------------------------------------------------------------------
echo "1. The tools and the versions this lab was written against"
# --------------------------------------------------------------------------

versions="$("${python_bin}" - <<'PY'
import platform
from importlib.metadata import version

print(f"python     {platform.python_version()}")
for name in ("numpy", "pandas", "pytest"):
    try:
        print(f"{name:<10} {version(name)}")
    except Exception as exc:  # pragma: no cover
        print(f"{name:<10} NOT INSTALLED ({exc})")
PY
)"
echo "${versions}"
echo

mismatch=0
while IFS= read -r line; do
  [ -z "${line}" ] && continue
  pkg="${line%%==*}"
  pinned="${line#*==}"
  installed="$("${python_bin}" -c "from importlib.metadata import version; print(version('${pkg}'))" 2>/dev/null || echo "MISSING")"
  if [ "${installed}" != "${pinned}" ]; then
    mismatch=1
    echo "  version mismatch: ${pkg} pinned ${pinned}, installed ${installed}"
  fi
done < "${lab_dir}/requirements/requirements.txt"
check "installed packages match requirements.txt exactly" "$( [ ${mismatch} -eq 0 ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "2. The measurements themselves, run over the synthetic populations"
# --------------------------------------------------------------------------

behaviour="$(cd "${lab_dir}/examples" && "${python_bin}" - <<'PY'
"""Drive ethics.py directly and print one machine-readable line per claim,
including the bias-versus-n ladder and the three fairness metrics."""

import ethics as et
import fixtures as fx

results = {}


def record(key, value):
    results[key] = value


# -- exercise 1: bias versus variance as n grows ----------------------------
ladder = et.bias_variance_ladder(sizes=(500, 5_000, 50_000), replicates=40, seed=138)
for _, row in ladder.iterrows():
    n = int(row["n"])
    record(f"bias_b_n{n}", f"{abs(row['bias_b']):.4f}")
    record(f"sd_b_n{n}", f"{row['sd_b']:.5f}")
    record(f"bias_a_n{n}", f"{abs(row['bias_a']):.4f}")
    record(f"rmse_n{n}", f"{row['in_sample_rmse']:.4f}")
    record(f"share_b_n{n}", f"{row['sample_share_b']:.4f}")
biases = [abs(v) for v in ladder["bias_b"]]
sds = list(ladder["sd_b"])
record("bias_b_spread", f"{max(biases) - min(biases):.4f}")
record("bias_b_flat", "yes" if max(biases) - min(biases) < 0.05 else "no")
record("sd_ratio_500_5000", f"{sds[0] / sds[1]:.3f}")
record("sd_ratio_5000_50000", f"{sds[1] / sds[2]:.3f}")
record("sd_ratio_total", f"{sds[0] / sds[2]:.3f}")
record("sd_falls_each_step", "yes" if all(a / b > 2.2 for a, b in zip(sds, sds[1:])) else "no")

# -- exercise 2: coverage mismatch ------------------------------------------
coverage = et.coverage_report(fx.SURVEY_COUNTS, fx.REFERENCE_SHARES)
record("coverage_tvd", f"{coverage['total_variation_distance']:.6f}")
record("coverage_ratio_west", f"{coverage['representation_ratio']['west']:.6f}")
record("coverage_flagged", ",".join(coverage["flagged"]))
record("coverage_worst_group", coverage["worst_group"])
clean = et.coverage_report(
    {"north": 400, "south": 350, "east": 150, "west": 100}, fx.REFERENCE_SHARES
)
record("coverage_clean_tvd", f"{clean['total_variation_distance']:.6f}")
record("coverage_clean_flagged_count", str(len(clean["flagged"])))

# -- exercise 3: the proxy gap ----------------------------------------------
proxy = et.proxy_gap_experiment(seed=138)
record("proxy_correlation", f"{proxy['proxy_correlation']:.4f}")
record("proxy_share_b_by_target", f"{proxy['selected_share_b_by_target']:.4f}")
record("proxy_share_b_by_proxy", f"{proxy['selected_share_b_by_proxy']:.4f}")
record(
    "proxy_b_need_ratio",
    f"{proxy['b_need_served_by_proxy'] / proxy['b_need_served_by_target']:.4f}",
)
record(
    "proxy_total_need_ratio",
    f"{proxy['need_served_by_proxy'] / proxy['need_served_by_target']:.4f}",
)

# -- exercise 4: aggregation bias -------------------------------------------
agg = et.aggregation_bias_experiment(n_per_group=2_000, seed=138)
record("agg_pooled_slope", f"{agg['pooled_slope']:.4f}")
record("agg_slope_a", f"{agg['per_group_fit']['A']['slope']:.4f}")
record("agg_slope_b", f"{agg['per_group_fit']['B']['slope']:.4f}")
record("agg_sign_flip", "yes" if agg["sign_flip"] else "no")
record("agg_pooled_rmse_a", f"{agg['pooled_rmse_by_group']['A']:.4f}")
record("agg_pooled_rmse_b", f"{agg['pooled_rmse_by_group']['B']:.4f}")
record("agg_own_rmse_a", f"{agg['per_group_rmse_by_group']['A']:.4f}")
record("agg_own_rmse_b", f"{agg['per_group_rmse_by_group']['B']:.4f}")
record("agg_pooled_worse_for_every_group", "yes" if agg["pooled_worse_for_every_group"] else "no")

# -- exercise 5: the fairness tension ---------------------------------------
calibration = et.calibration_by_score_bin()
record("fair_max_calibration_deviation", f"{float(calibration['deviation'].max()):.12f}")
rates = et.base_rates()
record("fair_base_rate_a", f"{rates['A']:.4f}")
record("fair_base_rate_b", f"{rates['B']:.4f}")
incompatibility = et.fairness_incompatibility()
for policy, result in incompatibility["results"].items():
    record(f"fair_{policy}_parity_gap", f"{result['parity_gap']:.4f}")
    record(f"fair_{policy}_eo_gap", f"{result['equal_opportunity_gap']:.4f}")
    record(f"fair_{policy}_precision_gap", f"{result['precision_gap']:.4f}")
    for group, metrics in result["per_group"].items():
        record(
            f"fair_{policy}_{group}",
            "sel={:.4f};tpr={:.4f};prec={:.4f}".format(
                metrics["selection_rate"],
                metrics["true_positive_rate"],
                metrics["precision"],
            ),
        )
record("fair_any_policy_satisfies_all", "yes" if incompatibility["any_policy_satisfies_all"] else "no")

# -- exercise 6: re-identification ------------------------------------------
register = et.synthetic_register(n=5_000, seed=138)
generalised = et.generalise_quasi_ids(register)
exact_unique = et.unique_row_count(register, et.QUASI_IDS)
coarse_unique = et.unique_row_count(generalised, et.GENERALISED_QUASI_IDS)
record("reid_rows", str(len(register)))
record("reid_unique_exact", str(exact_unique))
record("reid_unique_generalised", str(coarse_unique))
record("reid_unique_exact_fraction", f"{exact_unique / len(register):.4f}")
record("reid_reduction_factor", f"{exact_unique / coarse_unique:.1f}")

# -- exercise 7: k-anonymity and its limit ----------------------------------
record("k_level_generalised", str(et.k_anonymity_level(generalised, et.GENERALISED_QUASI_IDS)))
kept, suppressed = et.suppress_small_classes(generalised, et.GENERALISED_QUASI_IDS, 5)
record("k_level_after_suppression", str(et.k_anonymity_level(kept, et.GENERALISED_QUASI_IDS)))
record("k_rows_suppressed", str(suppressed))
record("k_rows_kept", str(len(kept)))
homogeneous = et.homogeneous_table()
record("k_homogeneous_table_level", str(et.k_anonymity_level(homogeneous, et.GENERALISED_QUASI_IDS)))
leaks = et.homogeneous_classes(homogeneous, et.GENERALISED_QUASI_IDS, "diagnosis", 4)
record("k_leaking_class_count", str(len(leaks)))
record("k_leaking_class_size", str(leaks[0]["size"]))
record("k_leaking_distinct_values", str(leaks[0]["distinct_sensitive_values"]))
record("k_disclosed_value", str(leaks[0]["disclosed_value"]))

# -- exercise 8: the datasheet contract -------------------------------------
bad = et.check_datasheet(fx.INCOMPLETE_DATASHEET)
good = et.check_datasheet(fx.COMPLETE_DATASHEET)
record("datasheet_incomplete_complete", "yes" if bad["complete"] else "no")
record("datasheet_missing", ",".join(bad["missing"]))
record("datasheet_missing_count", str(len(bad["missing"])))
record("datasheet_complete_complete", "yes" if good["complete"] else "no")
record("datasheet_required_field_count", str(len(et.DATASHEET_REQUIRED_FIELDS)))

# -- exercise 9: version drift ------------------------------------------------
v1, v2 = et.build_versions(n=2_000, seed=138)
diff = et.diff_versions(v1, v2)
record("drift_summary_identical", "yes" if diff["summary_identical"] else "no")
record("drift_mean_old", f"{diff['summary_old']['mean']:.10f}")
record("drift_mean_new", f"{diff['summary_new']['mean']:.10f}")
record("drift_std_old", f"{diff['summary_old']['std']:.10f}")
record("drift_std_new", f"{diff['summary_new']['std']:.10f}")
record("drift_share_b_old", f"{diff['group_shares_old']['B']:.4f}")
record("drift_share_b_new", f"{diff['group_shares_new']['B']:.4f}")
record("drift_composition_shift", f"{diff['composition_shift']:.4f}")
record("drift_summary_hides_the_change", "yes" if diff["summary_hides_the_change"] else "no")
record("drift_changelog_explains", "yes" if diff["changelog_explains"] else "no")
record("drift_new_changelog_entries", str(len(diff["changelog_new_entries"])))

for key, value in results.items():
    print(f"{key}={value}")
PY
)"
behaviour_status=$?
echo "${behaviour}"
echo

value_of() { echo "${behaviour}" | grep "^$1=" | cut -d= -f2-; }

check "the behaviour script ran without error" "$( [ ${behaviour_status} -eq 0 ] && echo yes || echo no )"

# -- exercise 1 --------------------------------------------------------------
check "bias for the under-represented group is flat across a hundredfold n: $(value_of bias_b_n500) -> $(value_of bias_b_n5000) -> $(value_of bias_b_n50000) (spread $(value_of bias_b_spread))" "$( [ "$(value_of bias_b_flat)" = yes ] && echo yes || echo no )"
check "every bias figure sits at the predicted 5.94, not near zero" "$( "${python_bin}" -c "import sys; vals=[float(v) for v in sys.argv[1:]]; sys.exit(0 if all(5.85 < v < 6.00 for v in vals) else 1)" "$(value_of bias_b_n500)" "$(value_of bias_b_n5000)" "$(value_of bias_b_n50000)" && echo yes || echo no )"
check "variance falls by more than 2.2x per tenfold: $(value_of sd_b_n500) -> $(value_of sd_b_n5000) -> $(value_of sd_b_n50000) (ratios $(value_of sd_ratio_500_5000), $(value_of sd_ratio_5000_50000))" "$( [ "$(value_of sd_falls_each_step)" = yes ] && echo yes || echo no )"
check "variance falls more than 5x overall (total ratio $(value_of sd_ratio_total)) while bias does not move at all" "$( "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) > 5.0 else 1)" "$(value_of sd_ratio_total)" && echo yes || echo no )"
check "the model looks fine at every n: in-sample RMSE $(value_of rmse_n500), $(value_of rmse_n5000), $(value_of rmse_n50000)" "$( "${python_bin}" -c "import sys; vals=[float(v) for v in sys.argv[1:]]; sys.exit(0 if all(1.0 < v < 1.4 for v in vals) else 1)" "$(value_of rmse_n500)" "$(value_of rmse_n5000)" "$(value_of rmse_n50000)" && echo yes || echo no )"
check "the group the frame does reach is fine: bias $(value_of bias_a_n500), $(value_of bias_a_n5000), $(value_of bias_a_n50000)" "$( "${python_bin}" -c "import sys; vals=[float(v) for v in sys.argv[1:]]; sys.exit(0 if all(v < 0.15 for v in vals) else 1)" "$(value_of bias_a_n500)" "$(value_of bias_a_n5000)" "$(value_of bias_a_n50000)" && echo yes || echo no )"
check "the sample's composition never converges on the population's: share_b stays near 0.01 ($(value_of share_b_n500), $(value_of share_b_n5000), $(value_of share_b_n50000))" "$( "${python_bin}" -c "import sys; vals=[float(v) for v in sys.argv[1:]]; sys.exit(0 if all(abs(v-0.01) < 0.005 for v in vals) else 1)" "$(value_of share_b_n500)" "$(value_of share_b_n5000)" "$(value_of share_b_n50000)" && echo yes || echo no )"

# -- exercise 2 --------------------------------------------------------------
check "coverage mismatch computed: total variation distance $(value_of coverage_tvd)" "$( [ "$(value_of coverage_tvd)" = "0.089583" ] && echo yes || echo no )"
check "the under-represented group is named: west appears at $(value_of coverage_ratio_west) of its population share" "$( "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) < 0.15 else 1)" "$(value_of coverage_ratio_west)" && echo yes || echo no )"
check "only the west is flagged (flagged: $(value_of coverage_flagged))" "$( [ "$(value_of coverage_flagged)" = "west" ] && [ "$(value_of coverage_worst_group)" = "west" ] && echo yes || echo no )"
check "a sample matching the reference scores exactly 0 and flags nobody" "$( [ "$(value_of coverage_clean_tvd)" = "0.000000" ] && [ "$(value_of coverage_clean_flagged_count)" = "0" ] && echo yes || echo no )"

# -- exercise 3 --------------------------------------------------------------
check "the proxy passes the usual check: correlation with the target is $(value_of proxy_correlation)" "$( "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) > 0.90 else 1)" "$(value_of proxy_correlation)" && echo yes || echo no )"
check "ranking by the target selects group B at $(value_of proxy_share_b_by_target); ranking by the proxy at $(value_of proxy_share_b_by_proxy)" "$( "${python_bin}" -c "import sys; t,p=float(sys.argv[1]),float(sys.argv[2]); sys.exit(0 if abs(t-0.25) < 0.05 and p < 0.10 else 1)" "$(value_of proxy_share_b_by_target)" "$(value_of proxy_share_b_by_proxy)" && echo yes || echo no )"
check "group B's true need served falls to $(value_of proxy_b_need_ratio) of what target-ranking would have served" "$( "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) < 0.35 else 1)" "$(value_of proxy_b_need_ratio)" && echo yes || echo no )"
check "and the aggregate barely moves -- total need served is $(value_of proxy_total_need_ratio) of the target-ranked total, which is why nobody notices" "$( "${python_bin}" -c "import sys; v=float(sys.argv[1]); sys.exit(0 if 0.97 < v < 1.0 else 1)" "$(value_of proxy_total_need_ratio)" && echo yes || echo no )"

# -- exercise 4 --------------------------------------------------------------
check "both subgroups trend downward (slopes $(value_of agg_slope_a), $(value_of agg_slope_b)) while the pooled fit trends up at $(value_of agg_pooled_slope)" "$( [ "$(value_of agg_sign_flip)" = yes ] && echo yes || echo no )"
check "pooled RMSE is worse for EVERY subgroup: A $(value_of agg_pooled_rmse_a) vs $(value_of agg_own_rmse_a), B $(value_of agg_pooled_rmse_b) vs $(value_of agg_own_rmse_b)" "$( [ "$(value_of agg_pooled_worse_for_every_group)" = yes ] && echo yes || echo no )"

# -- exercise 5 --------------------------------------------------------------
check "the score is perfectly calibrated for both groups (max deviation $(value_of fair_max_calibration_deviation))" "$( [ "$(value_of fair_max_calibration_deviation)" = "0.000000000000" ] && echo yes || echo no )"
check "the two groups have different base rates: A $(value_of fair_base_rate_a), B $(value_of fair_base_rate_b)" "$( [ "$(value_of fair_base_rate_a)" = "0.6600" ] && [ "$(value_of fair_base_rate_b)" = "0.3400" ] && echo yes || echo no )"
check "one threshold for all: parity gap $(value_of fair_single_threshold_parity_gap), equal-opportunity gap $(value_of fair_single_threshold_eo_gap), precision gap $(value_of fair_single_threshold_precision_gap)" "$( [ "$(value_of fair_single_threshold_parity_gap)" = "0.5000" ] && "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) > 0.30 and float(sys.argv[2]) > 0.10 else 1)" "$(value_of fair_single_threshold_eo_gap)" "$(value_of fair_single_threshold_precision_gap)" && echo yes || echo no )"
check "enforcing demographic parity closes its gap exactly ($(value_of fair_demographic_parity_parity_gap)) and opens the others ($(value_of fair_demographic_parity_eo_gap), $(value_of fair_demographic_parity_precision_gap))" "$( [ "$(value_of fair_demographic_parity_parity_gap)" = "0.0000" ] && "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) > 0.10 and float(sys.argv[2]) > 0.30 else 1)" "$(value_of fair_demographic_parity_eo_gap)" "$(value_of fair_demographic_parity_precision_gap)" && echo yes || echo no )"
check "enforcing equal opportunity closes its gap exactly ($(value_of fair_equal_opportunity_eo_gap)) and reopens the others ($(value_of fair_equal_opportunity_parity_gap), $(value_of fair_equal_opportunity_precision_gap))" "$( [ "$(value_of fair_equal_opportunity_eo_gap)" = "0.0000" ] && "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) > 0.10 and float(sys.argv[2]) > 0.25 else 1)" "$(value_of fair_equal_opportunity_parity_gap)" "$(value_of fair_equal_opportunity_precision_gap)" && echo yes || echo no )"
check "NO policy satisfies all three criteria -- the incompatibility, not a preference" "$( [ "$(value_of fair_any_policy_satisfies_all)" = no ] && echo yes || echo no )"

# -- exercise 6 --------------------------------------------------------------
check "$(value_of reid_unique_exact) of $(value_of reid_rows) rows are unique on birth year, postcode and sex ($(value_of reid_unique_exact_fraction) of the table)" "$( [ "$(value_of reid_unique_exact)" = "2723" ] && echo yes || echo no )"
check "coarsening ONE field cuts uniques to $(value_of reid_unique_generalised), a reduction factor of $(value_of reid_reduction_factor)" "$( [ "$(value_of reid_unique_generalised)" = "13" ] && "${python_bin}" -c "import sys; sys.exit(0 if float(sys.argv[1]) > 100 else 1)" "$(value_of reid_reduction_factor)" && echo yes || echo no )"

# -- exercise 7 --------------------------------------------------------------
check "generalisation alone does not buy k-anonymity (k is still $(value_of k_level_generalised))" "$( [ "$(value_of k_level_generalised)" = "1" ] && echo yes || echo no )"
check "suppressing classes under 5 achieves k=$(value_of k_level_after_suppression), at a cost of $(value_of k_rows_suppressed) rows ($(value_of k_rows_kept) kept)" "$( [ "$(value_of k_level_after_suppression)" -ge 5 ] && [ "$(value_of k_rows_suppressed)" = "794" ] && echo yes || echo no )"
check "a $(value_of k_homogeneous_table_level)-anonymous table still discloses: $(value_of k_leaking_class_count) class of size $(value_of k_leaking_class_size) with $(value_of k_leaking_distinct_values) distinct sensitive value, revealing '$(value_of k_disclosed_value)'" "$( [ "$(value_of k_homogeneous_table_level)" = "4" ] && [ "$(value_of k_leaking_class_count)" = "1" ] && [ "$(value_of k_leaking_distinct_values)" = "1" ] && [ "$(value_of k_disclosed_value)" = "diabetes" ] && echo yes || echo no )"

# -- exercise 8 --------------------------------------------------------------
check "an undocumented dataset fails the datasheet check, naming $(value_of datasheet_missing_count) fields: $(value_of datasheet_missing)" "$( [ "$(value_of datasheet_incomplete_complete)" = no ] && [ "$(value_of datasheet_missing_count)" = "8" ] && echo yes || echo no )"
check "a documented one passes, against all $(value_of datasheet_required_field_count) required provenance fields" "$( [ "$(value_of datasheet_complete_complete)" = yes ] && [ "$(value_of datasheet_required_field_count)" = "11" ] && echo yes || echo no )"

# -- exercise 9 --------------------------------------------------------------
check "two releases have identical summary statistics (mean $(value_of drift_mean_old) vs $(value_of drift_mean_new), std $(value_of drift_std_old) vs $(value_of drift_std_new))" "$( [ "$(value_of drift_summary_identical)" = yes ] && [ "$(value_of drift_mean_old)" = "$(value_of drift_mean_new)" ] && echo yes || echo no )"
check "and differ by $(value_of drift_composition_shift) in composition: group B goes from $(value_of drift_share_b_old) to $(value_of drift_share_b_new)" "$( [ "$(value_of drift_composition_shift)" = "0.2500" ] && [ "$(value_of drift_summary_hides_the_change)" = yes ] && echo yes || echo no )"
check "only the provenance record makes it visible: $(value_of drift_new_changelog_entries) new changelog entry naming the change" "$( [ "$(value_of drift_changelog_explains)" = yes ] && [ "$(value_of drift_new_changelog_entries)" = "1" ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "3. Reference suite -- examples/ must pass in full"
# --------------------------------------------------------------------------

examples_output="$(cd "${lab_dir}" && "${pytest_bin}" examples -q 2>&1)"
examples_status=$?
echo "${examples_output}" | tail -5
check "examples/ exits 0" "$( [ ${examples_status} -eq 0 ] && echo yes || echo no )"
examples_passed_line="$(echo "${examples_output}" | grep -E '^[0-9]+ passed' || true)"
check "examples/ reports 9 passed, 0 failed" "$( echo "${examples_passed_line}" | grep -qE '^9 passed' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "4. Exercise suite -- starter/ is all-skip on an untouched checkout"
# --------------------------------------------------------------------------

starter_output="$(cd "${lab_dir}" && "${pytest_bin}" starter -q 2>&1)"
starter_status=$?
echo "${starter_output}" | tail -5
check "starter/ (untouched) exits 0" "$( [ ${starter_status} -eq 0 ] && echo yes || echo no )"
check "starter/ (untouched) reports 9 skipped, 0 failed" "$( echo "${starter_output}" | grep -qE '^9 skipped' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "5. Never run 'pytest examples starter' in one invocation -- both"
echo "   directories define a module named test_ethics.py, and pytest"
echo "   collects by dotted module name. Documented, and run as two commands."
# --------------------------------------------------------------------------

combined_output="$(cd "${lab_dir}" && "${pytest_bin}" examples starter -q 2>&1)"
combined_status=$?
check "'pytest examples starter' aborts rather than silently passing" "$( [ ${combined_status} -ne 0 ] && echo yes || echo no )"
check "the collision is reported as an import file mismatch" "$( echo "${combined_output}" | grep -qi 'import file mismatch' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "6. Prove the suite can genuinely FAIL: solve every exercise in a"
echo "   scratch copy, confirm green, break one assertion on purpose,"
echo "   confirm a non-zero exit and a printed failure, then restore."
# --------------------------------------------------------------------------

scratch_dir="$(mktemp -d "${TMPDIR:-/tmp}/d138-scratch.XXXXXX")"
cleanup_scratch() { rm -rf "${scratch_dir}"; }
trap cleanup_scratch EXIT

for module in test_ethics.py ethics.py fixtures.py conftest.py; do
  cp "${lab_dir}/examples/${module}" "${scratch_dir}/${module}"
done

solved_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
solved_status=$?
check "scratch copy of the solved suite exits 0" "$( [ ${solved_status} -eq 0 ] && echo yes || echo no )"
check "scratch copy reports 9 passed" "$( echo "${solved_output}" | grep -qE '^9 passed' && echo yes || echo no )"

# Break exercise 6's exact uniqueness count on purpose.
sed -i.bak 's/assert exact == 2_723/assert exact == 9_999/' "${scratch_dir}/test_ethics.py"

broken_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
broken_status=$?
check "broken scratch copy exits non-zero" "$( [ ${broken_status} -ne 0 ] && echo yes || echo no )"
check "broken scratch copy prints a failure" "$( echo "${broken_output}" | grep -qiE 'failed|assert' && echo yes || echo no )"

mv "${scratch_dir}/test_ethics.py.bak" "${scratch_dir}/test_ethics.py"
restored_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
restored_status=$?
check "restored scratch copy exits 0 again" "$( [ ${restored_status} -eq 0 ] && echo yes || echo no )"
check "restored scratch copy reports 9 passed again" "$( echo "${restored_output}" | grep -qE '^9 passed' && echo yes || echo no )"

cleanup_scratch
trap - EXIT
echo

# --------------------------------------------------------------------------
echo "7. Offline, synthetic, and nothing left behind"
# --------------------------------------------------------------------------

network_hits="$(grep -rlE 'urllib|requests\.|socket|http://|https://' "${lab_dir}/examples" "${lab_dir}/starter" 2>/dev/null || true)"
check "no network call anywhere in the lab code -- everything is constructed locally" "$( [ -z "${network_hits}" ] && echo yes || echo no )"

seeded="$(grep -c 'default_rng' "${lab_dir}/examples/ethics.py" || true)"
check "every random draw goes through numpy default_rng with an explicit seed (${seeded} call sites)" "$( [ "${seeded}" -ge 4 ] && echo yes || echo no )"

find "${lab_dir}" -name '.venv' -prune -o -type d -name '__pycache__' -exec rm -rf {} + 2>/dev/null || true
find "${lab_dir}" -name '.venv' -prune -o -type d -name '.pytest_cache' -exec rm -rf {} + 2>/dev/null || true

stray="$(find "${lab_dir}" -name '.venv' -prune -o \( -type d -name '__pycache__' -print -o -type d -name '.pytest_cache' -print \) 2>/dev/null || true)"
check "no __pycache__ or .pytest_cache left behind" "$( [ -z "${stray}" ] && echo yes || echo no )"

leftover_tmp="$(find "${TMPDIR:-/tmp}" -maxdepth 1 -name 'd138-*' -print 2>/dev/null || true)"
check "no d138 temporary directory left in the system temp directory" "$( [ -z "${leftover_tmp}" ] && echo yes || echo no )"
echo

echo "-------------------------------------------------------------"
echo "${checks} checks, ${failures} failure(s)"
if [ "${failures}" -gt 0 ]; then
  exit 1
fi
exit 0