diff --git a/reports/tsi/README.md b/reports/tsi/README.md
index 01e59b7..33457a3 100644
--- a/reports/tsi/README.md
+++ b/reports/tsi/README.md
@@ -58,7 +58,7 @@
---
-> **The uncle-reference model analysed here.** Uncle references are *counting-only* — they never affect a block's validity ([§6.7](#s6-7)–[§6.8](#s6-8)) — and are deduplicated by **slot**, not by block ([§8.5](#s8-5)). The reference window is derived from a **window absorption parameter**, `w_u = W_abs/f` with `W_abs = 10` expected block-intervals and bound `W_abs ≤ 0.6·k`, which is the `W = 300 = 10/f` sizing rule of [§3.4](#s3-4). One structural constraint shapes the results: only the **first block of a fork** — one whose parent lies on the referencing chain — is countable, which is the price of verifying every reference from chain data alone. This report calls that the **countable** model, and measures it against an **unrestricted** baseline in which any orphan in the window is referenceable at any fork depth ([§2.1](#s2-1), [§3.2](#s3-2), [§9](#s9)). The result: **in the design regime `ρ < 1` the two are practically identical.** A dedicated 40-replicate **paired** sweep of the operating band (`δ_max` 1–5, `ρ ≈ 0.21–0.41`), in which both arms share the stake draw, the peering graph and every lottery outcome so the uncle rule is the only difference, puts every `U ≥ 1` cell of both models in 0.998–1.001 and pins the first-fork cost at **nothing measurable up to `δ_max` = 3 and ~0.1 % at `δ_max` = 4–5** ([§3.2a](#s3-2a)). Differences appear only at `ρ ≥ 1`: the restriction trims ~1.4 % at `ρ ≈ 1`, and under overload (`ρ ≈ 1.8`) recovery saturates near 0.95 where the unrestricted model reaches ~1.00. Deep forks are unrecoverable by construction, so past `ρ ≈ 1` the ceiling is set by the fork *depth* profile and **no uncle cap buys it back** — a second, independent reason for the [§8](#s8) "keep `ρ < 1`" rule. The closed-form accuracy with the measured effective utilization, `log(1−f)/log(1−f/q_u)`, predicts every countable cell to within 0.2 %.
+> **The uncle-reference model analysed here.** Uncle references are *counting-only* — they never affect a block's validity ([§6.7](#s6-7)–[§6.8](#s6-8)) — and are deduplicated by **slot**, not by block ([§8.5](#s8-5)). The reference window is derived from a **window absorption parameter**, `w_u = W_abs/f` with `W_abs = 10` expected block-intervals and bound `W_abs ≤ 0.6·k`, which is the `W = 300 = 10/f` sizing rule of [§3.4](#s3-4). One structural constraint shapes the results: only the **first block of a fork** — one whose parent lies on the referencing chain — is countable, which is the price of verifying every reference from chain data alone. This report calls that the **countable** model, and measures it against an **unrestricted** baseline in which any orphan in the window is referenceable at any fork depth ([§2.1](#s2-1), [§3.2](#s3-2), [§9](#s9)). The result: **in the design regime `ρ < 1` the two are practically identical.** A dedicated 40-replicate **paired** sweep of the operating band (`δ_max` 1–5, `ρ ≈ 0.21–0.41`), in which both arms share the stake draw, the peering graph and every lottery outcome so the uncle rule is the only difference, puts every `U ≥ 1` cell of both models in 0.998–1.001 and pins the first-fork cost at **nothing measurable up to `δ_max` = 3 and ~0.1 % at `δ_max` = 4–5** ([§3.2a](#s3-2a)). Measured against the unrestricted baseline under common random numbers, the cost grows monotonically with load — 0.13 % at `ρ ≈ 0.36`, 0.34 % at `ρ ≈ 0.56`, 1.0 % at `ρ ≈ 0.96` and 2.3 % under overload at `ρ ≈ 1.76`, where recovery saturates near 0.95 against the unrestricted model's ~1.00. Deep forks are unrecoverable by construction, so past `ρ ≈ 1` the ceiling is set by the fork *depth* profile and **no uncle cap buys it back** — a second, independent reason for the [§8](#s8) "keep `ρ < 1`" rule. The closed-form accuracy with the measured effective utilization, `log(1−f)/log(1−f/q_u)`, predicts every countable cell to within 0.2 %.
## 1. Executive summary
@@ -200,28 +200,30 @@ The U = 0 under-count **deepens with N** (more nodes → more concurrent proposa

-**What the first-fork restriction costs.** Running the delay × U grid under both referencing models on the same grid (`configs/countable-vs-old.yaml`, N = 1 000, blend, 3 hops, 5 replicates; [§9](#s9)) separates the regimes cleanly (`fig30`). Cells are replicate means ± SEM; `t` is the two-sample statistic for the countable−unrestricted gap, and cells with `t < 2` do **not** resolve a difference at this replicate count:
+**What the first-fork restriction costs.** The delay × U grid is run under both referencing models with **common random numbers** (`configs/countable-vs-old-paired.yaml`, 20 replicates, N = 1 000, blend, 3 hops; [§9](#s9)): a countable cell and its unrestricted twin share the stake draw, the peering graph and every lottery outcome, so the uncle rule is the only difference and each replicate is a matched pair. Cells are replicate means ± SEM; the gap and its `t` come from the 20 per-replicate differences.
-| `δ_max` (s) | ρ | U | countable | unrestricted | gap (t) |
-|---|---|---|---|---|---|
-| 4 | ≈ 0.36 | 1 / 2 / 4 | 0.9979 ± 0.0015 / 0.9982 ± 0.0014 / 0.9995 ± 0.0012 | 0.9981 ± 0.0006 / 0.9985 ± 0.0012 / 0.9992 ± 0.0013 | *none resolved* (0.1–0.2) |
-| 8 | ≈ 0.56 | 1 / 2 / 4 | 0.9965 ± 0.0015 / 0.9955 ± 0.0017 / 0.9963 ± 0.0015 | 0.9986 ± 0.0017 / 0.9993 ± 0.0016 / 0.9989 ± 0.0017 | *none resolved* (1.0–1.6) |
-| 16 | ≈ 0.96 | 1 | 0.9647 ± 0.0065 | 0.9610 ± 0.0043 | *not resolved* (0.5) |
-| 16 | ≈ 0.96 | 2 / 4 | 0.9862 ± 0.0010 / 0.9911 ± 0.0032 | 0.9999 ± 0.0005 / 1.0013 ± 0.0011 | −0.014 (12.1) / −0.010 (3.0) |
-| 32 | ≈ 1.76 | 1 | 0.608 ± 0.069 | 0.570 ± 0.042 | *not resolved* (0.5) |
-| 32 | ≈ 1.76 | 2 / 4 | 0.9343 ± 0.0015 / 0.9517 ± 0.0007 | 0.9774 ± 0.0064 / 1.0020 ± 0.0010 | −0.043 (6.5) / −0.050 (41.0) |
+| `δ_max` (s) | ρ | `U=0` (both) | `U=1` countable / unrestricted | `U=2` | `U=4` | pooled gap (t) |
+|---|---|---|---|---|---|---|
+| 4 | ≈ 0.36 | 0.674 | 0.9975 / 0.9989 | 0.9977 / 0.9991 | 0.9996 / 1.0005 | −0.0013 (4.0) |
+| 8 | ≈ 0.56 | 0.506 | 0.9962 / 0.9995 | 0.9960 / 0.9996 | 0.9958 / 0.9989 | −0.0034 (9.8) |
+| 16 | ≈ 0.96 | 0.298 | 0.9582 / 0.9610 | 0.9876 / 0.9998 | 0.9883 / 1.0002 | −0.0102 (17.5) |
+| 32 | ≈ 1.76 | 0.227 | 0.585 / 0.572 | 0.9394 / 0.9821 | 0.9539 / 1.0005 | −0.0228 (22.4) |
-**The `U = 0` rows are a negative control.** With no uncles the two models are identical by construction — no reference is ever taken, so any measured gap is pure between-run noise (the two models draw independent RNG streams by design, so the comparison is unpaired). That control reads **+0.010 (t = 0.8) at `δ_max = 4`** but **−0.230 (t = 2.1) at `δ_max = 32`**, where single replicates range from 0.05 to 0.60. Read the table against that floor: at `δ_max = 32` a gap must clear ~0.2 to mean anything on a single-replicate basis, which is exactly why the `U = 1` cells there are reported as unresolved while the tightly-clustered `U = 2` and `U = 4` cells are not.
+**The negative control is exact.** With `U = 0` no reference is ever taken, so under shared streams the two arms must be bit-identical rather than merely close: **all 80 replicate pairs differ by exactly 0.0**, at every delay.
-So: **at the operating loads (`ρ < 1`) no difference between the models is detectable at all** — the first-fork restriction costs nothing measurable where forks are shallow, and both models sit at ≥ 0.995 for every `U ≥ 1`. ([§3.2a](#s3-2a) puts a much tighter bound on this region.) Differences appear only at `ρ ≥ 1` and only at `U ≥ 2`: at `ρ ≈ 1` the restriction trims ~1.4 %, and in overload the countable recovery **saturates** — raising `U` 2 → 4 buys only 0.934 → 0.952 while the unrestricted model reaches ~1.00. The ceiling is the measured recovery rate `r` (0.92 at `δ_max = 32`, `U = 4`), set by the fork-*depth* profile rather than by capacity, so no uncle cap buys it back. The closed form `D̂/D = log(1−f)/log(1−f/q_u)` with the measured `q_u = q + (1−q)·r` predicts every countable cell to within 0.2 % (`fig31`).
+**The cost grows monotonically with load**, and is resolved at every delay: 0.13 % at `ρ ≈ 0.36`, 0.34 % at `ρ ≈ 0.56`, 1.0 % at `ρ ≈ 0.96`, 2.3 % at `ρ ≈ 1.76`. Eleven of the twelve `U ≥ 1` cells resolve individually (the exception is `U = 4` at `δ_max = 4`, `t = 1.07`), and the largest, `t = 29.0`, clears the Bonferroni threshold of 2.87 for twelve tests by an order of magnitude. At the design end the effect is small enough that it took pairing to see it at all — the earlier unpaired 5-replicate grid found nothing below `δ_max = 16` — and [§3.2a](#s3-2a) reaches the same `δ_max = 4` figure (−0.0013 at 20 replicates here, −0.0011 at 40 there) from an independent sweep.
-`fig32` shows `r` and the residual `1−r` per delay, rising from 0.3 % to 7.8 % of the waste as `δ_max` goes 4 → 32 at `U = 4`. That residual is the *combined* first-fork and capacity loss — the simulator cannot separate the two, because the restriction acts at selection (a deep orphan is simply never referenced) rather than at counting. At `U = 4`, where capacity is not binding for `ρ ≤ 1.8`, most of it is attributable to fork depth; the split is not measured directly.
+**At `U = 1` under overload the sign flips**, and this is the one place the countable rule *wins*: at `δ_max = 32` it reads 0.585 against 0.572, a gap of **+0.0127 (t = 7.6), positive in 19 of 20 pairs**. Both models have collapsed there — one uncle slot cannot drain `ρ ≈ 1.76` — but when capacity is the binding constraint, the countable rule's occupied-slot exclusion means its single reference always recovers a *new* slot, whereas the unrestricted rule deduplicates by block id and can spend its one reference on an orphan whose slot is already counted, recovering nothing. The measured recovery agrees: `q_u` = 0.591 countable against 0.579 unrestricted. So the slot-vs-block distinction of [§2.1](#s2-1) is worth most exactly where references are scarcest.
-
+**Everywhere else the restriction dominates**, and it bites hardest through the fork *depth* profile rather than capacity: raising `U` from 2 to 4 at `δ_max = 32` moves the countable model only 0.939 → 0.954 while the unrestricted model reaches ~1.00. The ceiling is the measured recovery rate `r` (0.92 at `δ_max = 32`, `U = 4`), so no uncle cap buys it back. The closed form `D̂/D = log(1−f)/log(1−f/q_u)` with the measured `q_u = q + (1−q)·r` predicts every countable cell to within 0.2 % (`fig31`).
+
+`fig32` shows `r` and the residual `1−r` per delay, rising from 0.3 % to 8.1 % of the waste as `δ_max` goes 4 → 32 at `U = 4`. That residual is the *combined* first-fork and capacity loss — the simulator cannot separate the two, because the restriction acts at selection (a deep orphan is simply never referenced) rather than at counting.
+
+

-
+
**Bootstrap is self-limiting — but only with uncles** (`fig1`; full scale, k = 2160, Blend, N = 1 000 and 5 000, genesis guesses 0.01×–2× the true stake; 0.1×–2× at N = 5 000). With uncles (U = 2, solid lines) the cold start is a non-event: whatever the guess, block production snaps back to the target `f` within ~2–2.5 epochs and the estimate lands on 1.0 by epoch 2 — identically at both network sizes. Without uncles (U = 0, dashed) the system also converges from any guess, **but to the wrong place**: the estimate settles below truth (mean ≈ 0.58× across guesses, ranging ≈0.48–0.85; ≈0.51× at N = 5 000) and the network then *chronically over-produces* blocks at ~1.85× the target rate (~2.0× at N = 5 000) — the [§6.2](#s6-2) load feedback in the flesh (a low estimate makes the lottery easier, extra blocks orphan, the count stays low). So uncle references are load-bearing from the very first epochs: they are what makes the bootstrap end *at the truth* rather than at a permanently overheated equilibrium.
@@ -879,7 +881,8 @@ Sweep studies are committed configs, run with `make ` (writes a dated `run
| soft uncle inclusion (reward share vs emergent `p_ref`) | `scripts/reward_mandate.py` | [§6.8](#s6-8), [§6.9](#s6-9) |
| U = 0 fluctuation series (zero delay, k ∈ {256, 1024, 2160}) | `scripts/appendix_fluct.py --run` (`runs/fluctuation_u0.parquet`) | [Appendix B](#sB) |
| CI smoke grid + analytic sanity checks | `configs/smoke.yaml`; `scripts/verify.py` (`make verify`) — validation only, no figures | — |
-| **countable vs unrestricted referencing** (accuracy over delay × U; measured `q_u`/recovery `r`) | `configs/countable-vs-old.yaml` run twice — default and with `--old`; figures + significance table by `scripts/plot_countable_vs_old.py` | [§2.1](#s2-1), [§3.2](#s3-2) |
+| **countable vs unrestricted referencing, PAIRED** (accuracy over delay × U; measured `q_u`/recovery `r`) | `configs/countable-vs-old-paired.yaml` run twice — default and with `--old`, 20 replicates, common random numbers; figures + paired-difference table by `scripts/plot_countable_vs_old.py` | [§2.1](#s2-1), [§3.2](#s3-2) |
+| countable vs unrestricted, unpaired (same grid, 5 replicates; superseded — retained as the power comparison) | `configs/countable-vs-old.yaml` run twice | [§3.2](#s3-2) |
| **fine delay band, PAIRED** (δ_max 1–5, 40 replicates, common random numbers — the primary design-regime measurement) | `configs/fine-delay-paired.yaml` run twice — default and with `--old`; paired-difference tests + figures by `scripts/plot_fine_delay.py` | [§3.2a](#s3-2a) |
| fine delay band, unpaired (same grid without `paired_streams`; superseded, retained as the power comparison) | `configs/fine-delay.yaml` run twice | [§3.2a](#s3-2a) |
| **window absorption sweep** (`W` in expected block-intervals, `w_u = W/f` derived) | `configs/absorption-window.yaml`; figure by `scripts/plot_countable_vs_old.py` | [§3.4](#s3-4) |
@@ -888,7 +891,7 @@ Sweep studies are committed configs, run with `make ` (writes a dated `run
Because the two models draw independent RNG streams, every countable-vs-unrestricted comparison is **unpaired**, and its resolution is set by the replicate spread rather than by the effect size. Each comparison sweep therefore includes a `U = 0` arm as a **negative control**: with no uncles the models are identical by construction, so the measured `U = 0` gap is a direct reading of the noise floor at that delay and replicate count. At `δ_max = 32` with 5 replicates that floor is ≈ 0.23 in `D̂/D` — larger than several real effects elsewhere in the grid — which is why [§3.2](#s3-2) reports a `t` statistic per cell and why the design regime is measured separately at 40 replicates ([§3.2a](#s3-2a)).
-All studies were **re-run on 2026-07-23/24 with the corrected slot-counting mechanism** ([§2.1](#s2-1)) and the early-stop optimisation; the resilient batch is `scripts/run_all_reruns.sh` (per-step log in `runs/rerun_status.log`). Canonical run directories (latest): fullscale N=5000/10000 = `2026-07-24_094519_fullscale`; fullscale N=1000/2000 = `2026-07-23_171803_fullscale-small`; uncle-window = `2026-07-24_001456`; window-uncles = `2026-07-24_014240`; block-rate = `2026-07-24_043943`; blend-hops-delay = `2026-07-24_064052`; window-scale = `2026-07-24_085234`; latency-shape = `2026-07-24_090014_expdist`; stake-tail = `2026-07-24_090044_pareto133`; heterogeneous-start = `2026-07-24_090114_default`; N-scaling = the `nscaling-{a,b}` + `nscaling32-{a,b}` runs; adversary grids = `runs/adversary_grid/`; jitter = `runs/jitter_grid/`; bootstrap = `runs/bootstrap_fullscale/`; fluctuation = `runs/fluctuation_u0.parquet`; fork-rate = `runs/fork_rate_vs_delay.parquet`; ρ-boundary = `2026-07-27_195627_rho-boundary`. Referencing-model studies: countable-vs-unrestricted = `2026-08-04_103536_cvo-countable` / `2026-08-04_104010_cvo-old`; window absorption = `2026-08-04_104633_absorption-window`; fine delay band (unpaired) = `2026-08-04_191441_fine-countable` / `2026-08-04_195353_fine-old`; fine delay band (paired, primary) = `2026-08-05_111109_fine-paired-countable` / `2026-08-05_113218_fine-paired-old`.
+All studies were **re-run on 2026-07-23/24 with the corrected slot-counting mechanism** ([§2.1](#s2-1)) and the early-stop optimisation; the resilient batch is `scripts/run_all_reruns.sh` (per-step log in `runs/rerun_status.log`). Canonical run directories (latest): fullscale N=5000/10000 = `2026-07-24_094519_fullscale`; fullscale N=1000/2000 = `2026-07-23_171803_fullscale-small`; uncle-window = `2026-07-24_001456`; window-uncles = `2026-07-24_014240`; block-rate = `2026-07-24_043943`; blend-hops-delay = `2026-07-24_064052`; window-scale = `2026-07-24_085234`; latency-shape = `2026-07-24_090014_expdist`; stake-tail = `2026-07-24_090044_pareto133`; heterogeneous-start = `2026-07-24_090114_default`; N-scaling = the `nscaling-{a,b}` + `nscaling32-{a,b}` runs; adversary grids = `runs/adversary_grid/`; jitter = `runs/jitter_grid/`; bootstrap = `runs/bootstrap_fullscale/`; fluctuation = `runs/fluctuation_u0.parquet`; fork-rate = `runs/fork_rate_vs_delay.parquet`; ρ-boundary = `2026-07-27_195627_rho-boundary`. Referencing-model studies: countable-vs-unrestricted (paired, primary) = `2026-08-05_120804_cvo-paired-countable` / `2026-08-05_123403_cvo-paired-old`; the superseded unpaired pair = `2026-08-04_103536_cvo-countable` / `2026-08-04_104010_cvo-old`; window absorption = `2026-08-04_104633_absorption-window`; fine delay band (unpaired) = `2026-08-04_191441_fine-countable` / `2026-08-04_195353_fine-old`; fine delay band (paired, primary) = `2026-08-05_111109_fine-paired-countable` / `2026-08-05_113218_fine-paired-old`.
Figures are embedded from [`report-figures/`](report-figures) and versioned here alongside the report; the simulator folder does **not** commit its own generated figures, so the copies checked in here are the figures of record. Figures are in `report-figures/` (`fig1`–`fig29`, plus [Appendix B](#sB)'s `figB1`–`figB2`; numbering is generation order, not order of appearance). Committed generators: `fig1` (bootstrap, k=2160) by `scripts/bootstrap_dynamics.py`; `fig2`,`fig4`,`fig5`,`fig17`–`fig22` by `scripts/regenerate_extra_figs.py` from the latest sweeps (`fig3` hops×delay×U grid by `scripts/hops_delay_grid.py`, `fig6` (block-rate `U_min` grid + ρ-collapse) rendered ad hoc from `runs/2026-07-24_043943_block-rate` with no committed generator; `fig26` deficit-vs-ρ by `scripts/rho_boundary_analysis.py`) (fullscale-derived `fig17`–`fig20` pool both sizes in that run, N = 5 000 and N = 10 000 — the generators filter on stake_dist/topology/degree/init_dest only, never on `n_nodes`); `fig8`,`fig9` by `scripts/adversary_figs.py` from `runs/adversary_grid/`; `fig10`–`fig12` by `scripts/dynamic_withhold.py`; `fig13`–`fig15` by `scripts/selfish_mining.py`/`selfish_rewards.py`/`reward_mandate.py`; `fig16` by `scripts/stake_vs_delay.py`; `fig23`–`fig24` by `scripts/nscaling_analysis.py`; `fig25` by `scripts/window_scale_analysis.py`; `fig27`–`fig28` by `scripts/reorg_depth.py` (fork rates via `--measure`; private-chain model `src/tsi_sim/reorg.py`); `fig29` by `scripts/churn.py`; `figB1`–`figB2` by `scripts/appendix_fluct.py`. `fig7` (feedback fixed-point) is an analytic overlay; `fig30`–`fig33` (countable-vs-unrestricted accuracy, `q_u`-prediction check, recovery rate, absorption-window sweep) by `scripts/plot_countable_vs_old.py` from the `cvo-countable`/`cvo-old`/`absorption-window` runs; `fig34`–`fig35` (design-regime accuracy and the model gap with 95 % CIs) by `scripts/plot_fine_delay.py` from the `fine-countable`/`fine-old` runs. Every figure type the per-node simulator generates appears in this report, and the fork-rate/reorg-depth study closes the previous reproducibility gap for the adversarial figures (`fig8`,`fig9` now have committed generators from `runs/adversary_grid/`).
diff --git a/reports/tsi/report-figures/fig30_countable_vs_old.png b/reports/tsi/report-figures/fig30_countable_vs_old.png
index b187b75..6d0bb0e 100644
Binary files a/reports/tsi/report-figures/fig30_countable_vs_old.png and b/reports/tsi/report-figures/fig30_countable_vs_old.png differ
diff --git a/reports/tsi/report-figures/fig31_countable_prediction.png b/reports/tsi/report-figures/fig31_countable_prediction.png
index 33b8c0f..3826c51 100644
Binary files a/reports/tsi/report-figures/fig31_countable_prediction.png and b/reports/tsi/report-figures/fig31_countable_prediction.png differ
diff --git a/reports/tsi/report-figures/fig32_countable_recovery.png b/reports/tsi/report-figures/fig32_countable_recovery.png
index 951ce4f..73a3e81 100644
Binary files a/reports/tsi/report-figures/fig32_countable_recovery.png and b/reports/tsi/report-figures/fig32_countable_recovery.png differ
diff --git a/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_countable_vs_old.py b/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_countable_vs_old.py
index db3273a..84e60a1 100644
--- a/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_countable_vs_old.py
+++ b/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_countable_vs_old.py
@@ -34,7 +34,14 @@ import numpy as np
import pandas as pd
from tsi_sim.plotting import style
-from tsi_sim.plotting.figures_pernode import equilibrium, recovery_rate, sem
+from tsi_sim.plotting.figures_pernode import (
+ Z95,
+ equilibrium,
+ paired_gaps,
+ pooled_by_delay,
+ recovery_rate,
+ sem,
+)
from tsi_sim.theory import expected_ratio, q_effective
DELAY = "blend_delay_max"
@@ -61,7 +68,8 @@ def _eq(df: pd.DataFrame, extra_cols: tuple[str, ...] = ()) -> pd.DataFrame:
return eq.groupby(keys, as_index=False).agg(**agg)
-def fig_accuracy_vs_delay(cnt: pd.DataFrame, old: pd.DataFrame) -> plt.Figure:
+def fig_accuracy_vs_delay(cnt: pd.DataFrame, old: pd.DataFrame,
+ paired: bool = False) -> plt.Figure:
"""Accuracy vs delay per U, with replicate SEM bars.
U=0 is the NEGATIVE CONTROL: with no uncles the two models are identical by
@@ -82,7 +90,8 @@ def fig_accuracy_vs_delay(cnt: pd.DataFrame, old: pd.DataFrame) -> plt.Figure:
ax.axhline(1.0, color="0.4", lw=0.8, ls=":")
ax.set_xlabel("max per-relay mixing delay (slots)")
ax.set_ylabel(r"equilibrium $\hat{D}/D_{true}$")
- ax.set_title("Accuracy vs delay: countable (solid) vs unrestricted (dashed) referencing")
+ ax.set_title("Accuracy vs delay: countable (solid) vs unrestricted (dashed)"
+ + (" — paired" if paired else ""))
ax.legend(ncol=2, fontsize="x-small")
return fig
@@ -171,20 +180,26 @@ def main() -> None:
cnt_raw, old_raw = _load(args.countable), _load(args.old)
f = float(cnt_raw["f"].iloc[0])
cnt, old = _eq(cnt_raw), _eq(old_raw)
+ # When both arms shared their RNG streams the comparison is PAIRED: the per-replicate
+ # difference cancels the shared variance, so it supersedes the unpaired two-sample test.
+ pg = paired_gaps(cnt_raw, old_raw)
prov = "tsi-sim-pernode countable-vs-old.yaml (+--old) / absorption-window.yaml"
written = []
- written += style.save(fig_accuracy_vs_delay(cnt, old), out / "cvo_accuracy_vs_delay", prov)
+ written += style.save(fig_accuracy_vs_delay(cnt, old, pg is not None),
+ out / "cvo_accuracy_vs_delay", prov)
written += style.save(fig_prediction_vs_sim(cnt, f), out / "cvo_prediction_vs_sim", prov)
written += style.save(fig_recovery_vs_delay(cnt), out / "cvo_recovery_vs_delay", prov)
written += style.save(fig_absorption_window(_load(args.absorption)),
out / "absorption_window", prov)
- # Headline numbers for the report. Every countable-vs-old gap is printed with the
- # two-sample t = |diff| / SE(diff) over replicates; |t| < 2 means the cell does NOT
- # resolve a model difference at this replicate count and must not be read as one.
- # The U=0 rows are the negative control (models identical by construction, so their
- # t is a pure noise reading).
- print(f"{'cell':<14} {'countable':>17} {'old':>17} {'diff':>9} {'t':>6} verdict")
+ # Headline numbers for the report. The gap is printed with t = |gap| / SE(gap); |t| < 2
+ # means the cell does NOT resolve a model difference at this replicate count and must not
+ # be read as one. Under paired_streams the SE is the per-replicate paired one (the shared
+ # variance cancels); otherwise it is the unpaired two-sample SE. The U=0 rows are the
+ # negative control — identical models, so paired they must be EXACTLY zero.
+ kind = "PAIRED (common random numbers)" if pg is not None else "unpaired two-sample"
+ print(f"\ncountable vs unrestricted — {kind}, {int(cnt.n_rep.min())} replicates/arm")
+ print(f"{'cell':<14} {'countable':>17} {'unrestricted':>17} {'gap':>9} {'t':>6} verdict")
for u in sorted(cnt["max_uncles"].unique()):
for _, row in cnt[cnt.max_uncles == u].sort_values(DELAY).iterrows():
q, qu = row.mean_q, row.mean_q_eff
@@ -193,10 +208,15 @@ def main() -> None:
if not len(o):
continue
orow = o.iloc[0]
- diff = row.mean_ratio - orow.mean_ratio
- se = float(np.hypot(row.sem_ratio, orow.sem_ratio))
- t = abs(diff) / se if se > 0 else float("inf")
- verdict = ("CONTROL (must be 0)" if u == 0 else
+ if pg is not None:
+ pr = pg[(pg.max_uncles == u) & (pg[DELAY] == row[DELAY])].iloc[0]
+ diff, se, t = float(pr.gap), float(pr.se), float(pr.t)
+ else:
+ diff = row.mean_ratio - orow.mean_ratio
+ se = float(np.hypot(row.sem_ratio, orow.sem_ratio))
+ t = abs(diff) / se if se > 0 else float("inf")
+ verdict = ("CONTROL (exactly 0)" if u == 0 and diff == 0.0 else
+ "CONTROL" if u == 0 else
"resolved" if t >= 2 else "NOT RESOLVED (noise)")
print(f"U={u} delay={row[DELAY]:>5g} "
f"{row.mean_ratio:.4f}+-{row.sem_ratio:.4f} "
@@ -205,6 +225,16 @@ def main() -> None:
print(f"{'':>14} q={q:.4f} q_u={qu:.4f} r={r:.4f} "
f"pred={float(expected_ratio(f, q_effective(q, r))):.4f} "
f"n_rep={int(row.n_rep)}")
+ if pg is not None:
+ print("\npooled over U>=1 per delay (the caps measure the same difference):")
+ for _, r_ in pooled_by_delay(pg).iterrows():
+ mark = " <-- resolved" if r_.t >= 2 else ""
+ print(f" delay={r_[DELAY]:>5g}: {r_.gap:+.5f} +-{Z95 * r_.se:.5f} "
+ f"t={r_.t:6.2f}{mark}")
+ ctl = pg[pg.max_uncles == 0]
+ ex, tot = int(ctl.n_zero.sum()), int(ctl.n_pair.sum())
+ print(f" U=0 control must be exactly 0 under pairing: {ex}/{tot} pairs "
+ f"-> {'PASSES' if ex == tot else 'FAILS'}")
print(f"wrote {len(written)} files -> {out}")
diff --git a/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_fine_delay.py b/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_fine_delay.py
index 3750fb5..3c86a06 100644
--- a/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_fine_delay.py
+++ b/tools/simulators/tsi/tsi-sim-pernode/scripts/plot_fine_delay.py
@@ -28,52 +28,29 @@ from __future__ import annotations
import argparse
from pathlib import Path
+from statistics import NormalDist
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from tsi_sim.plotting import style
-from tsi_sim.plotting.figures_pernode import equilibrium, rho_for, sem
+from tsi_sim.plotting.figures_pernode import (
+ Z95,
+ equilibrium,
+ paired_gaps,
+ pooled_by_delay,
+ rho_for,
+ sem,
+)
DELAY = "blend_delay_max"
-# Normal approximation: with 40 replicates per arm the t-quantile is within ~2% of 1.96,
-# and the replicate spread itself is the dominant uncertainty, so 1.96 is precise enough.
-Z95 = 1.96
def _load(run_dir: str | Path) -> pd.DataFrame:
return pd.read_parquet(Path(run_dir) / "results.parquet")
-def paired_gaps(cnt_raw: pd.DataFrame, old_raw: pd.DataFrame) -> pd.DataFrame | None:
- """Per-replicate differences, when both arms were run with ``paired_streams``.
-
- Under common random numbers replicate *i* of each arm shares the stake draw, the peering
- graph and the lottery outcomes, so ``d_i = countable_i - unrestricted_i`` is a PAIRED
- observation and the shared variance cancels. The test is then a one-sample t on the d_i,
- which is what makes a sub-0.1 % effect reachable per cell instead of only after pooling.
-
- Returns None when the runs are not paired, so the caller falls back to the unpaired test.
- """
- if not (cnt_raw.get("paired_streams", pd.Series([False])).all()
- and old_raw.get("paired_streams", pd.Series([False])).all()):
- return None
- c, o = equilibrium(cnt_raw), equilibrium(old_raw)
- keys = [DELAY, "max_uncles", "replicate"]
- m = c[[*keys, "mean_ratio"]].merge(o[[*keys, "mean_ratio"]], on=keys,
- suffixes=("_c", "_o"))
- m["d"] = m.mean_ratio_c - m.mean_ratio_o
- rows = []
- for (dl, u), s in m.groupby([DELAY, "max_uncles"]):
- d = s.d.to_numpy()
- se = sem(d)
- rows.append({DELAY: dl, "max_uncles": u, "gap": float(d.mean()), "se": se,
- "ci95": Z95 * se, "t": abs(d.mean()) / se if se > 0 else np.inf,
- "n_pair": len(d), "n_zero": int((d == 0.0).sum())})
- return pd.DataFrame(rows).sort_values(["max_uncles", DELAY])
-
-
def _cells(df: pd.DataFrame) -> pd.DataFrame:
"""Per (delay, U): replicate mean, SEM and count of the equilibrium accuracy."""
return equilibrium(df).groupby([DELAY, "max_uncles"], as_index=False).agg(
@@ -127,24 +104,6 @@ def fig_accuracy(cnt: pd.DataFrame, old: pd.DataFrame) -> plt.Figure:
return fig
-def pooled_by_delay(g: pd.DataFrame) -> pd.DataFrame:
- """Inverse-variance pooled gap across the U >= 1 arms, per delay.
-
- Individual cells are underpowered against a sub-0.1 % effect even at 40 replicates, but
- the three uncle caps are independent measurements of the same underlying difference, so
- pooling them buys back a factor of ~sqrt(3) and is what actually resolves the trend.
- """
- u = g[g.max_uncles > 0]
- rows = []
- for d, s in u.groupby(DELAY):
- w = 1.0 / s.se.to_numpy() ** 2
- p = float((s.gap.to_numpy() * w).sum() / w.sum())
- e = float(np.sqrt(1.0 / w.sum()))
- rows.append({DELAY: d, "gap": p, "se": e, "ci95": Z95 * e,
- "t": abs(p) / e if e > 0 else np.inf})
- return pd.DataFrame(rows).sort_values(DELAY)
-
-
def fig_gap(g: pd.DataFrame) -> plt.Figure:
"""The countable − unrestricted gap with 95% CIs, zoomed to the U >= 1 scale.
@@ -232,9 +191,10 @@ def main() -> None:
print(f"\nn_rep = {int(g.n_rep_c.min())}/{int(g.n_rep_o.min())} per arm")
print(f"widest 95% CI half-width at U>=1: +-{worst.ci95.max():.4f} "
f"({100 * worst.ci95.max():.2f} pp)")
- # Per-cell significance must be read against the number of cells tested: with 15 cells,
- # ~0.75 are expected to clear t=2 by chance alone, so quote the Bonferroni threshold.
- bonf = 2.935 if len(worst) == 15 else float("nan")
+ # Per-cell significance must be read against the number of cells tested — some are
+ # expected to clear t=2 by chance alone — so quote the Bonferroni threshold for the grid
+ # actually run, not a constant baked in for one grid size.
+ bonf = NormalDist().inv_cdf(1 - 0.05 / (2 * len(worst)))
print(f"per-cell: {int((worst.t >= 2).sum())}/{len(worst)} cells with t>=2 "
f"(expected by chance {0.05 * len(worst):.2f}); max t = {worst.t.max():.2f} "
f"vs Bonferroni threshold {bonf:.3f}")
@@ -264,8 +224,10 @@ def main() -> None:
print(f" U=0 control under pairing must be EXACTLY zero: {exact}/{tot} pairs are 0.0"
f" -> {'PASSES' if exact == tot else 'FAILS — streams are not shared'}")
res = pg[(pg.max_uncles > 0) & (pg.t >= 2)]
- print(f" per-cell resolved at |t|>=2: {len(res)}/{len(pg[pg.max_uncles>0])}"
- f" (unpaired: {int((unpaired[unpaired.max_uncles>0].t>=2).sum())}/15)")
+ n_cells = len(pg[pg.max_uncles > 0])
+ n_unpaired = int((unpaired[unpaired.max_uncles > 0].t >= 2).sum())
+ print(f" per-cell resolved at |t|>=2: {len(res)}/{n_cells}"
+ f" (same data, unpaired test: {n_unpaired}/{n_cells})")
ctl = g[g.max_uncles == 0]
if len(ctl):
diff --git a/tools/simulators/tsi/tsi-sim-pernode/src/tsi_sim/plotting/figures_pernode.py b/tools/simulators/tsi/tsi-sim-pernode/src/tsi_sim/plotting/figures_pernode.py
index 7843d86..92398bf 100644
--- a/tools/simulators/tsi/tsi-sim-pernode/src/tsi_sim/plotting/figures_pernode.py
+++ b/tools/simulators/tsi/tsi-sim-pernode/src/tsi_sim/plotting/figures_pernode.py
@@ -83,6 +83,58 @@ def rho_for(df: pd.DataFrame, delay) -> np.ndarray:
return f * (hops * np.asarray(delay, dtype=float) / 2.0 + (hops + 1) * ell)
+DELAY = "blend_delay_max"
+# Normal approximation: at 20+ replicates the t-quantile is within a few percent of
+# 1.96, and the replicate spread dominates, so 1.96 is precise enough for these CIs.
+Z95 = 1.96
+
+
+def paired_gaps(cnt_raw: pd.DataFrame, old_raw: pd.DataFrame) -> pd.DataFrame | None:
+ """Per-replicate differences, when both arms were run with ``paired_streams``.
+
+ Under common random numbers replicate *i* of each arm shares the stake draw, the peering
+ graph and the lottery outcomes, so ``d_i = countable_i - unrestricted_i`` is a PAIRED
+ observation and the shared variance cancels. The test is then a one-sample t on the d_i,
+ which is what makes a sub-0.1 % effect reachable per cell instead of only after pooling.
+
+ Returns None when the runs are not paired, so the caller falls back to the unpaired test.
+ """
+ if not (cnt_raw.get("paired_streams", pd.Series([False])).all()
+ and old_raw.get("paired_streams", pd.Series([False])).all()):
+ return None
+ c, o = equilibrium(cnt_raw), equilibrium(old_raw)
+ keys = ["blend_delay_max", "max_uncles", "replicate"]
+ m = c[[*keys, "mean_ratio"]].merge(o[[*keys, "mean_ratio"]], on=keys,
+ suffixes=("_c", "_o"))
+ m["d"] = m.mean_ratio_c - m.mean_ratio_o
+ rows = []
+ for (dl, u), s in m.groupby(["blend_delay_max", "max_uncles"]):
+ d = s.d.to_numpy()
+ se = sem(d)
+ rows.append({"blend_delay_max": dl, "max_uncles": u, "gap": float(d.mean()), "se": se,
+ "ci95": Z95 * se, "t": abs(d.mean()) / se if se > 0 else np.inf,
+ "n_pair": len(d), "n_zero": int((d == 0.0).sum())})
+ return pd.DataFrame(rows).sort_values(["max_uncles", "blend_delay_max"])
+
+
+def pooled_by_delay(g: pd.DataFrame) -> pd.DataFrame:
+ """Inverse-variance pooled gap across the U >= 1 arms, per delay.
+
+ Individual cells are underpowered against a sub-0.1 % effect even at 40 replicates, but
+ the three uncle caps are independent measurements of the same underlying difference, so
+ pooling them buys back a factor of ~sqrt(3) and is what actually resolves the trend.
+ """
+ u = g[g.max_uncles > 0]
+ rows = []
+ for d, s in u.groupby(DELAY):
+ w = 1.0 / s.se.to_numpy() ** 2
+ p = float((s.gap.to_numpy() * w).sum() / w.sum())
+ e = float(np.sqrt(1.0 / w.sum()))
+ rows.append({DELAY: d, "gap": p, "se": e, "ci95": Z95 * e,
+ "t": abs(p) / e if e > 0 else np.inf})
+ return pd.DataFrame(rows).sort_values(DELAY)
+
+
def _lat_axis(topo: str) -> tuple[str, str]:
"""(dataframe column, axis label) for the dominant latency knob of a graph topology."""
if topo == "blend":