Correct the #982 receipts: real pool at frame weight vs PUF-clone half ×2

changelog.d/982-real-pool-bands.fixed.md

@@ -0,0 +1 @@
+Corrected the #982 experiment README: it now labels its PUF-clone-half ×2 tables as distribution checks, and a new section scores the real pool at frame weight, which before calibration reaches about half of SOI at $5M+ proxy AGI after #1033.

experiments/982-puf-self-prediction/README.md

@@ -1,10 +1,13 @@
 # PUF imputation self-prediction (#982): receipts
 
-Issue: PolicyEngine/microcosm#982. Branch `us-982-puf-demographic-rank-predictors`.
-Scripts in this folder: `compare_predictors.py` (fit and predict per predictor
-variant) and `score_variants.py` (score against SOI). Scored output:
-`results/scores_8_trees.json`; input counts and per-run metadata (trees, rows,
-peak memory): `results/inputs_and_runs.json`.
+Issue: PolicyEngine/microcosm#982. Fix merged in #1033 (branch
+`us-982-puf-demographic-rank-predictors`). Scripts in this folder:
+`compare_predictors.py` (fit and predict per predictor variant),
+`score_variants.py` (score the PUF-clone half ×2 against SOI) and
+`actual_pool_bands.py` (score the real pool, both halves at frame weight).
+Scored output: `results/scores_8_trees.json` and
+`results/actual_pool_bands.json`; input counts and per-run metadata (trees,
+rows, peak memory): `results/inputs_and_runs.json`.
 
 ## The defect
 
@@ -16,7 +19,9 @@ variant below), PUF-clone wages have rank correlation 1.000 with the same
 unit's survey wages; among units with positive survey wages, 99.98% land within
 10% of their survey value; and no clone exceeds the survey maximum.
 Survey topcodes and underreporting therefore carried into the PUF half, which
