Two-universe calibration (ONS frame + HMRC registered subset), full regeneration, and paper revision

.github/workflows/paper.yml

@@ -66,6 +66,13 @@ jobs:
           cd paper
           quarto render index.qmd --to pdf
 
+      - name: Upload rendered PDF
+        uses: actions/upload-artifact@v4
+        with:
+          name: paper-pdf
+          path: paper/out/*.pdf
+          if-no-files-found: error
+
       - name: Render paper HTML
         run: |
           cd paper

Makefile

@@ -1,4 +1,4 @@
-.PHONY: sync check test lint paper
+.PHONY: sync check test lint paper etl reproduce reproduce-heavy figures-sync
 
 sync:
 	uv sync --extra dev --locked
@@ -13,3 +13,39 @@ check: test lint
 
 paper:
 	cd paper && TEXMFVAR=/tmp/firm-microsim-tex-cache latexmk -lualatex -interaction=nonstopmode -halt-on-error main.tex
+
+# ---------------------------------------------------------------------------
+# Reproduction (issue #40). `reproduce` regenerates both synthetic vintages and
+# every checked results/*.txt and figure the paper reads, then syncs figures
+# into paper/figures/. `reproduce-heavy` adds the two multi-build artifacts
+# (seed sensitivity = two extra full builds; placebo B = one extra build).
+# The PDF is rendered by the Paper workflow in CI (quarto is not required
+# locally); see paper/README.md.
+# ---------------------------------------------------------------------------
+etl:
+	uv run python scripts/etl_ons_tables.py
+
+reproduce: etl
+	uv run firm-microsim --seed 42
+	uv run firm-microsim-static
+	uv run python analysis/static_results_dump.py
+	uv run firm-microsim-bunching --figures-only
+	uv run python scripts/bunching_inference.py
+	uv run firm-microsim-notch
+	uv run firm-microsim-dynamic
+	uv run firm-microsim-reform-menu
+	uv run firm-microsim-dominated-region
+	uv run firm-microsim-verify-optimum
+	uv run firm-microsim-formulation-a-optima
+	uv run python analysis/recovery_bunching.py
+	$(MAKE) figures-sync
+
+reproduce-heavy: reproduce
+	uv run python scripts/seed_sensitivity.py --seeds 42 7 99
+	uv run firm-microsim-placebo
+
+figures-sync:
+	@for f in paper/figures/*.png; do \
+	  b=$$(basename $$f); \
+	  if [ -f results/$$b ]; then cp results/$$b paper/figures/$$b; else echo "no results/$$b"; fi; \
+	done

README.md

@@ -41,15 +41,27 @@ threshold:
    and intermediate inputs. These within-band values are explicit modelling
    assumptions, not recovered administrative microdata.
 2. **Calibrate firm weights** by multi-objective optimisation (Adam, symmetric
-   relative-error loss) so weighted totals reproduce the official targets — HMRC
-   VAT-registered counts by turnover band and by sector, ONS employment-band
-   totals, and HMRC VAT-liability totals — with turnover bands weighted most
-   heavily. VAT registration is then assigned: mandatory above the threshold,
-   voluntary below at the HMRC-calibrated rate.
-
-The result is ~2.94M firm rows weighted to ~2.5M UK firms. Because the population
-is calibrated **to** the HMRC aggregates, agreement with them is an internal
-consistency check, not external validation.
+   relative-error loss) on **two declared universes** (issue #37):
+   - the **ONS VAT/PAYE enterprise frame** (~2.72M enterprises): the population
+     total, the employment-band totals and the OBR near-threshold shape;
+   - the **HMRC VAT-registered subset** (~2.18M traders): trader counts by
+     turnover band and by trade sector, and net VAT liability by band. Each
+     frame firm enters these rows with a *registration propensity* — the HMRC
+     count in its turnover band divided by the frame's mass in that band
+     (0.89 below the threshold; 0.51–1.00 above it up to £10m, the remainder
+     being PAYE-only or exempt-sector enterprises outside the VAT net; 0.64 in
+     the >£10m band, where the shortfall is mostly the uniform £5m–£50m draw for
+     the open ONS "5000+" band placing too many rows above £10m — issue #40). HMRC
+     negative/zero-turnover traders are appended before calibration as an
+     out-of-frame stratum.
+
+   VAT **scope** and **registration** flags are then assigned by seeded weighted
+   selection per band so registered totals match HMRC to within one weight.
+
+The result is ~2.94M firm rows (2.72M frame draws plus 0.22M appended traders)
+weighted to ~2.72M frame enterprises, of which ~2.18M are VAT-registered.
+Because the population is calibrated **to** the HMRC aggregates, agreement
+with them is an internal consistency check, not external validation.
 
 The official target surface is also being mirrored into PolicyEngine Ledger and
 Populace. This repository keeps the paper's archived CSV inputs and generator for
@@ -102,7 +114,15 @@ df, report = firm_microsim.generate(return_report=True)
 
 Output is written to `data/synthetic/synthetic_firms.csv`
 (`sic_code, annual_turnover_k, annual_input_k, vat_liability_k, employment,
-weight, vat_registered`).
+weight, vat_scope, vat_registered, in_frame`). `in_frame` marks ONS-frame
+enterprises; `vat_scope` marks firms in the VAT net (registered above the
+threshold, registrable below it); `vat_registered` marks HMRC-count-matched
+traders. Static and behavioural costings use in-scope firms only.
+
+The processed ONS tables are rebuilt from the raw workbooks by
+`scripts/etl_ons_tables.py` (Table 8 for turnover sizebands, Table 3 for
+employment sizebands, both counting enterprises); `--check` fails if a checked
+CSV differs from the extract, and `tests/test_data_universes.py` runs it.
 
 ## Paper
 
@@ -143,7 +163,8 @@ the manuscript without changing the computational results.
 
 The population is calibrated to **five** official ONS + HMRC target groups; the
 validator scores each dimension as
-`accuracy = max(0, 1 − |synthetic − target| / |target|)`. The displayed error is
+`accuracy = max(0, 1 − |synthetic − target| / |target|)`, HMRC dimensions on
+VAT-registered rows and ONS dimensions on frame rows. The displayed error is
 the clipped complement of that score, not a signed relative error. **Overall**
 is the simple mean over the five calibrated dimensions below.
 Reproduce with:
@@ -154,12 +175,12 @@ firm-microsim-report
 
 | Calibrated dimension | 85k (2023-24) | 90k (2024-25) |
 | --- | ---: | ---: |
-| HMRC turnover bands | 97.3% | 96.6% |
-| ONS population | 88.3% | 91.1% |
-| Employment bands | 92.1% | 92.3% |
-| Sector distribution | 90.1% | 90.4% |
-| VAT liability by band (6 calibrated bands) | 96.5% | 97.1% |
-| **Overall (5 calibrated dimensions)** | **92.9%** | **93.5%** |
+| HMRC turnover bands (registered rows; hold by construction) | 100.0% | 100.0% |
+| ONS population (frame rows) | 99.7% | 99.8% |
+| Employment bands (frame rows) | 99.8% | 99.9% |
+| Sector distribution (registered rows) | 93.3% | 92.2% |
+| VAT liability by band (6 calibrated bands, registered rows) | 97.1% | 98.0% |
+| **Overall (5 calibrated dimensions)** | **98.0%** | **98.0%** |
 
 **VAT liability by *sector*** is **not** a calibration target — it is reported as
 an informational diagnostic only, and neither is the **below-threshold
@@ -190,8 +211,9 @@ remains a true total. A vintage builds in ~15 seconds instead of ~13 minutes;
 headline aggregates reproduce the full build within ~0.3% and local bunching
 statistics within ~5%. Use for development only — release artifacts are
 full-size. Generator-seed sensitivity of the full build is recorded in
-`results/seed_sensitivity.txt` (E ±111, raise ±£1.4m / taper ±£0.3m across
-seeds; reproduce with `scripts/seed_sensitivity.py`).
+`results/seed_sensitivity.txt` (E ±166, raise ±£1.1m / taper ±£1.1m, base
+±£0.45bn across seeds — the scope flag is a seeded draw; reproduce with
+`scripts/seed_sensitivity.py`).
 
 ## Populace/Ledger migration check
 
@@ -229,7 +251,7 @@ targets and the paper's processed 2024-25 numeric inputs: six normalized source
 tables checked, zero mismatches, max numeric difference 0. It does **not** exactly
 replicate the paper's generated synthetic population: Populace's shared optimizer
 landed at 93.8% overall accuracy under its own validator versus the paper's
-then-90.5% (2024-25 scores 93.5% on the corrected build,
+then-90.5% (2024-25 scores 98.0% on the current two-universe build,
 `results/calibration_accuracy.txt`), but that overall pair is **not
 like-for-like**: HMRC turnover-band accuracy uses different band sets, and sector
 distribution reflects different calibration-target definitions. The directly
@@ -281,24 +303,27 @@ firm-microsim-static          # -> results/{vat_threshold_revenue_impact,revenue
 
 - `vat_threshold_revenue_impact.png` — the £85k→£90k anchor reform vs HMRC's
   published costing, by fiscal year. **Built on the £85k / 2023-24 vintage** —
-  the pre-reform basis HMRC actually had at the 6 March 2024 costing (the
-  threshold was still £85k until 1 April 2024). Full-deregistration model
-  −357/−364/−224/−78/+119 vs HMRC −150/−185/−125/−50/+65 £m (43% voluntary
-  retention: −203/−207/−128/−44/+68); both turn positive by 2028-29. The gap
-  shows that registration dynamics and population scope matter. See
-  `results/static_sweep.txt`.
+  the pre-reform basis HMRC actually had at the 6 March 2024 costing. Each year
+  differences revenue under the £90k policy and under the counterfactual
+  threshold path (85/85/87/89/92k), releasing in-scope firms only below the
+  £88k deregistration threshold: −250/−234/−200/−65/+123 vs HMRC
+  −150/−185/−125/−50/+65 £m; with 43% voluntary retention −142/−134/−114/−37/+123.
+  HMRC's figures lie between the two conventions in every release year. See
+  `results/static_sweep.txt` for the whole-band and fixed-preference
+  sensitivities.
 - `revenue_impact_2025_26.png` / `firms_impact_2025_26.png` — the forward static
   sweep of registration thresholds (£70k–£120k) vs the current £90k baseline,
   **on the £90k / 2024-25 vintage**.
 
-**Two vintages, two exercises.** The anchor reform uses the £85k vintage, where
-the affected `[85,90)k` band sits *above* the £85k registration threshold and is
-cleanly populated with registered firms — so a **simple band-sum** suffices (no
-de-bunching). The forward sweep uses the current £90k vintage; there the
-`[85,90)k` firms are *below* threshold and the calibration concentrates weight
-on them, so the sweep instead fits the clean above-threshold firm/liability
-profile and extrapolates it across the threshold
-(`StaticVATModel._counterfactual_bins`, unaged turnover scaled to the fiscal
-year by a nominal-growth factor). Revenue and the anchor reform match the paper
-and HMRC closely; the forward-sweep firm-count magnitudes run low because the
-regenerated population has a lower near-threshold VAT-paying-firm density.
+**Two vintages, two exercises, one convention.** The anchor reform uses the
+£85k vintage, where the affected `[85,90)k` band sits *above* the £85k
+registration threshold and is cleanly populated with in-scope registered firms,
+so a direct band-sum suffices. The forward sweep uses the current £90k vintage
+and is likewise a direct mechanical reclassification of in-scope firms
+(`StaticVATModel.threshold_sweep`); voluntary registrants below the data-year
+threshold stay registered, firms aged across it are released by a rise unless
+protected by the £2k deregistration gap, and out-of-scope enterprises never
+remit. `results/static_sweep.txt` decomposes each anchor year's release. In both exercises turnover **and** liability are aged to the
+fiscal year by the same nominal-growth factor (relative to each vintage's own
+data year), so band membership is evaluated on aged turnover — the fiscal-drag
+convention. `results/static_sweep.txt` holds the machine-readable table.

analysis/static_results_dump.py

@@ -11,7 +11,7 @@
 from __future__ import annotations
 
 from firm_microsim.config import RESULTS_DIR
-from firm_microsim.static.model import StaticVATModel, SWEEP_THRESHOLDS
+from firm_microsim.static.model import FISCAL_YEARS, StaticVATModel, SWEEP_THRESHOLDS
 
 LLAT_VOLUNTARY_SHARE = 0.43  # Liu-Lockwood-Almunia-Tam (2021): ~43% below-threshold
 
@@ -33,18 +33,47 @@ def main() -> None:
     anchor = anchor_model.anchor_reform()
     W(anchor.to_string(index=False))
     W("")
+    W("Release decomposition (headline; liability GBPm of firms registered under the")
+    W("baseline threshold but not under GBP90k, by data-year status):")
+    for fy in FISCAL_YEARS:
+        df = anchor_model._aged(anchor_model._growth(fy["year"]))
+        base_t, pol_t = float(fy["baseline"]), float(fy["policy"])
+        rel = anchor_model._registered(df, base_t) & ~anchor_model._registered(df, pol_t)
+        lw = df["liab"] * df["weight"]
+        vol = float(lw[rel & df["voluntary"]].sum()) / 1e6
+        man = float(lw[rel & df["mandatory"]].sum()) / 1e6
+        nev = float(lw[rel & ~df["voluntary"] & ~df["mandatory"]].sum()) / 1e6
+        W(f"  {fy['year']}: released {float(df['weight'][rel].sum()):,.0f} firms; "
+          f"voluntary-at-data-year {vol:.1f}m, mandatory-at-data-year {man:.1f}m, "
+          f"not-registered-at-data-year {nev:.1f}m")
+    W("")
+    W("Deregistration-threshold sensitivity: the headline releases registered")
+    W("firms only below the GBP88k deregistration threshold ([85k, 88k) in the")
+    W("raise years). Releasing the whole [baseline, 90k) band instead (gap = 0):")
+    for _, row in anchor_model.anchor_reform(gap=0.0).iterrows():
+        W(f"  {row['year']}: whole-band release {float(row['policyengine_impact_m']):+,.1f}m")
+    W("")
     W("Voluntary-retention sensitivity (anchor, per year): headline assumes")
     W("every released firm deregisters (full liability lost). If the Liu et")
     W(f"al. (2021) voluntary share ({LLAT_VOLUNTARY_SHARE:.0%}) of released-firm liability is")
     W("retained, the impact scales accordingly:")
-    for _, row in anchor.iterrows():
-        impact = float(row["policyengine_impact_m"])
-        W(f"  {row['year']}: headline {impact:+,.1f}m -> retention-adjusted "
-          f"{impact * (1 - LLAT_VOLUNTARY_SHARE):+,.1f}m "
+    ret = anchor_model.anchor_reform(retention=LLAT_VOLUNTARY_SHARE)
+    for (_, row), (_, r2) in zip(anchor.iterrows(), ret.iterrows()):
+        W(f"  {row['year']}: headline {float(row['policyengine_impact_m']):+,.1f}m -> "
+          f"retention-adjusted {float(r2['policyengine_impact_m']):+,.1f}m "
           f"(HMRC {float(row['hmrc_impact_m']):+,.0f}m)")
     W("")
+    W("Fixed-preference sensitivity: baseline voluntary registrants keep their")
+    W("registration wherever the threshold moves (treats the frame's ~89% below-")
+    W("threshold registered share as revealed preference; frame-selection bias")
+    W("makes this a lower bound on the revenue loss):")
+    for _, row in anchor_model.anchor_reform(retain_voluntary=True).iterrows():
+        W(f"  {row['year']}: fixed-preference {float(row['policyengine_impact_m']):+,.1f}m")
+    W("")
     W("Threshold sweep (2024-25 vintage, GBP 90k baseline, 2025-26 fiscal year)")
-    W("method: direct mechanical reclassification on the calibrated population")
+    W("method: direct mechanical reclassification of in-scope VAT firms; voluntary")
+    W("registrants below the data-year threshold held registered, firms aged across")
+    W("it released unless gap-protected; turnover and liability aged together")
     W("-" * 74)
     sweep = sweep_model.threshold_sweep(year="2025-26")
     W(sweep.to_string(index=False))
@@ -57,10 +86,10 @@ def main() -> None:
     W("-" * 74)
     for vintage, model in (("2023-24", anchor_model), ("2024-25", sweep_model)):
         df = model._aged(1.0)
-        thr = 85_000.0 if vintage == "2023-24" else 90_000.0
-        reg = df["turnover"] >= thr
-        base = float((df.loc[reg, "liab"] * df.loc[reg, "weight"]).sum()) / 1e9
-        W(f"  {vintage}: registered base (>= threshold, unaged) = {base:.1f}bn")
+        base = model._mandatory_base(df, model.data_threshold) / 1e9
+        vol = (model._revenue(df, model.data_threshold) - base * 1e9) / 1e9
+        W(f"  {vintage}: in-scope base (>= threshold, unaged) = {base:.1f}bn; "
+          f"model-implied voluntary below-threshold remittance = {vol:.1f}bn (not calibrated)")
     W("")
     W(f"Sweep thresholds: {SWEEP_THRESHOLDS}")
 

data/README.md

@@ -49,7 +49,19 @@ threshold automatically. The default is `2023-24`.
 
 The raw workbooks are the authoritative inputs. Everything in `processed/` is a
 faithful extract of tables inside those workbooks (no modelling applied), kept as
-CSV so the package can read them without an Excel engine.
+CSV so the package can read them without an Excel engine. The two ONS tables per
+vintage are produced by `scripts/etl_ons_tables.py` (run `--check` to verify the
+checked CSVs against the workbooks); the HMRC tables are manual extracts of the
+named sheets described below.
+
+### Statistical units and universes
+
+Both ONS tables count **enterprises** (Table 8: turnover sizebands; Table 3:
+employment sizebands), so their per-SIC totals agree to ONS disclosure rounding.
+An earlier 2023-24 employment CSV reproduced Table 18 — **local units**, i.e.
+individual sites (3,173,660 against 2,724,775 enterprises) — and was replaced in
+issue #37. The generator treats the ONS frame as the modelled population and the
+HMRC trader tables as its VAT-registered subset (see the top-level README).
 
 ### 2024-25 band re-aggregation note
 
@@ -98,8 +110,9 @@ liability cells are negative (net repayment positions).
 
 ### `ons_firm_turnover.csv`
 
-Count of enterprises by SIC sector and **turnover band (£ thousands)**, from the
-ONS UK Business 2024 workbook.
+Count of enterprises by SIC division and **turnover band (£ thousands)**, from
+**Table 8** (United Kingdom columns) of the ONS UK Business workbook for the
+vintage.
 
 - Columns: `SIC Code`, `Description`, then one count column per turnover band,
   then `Total`.
@@ -115,8 +128,9 @@ ONS UK Business 2024 workbook.
 
 ### `ons_firm_employment.csv`
 
-Count of enterprises by SIC sector and **employment size band (number of
-employees)**, from the ONS UK Business 2024 workbook.
+Count of enterprises by SIC division and **employment size band (number of
+employees)**, from **Table 3** (United Kingdom columns) of the ONS UK Business
+workbook for the vintage.
 
 - Columns: `SIC Code`, `Description`, then one count column per employment band,
   then `Total`.