1"""F1-conf-01: 12-month CPIAUCSL growth from two ALFRED vintages.2Run from workspace root: ~/wp29/venv/bin/python out/analysis.py3"""4import json, os5import numpy as np6import pandas as pd78D, O = "data", "out"9V16, V26, ASM = "cpi_vintage_2016-07-15.csv", "cpi_vintage_2026-10-28.csv", "cpi_assembled.csv"10W0, W1 = "2001-01-01", "2026-09-01"1112def load(f):13 d = pd.read_csv(os.path.join(D, f), dtype={"date": str})14 return d1516a, b, asm = load(V16), load(V26), load(ASM)17n_rows = {V16: len(a), V26: len(b), ASM: len(asm)}18issues = []1920# --- validation of vintage files21for name, d in ((V16, a), (V26, b)):22 assert not d.date.duplicated().any(), name23 dt = pd.to_datetime(d.date)24 assert (dt.diff().dropna().dt.days <= 31).all(), f"gap in {name}"25missing_src = b[b.cpi.isna()].date.tolist()26for m in missing_src:27 issues.append({"issue": f"Vintage 2026-10-28 reports no CPI value for {m[:7]} (empty cell; the source publishes no observation for that month).",28 "handling": "Not interpolated or imputed (rule). Not back-filled from the older vintage, which does not cover that month anyway. g_t for that month is reported missing; g_{t+12} (2026-10) lies outside the window."})2930# --- splice: most recent vintage covering the month31a["vintage"] = "2016-07-15"; b["vintage"] = "2026-10-28"32b_cov = b.dropna(subset=["cpi"]) # month 'available' in vintage only if it has a value33lvl = pd.concat([a[~a.date.isin(b.date)], b]).sort_values("date").reset_index(drop=True)34# for a month listed but empty in the newest vintage, the older vintage has no value either (checked):35for m in missing_src:36 assert m not in set(a.date)37lvl = lvl.set_index(pd.to_datetime(lvl.date))3839overlap = a.merge(b, on="date", suffixes=("_16", "_26"))40diff = (overlap.cpi_16 - overlap.cpi_26).abs()41issues.append({"issue": f"The two vintages overlap for {len(overlap)} months ({overlap.date.min()[:7]}..{overlap.date.max()[:7]}) and disagree (revised seasonal factors): max abs level difference {diff.max():.3f} index points.",42 "handling": "Per the availability rule the 2026-10-28 vintage is used for every month it covers (2014-01 onward); the 2016-07-15 vintage only for months up to 2013-12."})43issues.append({"issue": "For 2014-01..2014-12 the numerator P_t comes from the 2026 vintage while the lag P_{t-12} (2013) comes from the 2016 vintage, so these 12 growth rates mix two vintages (seasonal-factor vintage mismatch).",44 "handling": "Computed as specified by the rule; magnitude quantified as a sensitivity (2016-vintage-only growth for 2014 months) in the report and in the 'sensitivity' block of results.json."})4546# --- audit of the hand-assembled file (not used for estimation)47m = asm.merge(lvl.reset_index(drop=True)[["date", "cpi", "vintage"]], on="date", how="outer", suffixes=("_asm", "_ref"), indicator=True)48only_asm = m[m._merge == "left_only"].date.tolist()49only_ref = m[m._merge == "right_only"].date.tolist()50both = m[m._merge == "both"]51mism = both[~((both.cpi_asm == both.cpi_ref) | (both.cpi_asm.isna() & both.cpi_ref.isna()))]52wrongv = both[both.vintage_asm != both.vintage_ref]53if len(wrongv):54 issues.append({"issue": f"cpi_assembled.csv uses the 2016-07-15 vintage for {len(wrongv)} months ({wrongv.date.min()[:7]}..{wrongv.date.max()[:7]}) that are covered by the more recent 2026-10-28 vintage.",55 "handling": "Assembled file not used for estimation; levels taken directly from the two vintage files."})56if only_ref or only_asm:57 shifted = both[(both.date > "2016-06-01")]58 issues.append({"issue": f"cpi_assembled.csv is missing month(s) {[d[:7] for d in only_ref]} and contains month(s) {[d[:7] for d in only_asm]} absent from the sources; from 2016-07 onward the 2026-vintage values are dated one month late (e.g. 2026-09 value 336.298 appears as 2026-10). {len(mism)} common months have values differing from the authoritative splice.",59 "handling": "Assembled file treated as unreliable and used only for this audit; all estimates use the vintage files."})6061# --- growth62full = pd.date_range("1990-01-01", W1, freq="MS")63P = lvl.cpi.reindex(full)64V = lvl.vintage.reindex(full)65lag = P.shift(12)66g = 100 * (P / lag - 1)67win = pd.date_range(W0, W1, freq="MS")68df = pd.DataFrame({"date": win.strftime("%Y-%m-%d"), "cpi": P.reindex(win).values,69 "cpi_lag12": lag.reindex(win).values, "vintage_used": V.reindex(win).values,70 "growth_12m": g.reindex(win).values})71n0 = len(df)72s1 = df[df.cpi.notna()]73s2 = s1[s1.cpi_lag12.notna()]74an = s2.reset_index(drop=True)75ex = pd.DataFrame([76 {"step": 1, "rule": "P_t missing in all sources", "n_before": n0, "n_after": len(s1), "n_removed": n0 - len(s1),77 "rationale": "No index level for month t in the most recent vintage covering it (" + ", ".join(d[:7] for d in df[df.cpi.isna()].date) + "); imputation not allowed."},78 {"step": 2, "rule": "P_{t-12} missing in all sources", "n_before": len(s1), "n_after": len(s2), "n_removed": len(s1) - len(s2),79 "rationale": "Lagged index level required for g_t; imputation not allowed."},80])81an_out = an.copy()82an_out["vintage_used"] = an_out.vintage_used83an_out.to_csv(os.path.join(O, "analytic.csv"), index=False, float_format="%.10g")84ex.to_csv(os.path.join(O, "exclusions.csv"), index=False)8586gs = an.set_index("date").growth_12m87def gm(ym):88 k = ym + "-01"89 return float(round(gs[k], 6)) if k in gs.index else None90mean = lambda s: float(round(s.mean(), 6))91missing_months = df[df.growth_12m.isna()].date.str[:7].tolist()9293# sensitivity: 2014 growth from 2016 vintage alone (same-vintage numerator and lag)94Pa = a.set_index(pd.to_datetime(a.date)).cpi95ga = 100 * (Pa / Pa.shift(12) - 1)96m14 = pd.date_range("2014-01-01", "2016-06-01", freq="MS")97sens = pd.DataFrame({"spliced": g.reindex(m14).values, "v2016_only": ga.reindex(m14).values}, index=m14)98sens["diff"] = sens.spliced - sens.v2016_only99alt = gs.copy(); alt.index = pd.to_datetime(alt.index)100alt.loc[m14[:12]] = ga.reindex(m14[:12]).values101102res = {103 "task_id": "F1-conf-01",104 "first_growth_month": an.date.iloc[0][:7],105 "last_month": lvl.cpi.dropna().index.max().strftime("%Y-%m"),106 "n_months_growth": int(len(an)),107 "n_months_missing_growth": int(n0 - len(an)),108 "n_analysed": int(len(an)),109 "mean_growth_12m": mean(gs),110 "mean_growth_12m_2001_2013": mean(gs[gs.index <= "2013-12-01"]),111 "mean_growth_12m_2014_2026": mean(gs[gs.index >= "2014-01-01"]),112 "growth_2014_06": gm("2014-06"), "growth_2016_07": gm("2016-07"), "growth_2020_05": gm("2020-05"),113 "growth_2022_06": gm("2022-06"), "growth_2026_09": gm("2026-09"),114 "max_growth_12m": float(round(gs.max(), 6)), "max_growth_month": gs.idxmax()[:7],115 "min_growth_12m": float(round(gs.min(), 6)), "min_growth_month": gs.idxmin()[:7],116 "feasibility": "partial" if n0 - len(an) else "complete",117 "feasibility_reason": (f"g_t not computable for {missing_months}: the source (vintage 2026-10-28) publishes no CPI value for 2025-10; "118 "imputation is prohibited, so all window statistics use the remaining months.") if n0 - len(an) else None,119 "data_issues": issues,120 "n_input_rows": n_rows,121 "random_seed": None,122 "missing_growth_months": missing_months,123 "sensitivity": {124 "growth_2014_mixed_minus_2016vintage_mean_abs_pp": float(round(sens["diff"].iloc[:12].abs().mean(), 6)),125 "growth_2014_mixed_minus_2016vintage_max_abs_pp": float(round(sens["diff"].iloc[:12].abs().max(), 6)),126 "growth_2015_2016H1_vintage_revision_max_abs_pp": float(round(sens["diff"].iloc[12:].abs().max(), 6)),127 "mean_growth_12m_2014_2026_using_2016vintage_for_2014": mean(alt[alt.index >= "2014-01-01"]),128 "growth_2014_06_2016vintage_only": float(round(ga["2014-06-01"], 6)),129 },130}131with open(os.path.join(O, "results.json"), "w") as f:132 json.dump(res, f, indent=2)133print(json.dumps(res, indent=2))134