-is where the missing $5M+ tail came from (#958).
+is why that half lacked the $5M+ tail (#958). The survey half, which carries
+the other half of the weight, has no units there at all; see
+[Real pool at frame weight](#real-pool-at-frame-weight).
 
 ## The fix
 
@@ -39,12 +44,19 @@ middle.
 
 All numbers below are pre-calibration, from one fit on the saved 2026-09-12
 donor frame (211,677 PUF tax units) predicted onto the 231,007 saved PUF-clone
-recipients, weights = PUF-clone half x2. AGI is a proxy: the sum of the 15
-imputed income items. Run with 8 trees per variant (memory was constrained);
-every variant uses the same configuration. SOI references: Table 1.1 TY2023
-counts and AGI; Table 1.4 shares.
+recipients, weights = PUF-clone half ×2 ("as if this half were the whole
+population"). AGI is a proxy: the sum of the 15 imputed income items. Run with
+8 trees per variant (memory was constrained); every variant uses the same
+configuration. SOI references: Table 1.1 TY2023 counts and AGI; Table 1.4
+shares.
 
-### Returns by size of proxy AGI
+The ×2 weighting compares the variants' distributions with one another and
+with SOI's shape. It is not a level comparison with SOI: in the real pool the
+survey half carries the other 50% of the weight and has no units at $5M+ proxy
+AGI. For levels, see
+[Real pool at frame weight](#real-pool-at-frame-weight).
+
+### Returns by size of proxy AGI: PUF-clone half ×2 (distribution check, not a level comparison)
 
 | Band | Old design | Rank only | Rank + flag (chosen) | Within-group rank + flag | SOI 1.1 |
 |---|---|---|---|---|---|
@@ -54,10 +66,10 @@ counts and AGI; Table 1.4 shares.
 | $5M-10M | 23,600 | 47,582 | 48,629 | 49,569 | 49,262 |
 | $10M+ | 8,186 | 31,404 | 33,920 | 29,041 | 30,382 |
 
-AGI in $10M+: old $185B, chosen $1,037B, SOI $908B. Capital gains share of
-$10M+ AGI: old 0.1%, chosen 42.6%, SOI Table 1.4 39.5%.
+AGI in $10M+ (×2): old $185B, chosen $1,037B, SOI $908B. Capital gains share
+of $10M+ AGI: old 0.1%, chosen 42.6%, SOI Table 1.4 39.5%.
 
-### Wages: self-prediction and participation
+### Wage self-prediction and participation: PUF-clone half ×2 (distribution check, not a level comparison)
 
 | Measure | Old | Rank only | Rank + flag | Within-group |
 |---|---|---|---|---|
@@ -74,16 +86,67 @@ ranks level with low but positive PUF incomes; rank alone handed positive wages
 to 47.9% (weighted) of survey units with no wages or self-employment income. Why pooled rather than within-group ranks: ranking
 earners only among earners matched survey earners to a PUF earner distribution
 that includes many small separate returns (for example dependents filing their
-own), cutting imputed wages to $7.93T; the pooled rank keeps $9.68T with the
-same participation fidelity.
+own), cutting imputed wages to $7.93T (×2); the pooled rank keeps $9.68T (×2)
+with the same participation fidelity.
+
+## Real pool at frame weight
+
+The PUF-detail expansion (`microcosm.build.us_runtime.puf_support`) splits
+each unit's weight evenly between two clone channels. Clone index 0, the
+native survey channel, keeps the survey values; clone index 1
+(`PUF_TAX_DETAIL_CLONE_INDEX`) receives the imputed ones.
+`actual_pool_bands.py` combines both halves at frame weight and writes
+`actual_pool_bands.json` into its input directory;
+`results/actual_pool_bands.json` holds the recorded copy. The figures come
+from the same fit as above (8 trees, proxy AGI), before calibration and before
+main's capital-gains tail stage (#567). The fit also predates #1251 (the
+donor's partnership and S-corporation income moved to Schedule E line 32) and
+#1259 (a bound on each QRF draw), and these figures reflect neither. "Old
+design" is the `current` variant; "#1033" is `demographic_midrank6_earn`, the
+rank + flag design that #1033 merged.
+
+Clone 0 carries 50.0% of the tax-unit weight (`clone0_weight_share`) and has
+no units at $5M+ proxy AGI in either variant: in both $5M+ bands,
+`clone1_returns` equals `returns`. Its proxy AGI sums the 11 of the 15 items
+that its person table carries; the 15-item list's estate, farm, partnership
+and S-corporation columns are absent (`clone0_items_used`).
+
+| Band | Old design | #1033 | SOI Table 1.1 |
+|---|---|---|---|
+| $1M–1.5M returns | 701k (1.90×) | 483k (1.31×) | 369k |
+| $5M–10M returns | 11.8k (0.24×) | 24.3k (0.49×) | 49.3k |
+| $10M+ returns | 4.1k (0.13×) | 17.0k (0.56×) | 30.4k |
+| $10M+ AGI | $92B (0.10×) | $519B (0.57×) | $908B |
+| Colorado income above $1M | $20.2B (0.64×; largest record 28%) | $66.3B (2.1×; largest record 81%) | $31.7B (Historic Table 2, via #940) |
+
+"Largest record" is the share of the state's income above $1M carried by its
+single largest contributor.
+
+#1033 nearly quadruples (3.9×) the pool's proxy AGI at $5M+ ($174B to $675B)
+and lifts its $5M+ returns from 15.9k to 41.3k. On its own, though, it
+reaches about half of SOI there (52% of the 79.6k returns and 54% of the $1,244B), not
+SOI. At $5M+ the real-pool figures are exactly half of the ×2 ones (for
+example 48,629 to 24,315 returns at $5M–10M), because only the PUF-clone half
+has units there. The bands in between ($1.5M–5M) sit at 0.59–0.84× SOI in
+returns and 0.52–0.83× in AGI across the two designs. At $1M–1.5M the survey half alone holds 284k returns, 77% of
+SOI's count, in both designs. All five bands, with AGI and record counts, are
+in `results/actual_pool_bands.json`.
+
+These figures correct the ×2 level comparisons first posted on #982, #940 and
+#1033; each of those threads carries a correction.
 
 ## Known limits
 
-- Colorado (#940) is not fixed by this change. The chosen variant has 15
-  Colorado records at $1M+ proxy AGI; the largest (survey six-item total
-  $1.1M, weight 867) draws $125.4M of proxy AGI and carries 84% of the state's
-  income above $1M (`colorado_1m_plus` in the scores), so state top-tail totals
-  are noisy until state x AGI-band amount targets (#940) enter calibration.
+- Colorado (#940) is not fixed by this change. In the PUF-clone half the
+  chosen variant has 15 Colorado records at $1M+ proxy AGI. The largest
+  (survey six-item total $1.1M) draws $125.4M of proxy AGI at frame weight
+  434, shown as 867 in the ×2 scores. In the real pool at frame weight,
+  Colorado's income above $1M is $66.3B, 2.1× SOI's $31.7B, and that one
+  record carries 81% of it. The ×2 scores put the total at $128.5B, about 4×
+  SOI, with that record at 84% (`colorado_1m_plus`); those are not levels.
+  #1102 later added state × AGI-band return and amount targets (#940) to the
+  US calibration target surface; this pre-calibration run predates them and
+  shows nothing about a build that binds them.
 - All variants ran at 8 trees because memory was constrained; production fits
   more trees. The comparison is like for like across variants.
 - The within-10% and wage-total rows are pre-calibration; calibration targets

experiments/982-puf-self-prediction/actual_pool_bands.py

@@ -0,0 +1,154 @@
+"""Real-pool proxy-AGI bands at frame weight (microcosm#982 correction).
+
+Usage::
+
+    python actual_pool_bands.py <dir with imputed_<variant>.parquet>
+
+``score_variants.py`` weights the PUF-clone half x2 ("as if this half were the
+whole population"). The real pool is the survey half (clone index 0, survey
+values, half the weight) plus the PUF-clone half (clone index 1, imputed
+values, the other half), before main's capital-gains tail stage. This script
+combines both halves at frame weight on the same proxy AGI and compares the
+bands with SOI Table 1.1. The survey half sums whichever of the 15 proxy-AGI
+items its person table carries; ``clone0_items_used`` in the output lists them.
+
+Writes ``actual_pool_bands.json`` into the input directory and prints a
+summary. Inputs are restricted build artifacts on the author's machine; the
+checkpoint path is recorded as provenance, not as a portable interface.
+"""
+
+from __future__ import annotations
+
+import json
+import os
+import sys
+import warnings
+
+import numpy as np
+import pandas as pd
+from score_variants import COLORADO_SOI_INCOME_ABOVE_1M_BN, INCOME, SOI_TABLE_1_1
+
+from microcosm.build.frame_checkpoint import load_frame_checkpoint
+
+warnings.filterwarnings("ignore")
+
+CHECKPOINTS = os.path.expanduser(
+    "~/spm-rebuild-20260908/rollout/codex-continuation-20260912/"
+    "genuine-build-76325eb7-20260912/base-checkpoints"
+)
+VARIANTS = ("current", "demographic_midrank6_earn")
+COLORADO_FIPS = 8
+
+
+def survey_half(frame) -> tuple[pd.DataFrame, float, list[str]]:
+    """Clone index 0 at frame weight: survey proxy AGI, weight and state per unit.
+
+    Also returns clone 0's share of all tax-unit weight and the income items
+    found on the person table.
+    """
+
+    person, tax_unit = frame.table("person"), frame.table("tax_unit")
+    household = frame.table("household").copy()
+    household["w"] = frame.resolve_weights("household").values
+    clone = tax_unit["tax_unit_support_clone_index"].to_numpy()
+    have = [column for column in INCOME if column in person.columns]
+    survey = (
+        person[have]
+        .astype(float)
+        .sum(axis=1)
+        .groupby(person["person_tax_unit_id"])
+        .sum()
+    )
+    unit_household = (
+        person.groupby("person_tax_unit_id")["person_household_id"]
+        .first()
+        .reindex(tax_unit["tax_unit_id"])
+        .to_numpy()
+    )
+    by_household = household.set_index("household_id")
+    weight = by_household["w"].reindex(unit_household).to_numpy()
+    state = by_household["state_fips"].reindex(unit_household).to_numpy()
+    zero = clone == 0
+    half = pd.DataFrame(
+        {
+            "agi": survey.reindex(tax_unit["tax_unit_id"][zero]).fillna(0).to_numpy(),
+            "w": weight[zero],
+            "state": state[zero],
+            "clone": 0,
+        }
+    )
+    return half, float(weight[zero].sum() / weight.sum()), have
+
+
+def pool_report(pool: pd.DataFrame) -> dict[str, object]:
+    """SOI Table 1.1 bands and Colorado's income above $1M for one pool."""
+
+    bands = {}
+    for band, (low, high, returns, soi_agi) in SOI_TABLE_1_1.items():
+        mask = (pool.agi >= low) & (pool.agi < high)
+        bands[band] = {
+            "returns": float(pool.w[mask].sum()),
+            "soi_returns": returns,
+            "agi_bn": float((pool.w * pool.agi)[mask].sum() / 1e9),
+            "soi_agi_bn": soi_agi,
+            "records": int(mask.sum()),
+            "clone1_returns": float(pool.w[mask & (pool.clone == 1)].sum()),
+        }
+    colorado = pool.state == COLORADO_FIPS
+    above = pool.w * np.maximum(pool.agi - 1e6, 0)
+    return {
+        "bands": bands,
+        "colorado_income_above_1m_bn": float(above[colorado].sum() / 1e9),
+        "colorado_soi_bn": COLORADO_SOI_INCOME_ABOVE_1M_BN,
+        "colorado_largest_record_share": float(
+            above[colorado].max() / above[colorado].sum()
+        ),
+    }
+
+
+def main() -> None:
+    directory = sys.argv[1]
+    frame = load_frame_checkpoint(
+        f"{CHECKPOINTS}/002_clone_feature_extraction.frame.h5"
+    ).frame
+    survey, survey_share, items = survey_half(frame)
+    out = {
+        "clone0_weight_share": survey_share,
+        "clone0_items_used": items,
+        "variants": {},
+    }
+    for variant in VARIANTS:
+        imputed = pd.read_parquet(f"{directory}/imputed_{variant}.parquet")
+        clone1 = pd.DataFrame(
+            {
+                "agi": imputed[INCOME].sum(axis=1).to_numpy(),
+                "w": imputed["w"].to_numpy(),
+                "state": imputed["state_fips"].to_numpy(),
+                "clone": 1,
+            }
+        )
+        pool = pd.concat([survey, clone1], ignore_index=True)
+        out["variants"][variant] = pool_report(pool)
+    with open(f"{directory}/actual_pool_bands.json", "w") as handle:
+        json.dump(out, handle, indent=1)
+
+    print(f"clone0 weight share {out['clone0_weight_share']:.3f}")
+    for variant, report in out["variants"].items():
+        print("==", variant)
+        for band, x in report["bands"].items():
+            print(
+                f"  {band:12} returns {x['returns']:>9,.0f} vs SOI "
+                f"{x['soi_returns']:>8,} ({x['returns'] / x['soi_returns']:.2f}x)  "
+                f"AGI ${x['agi_bn']:6.0f}B vs ${x['soi_agi_bn']:4.0f}B  "
+                f"records {x['records']:5}  clone1 {x['clone1_returns']:>9,.0f}"
+            )
+        print(
+            "  Colorado income above $1M "
+            f"${report['colorado_income_above_1m_bn']:.1f}B vs SOI "
+            f"${report['colorado_soi_bn']:.1f}B; largest record share "
+            f"{report['colorado_largest_record_share']:.0%}"
+        )
+
+
+if __name__ == "__main__":
+    main()