Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Assistant: claude-code Assistant-Model: opus Assistant-Process: 14913@bnt-lap001 Assistant-Session: add7ef69-b901-4756-9a67-34cf1125d96f
471 lines
23 KiB
Python
471 lines
23 KiB
Python
"""Model matching and regime flipping (ASSEMBLY-WP-0012-T01/T02/T03/T06).
|
||
|
||
Candidates: band anchor (B1-B4) x legitimacy loop (off/on), each from the
|
||
jurisdiction's calibrated state (calibrate_data.start_state). Method:
|
||
research/ModelMatching.md.
|
||
|
||
Fit per candidate: the adoption indicator is taken as proportional to K_Omega,
|
||
so two free values are fitted, a time offset t0 and a scale s
|
||
(log y = log s + log K_Omega(t0 + tau) + noise, noise sd = the series'
|
||
log_sigma). The calibrated K_Omega level range at the aligned time
|
||
(calibration_targets.json) is scored as well, so shape and level both count.
|
||
`sentiment_ll` scores a sentiment slope, but it is not used yet: the Pew
|
||
points (2021, 2026) precede the model start (checkpoint 0). It is switched on
|
||
once the sentiment series has points after checkpoint 0. All candidates have the same
|
||
number of fitted values, so AIC weights reduce to likelihood weights.
|
||
|
||
Flip: the best candidate changes and the new best has weight >= FLIP_WEIGHT
|
||
in two consecutive evaluations (research/ModelMatching.md § 3).
|
||
|
||
Usage: python3 -B simulations/gadm02/matching.py # match observed series
|
||
python3 -B simulations/gadm02/matching.py --synthetic # recovery and flip test (single-state filter)
|
||
python3 -B simulations/gadm02/matching.py --imm # IMM filter, several indicators (WP-0017)
|
||
"""
|
||
|
||
from __future__ import annotations
|
||
|
||
import json
|
||
import math
|
||
import random
|
||
import sys
|
||
from pathlib import Path
|
||
|
||
sys.dont_write_bytecode = True
|
||
|
||
import calibrate_data as cd # noqa: E402
|
||
import gadm04 as g # noqa: E402
|
||
import grounding as gr # noqa: E402
|
||
|
||
HERE = Path(__file__).resolve().parent
|
||
OUT = HERE / "results"
|
||
SERIES = json.loads((HERE / "observed_series.json").read_text())["series"]
|
||
CANDIDATES = [(a, lp) for a in ("B1", "B2", "B3", "B4") for lp in ("off", "on")]
|
||
STEP = 0.1
|
||
T_PATH = 30.0
|
||
FLIP_WEIGHT = 0.75
|
||
SWITCH = 0.02 # prior probability per quarter that the regime changes (filter)
|
||
_CACHE: dict = {}
|
||
|
||
|
||
def path(jur: str, anchor: str, loop: str) -> list[dict]:
|
||
key = (jur, anchor, loop)
|
||
if key not in _CACHE:
|
||
p, s = cd.start_state(jur, anchor)
|
||
_CACHE[key] = g.simulate({**p, "sent": 1.0 if loop == "on" else 0.0}, s, T=T_PATH, record_every=1)
|
||
return _CACHE[key]
|
||
|
||
|
||
def interp(rows: list[dict], t: float, key: str) -> float:
|
||
x = t / STEP
|
||
i = min(int(x), len(rows) - 2)
|
||
w = x - i
|
||
return (1 - w) * rows[i][key] + w * rows[i + 1][key]
|
||
|
||
|
||
def fit_series(rows: list[dict], points: list, log_sigma: float,
|
||
fixed: tuple[float, float] | None = None) -> tuple[float, float, float]:
|
||
"""(log-likelihood, t0, log s) for an indicator proportional to K_Omega.
|
||
|
||
fixed = (t0, log s): no fitting, score the points at this alignment and scale. Used after the
|
||
first evaluation: the calendar fixes how t0 advances and the measurement fixes the scale."""
|
||
if fixed is not None:
|
||
t0, ls = fixed
|
||
if t0 + points[-1][0] > T_PATH - STEP:
|
||
return -1e9, t0, ls
|
||
sse = sum((math.log(y) - math.log(max(1e-9, interp(rows, t0 + tau, "Ko"))) - ls) ** 2 for tau, y in points)
|
||
return -0.5 * sse / log_sigma ** 2, t0, ls
|
||
span = points[-1][0]
|
||
best = None
|
||
t0 = 0.0
|
||
while t0 + span <= T_PATH - 1:
|
||
r = [math.log(y) - math.log(max(1e-9, interp(rows, t0 + tau, "Ko"))) for tau, y in points]
|
||
ls = sum(r) / len(r)
|
||
sse = sum((x - ls) ** 2 for x in r)
|
||
if best is None or sse < best[0]:
|
||
best = (sse, t0, ls)
|
||
t0 += 0.05
|
||
sse, t0, ls = best
|
||
return -0.5 * sse / log_sigma ** 2, t0, ls
|
||
|
||
|
||
def sentiment_ll(rows: list[dict], t0: float, series: dict) -> float:
|
||
(ta, ya), (tb, yb) = series["points"][0], series["points"][-1]
|
||
obs = (yb - ya) / (tb - ta)
|
||
span = tb - ta
|
||
a = max(0.0, t0 - span) # the sentiment window ends at the adoption alignment
|
||
mod = (interp(rows, a + span, "Sig") - interp(rows, a, "Sig")) / span
|
||
return -0.5 * ((mod - obs) * span / series["sigma"]) ** 2
|
||
|
||
|
||
def score(jur: str, rows_for: callable, adoption: dict, sentiment: dict | None = None,
|
||
level: dict | None = None, fixed: dict | None = None) -> dict:
|
||
out = {}
|
||
for a, lp in CANDIDATES:
|
||
rows = rows_for(jur, a, lp)
|
||
ll, t0, ls = fit_series(rows, adoption["points"], adoption["log_sigma"],
|
||
fixed=fixed.get((a, lp)) if fixed else None)
|
||
if level is not None: # K_Omega level at the end of the observed window
|
||
ko = interp(rows, t0 + adoption["points"][-1][0], "Ko")
|
||
ll += -0.5 * cd.z(ko, level) ** 2
|
||
if sentiment is not None:
|
||
ll += sentiment_ll(rows, t0, sentiment)
|
||
out[(a, lp)] = {"ll": ll, "t0": t0, "s": math.exp(ls), "ls": ls}
|
||
mx = max(v["ll"] for v in out.values())
|
||
z = sum(math.exp(v["ll"] - mx) for v in out.values())
|
||
for v in out.values():
|
||
v["w"] = math.exp(v["ll"] - mx) / z
|
||
return out
|
||
|
||
|
||
def best(res: dict) -> tuple:
|
||
return max(res, key=lambda k: res[k]["w"])
|
||
|
||
|
||
def marginal(res: dict, part: int) -> dict:
|
||
"""Weights summed over the other part of the candidate (0 = anchor, 1 = loop state)."""
|
||
out: dict = {}
|
||
for c, v in res.items():
|
||
out[(c[part],)] = out.get((c[part],), 0.0) + v["w"]
|
||
return {k: {"w": w} for k, w in out.items()}
|
||
|
||
|
||
def flips(evaluations: list[tuple[float, dict]], by: int | None = None) -> list[tuple[float, tuple, tuple]]:
|
||
"""Flip events: best changes and the new best holds weight >= FLIP_WEIGHT twice in a row.
|
||
|
||
by = 0 judges the anchor (weights summed over the loop state), by = 1 the loop state; None the full
|
||
candidate. Judging the anchor matters when the loop is not identifiable: its weight would otherwise
|
||
split between two indistinguishable candidates and never reach FLIP_WEIGHT."""
|
||
if by is not None:
|
||
evaluations = [(t, marginal(res, by)) for t, res in evaluations]
|
||
events, current, pending = [], None, None
|
||
for t, res in evaluations:
|
||
b, wb = best(res), res[best(res)]["w"]
|
||
if current is None:
|
||
if wb >= FLIP_WEIGHT:
|
||
current = b
|
||
continue
|
||
if b != current and wb >= FLIP_WEIGHT:
|
||
if pending == b:
|
||
events.append((t, current, b))
|
||
current, pending = b, None
|
||
else:
|
||
pending = b
|
||
else:
|
||
pending = None
|
||
return events
|
||
|
||
|
||
def fmt(c: tuple) -> str:
|
||
return f"{c[0]}/loop {c[1]}"
|
||
|
||
|
||
def main_observed() -> None:
|
||
lines = [f"# {g.MODEL_VERSION}: model matching on observed series", "",
|
||
"Generated by `simulations/gadm02/matching.py` (ASSEMBLY-WP-0012). Method:",
|
||
"`research/ModelMatching.md`. Data: `observed_series.json`. Weights are likelihood (AIC) weights",
|
||
"over the eight candidates; all candidates fit the same two values (t0, s). Each score adds the",
|
||
"calibrated K_Ω level range at the end of the window. Sentiment is not scored yet (data precede",
|
||
"checkpoint 0).", ""]
|
||
for jur, adop, sent in (("US", "US_adoption", None), ("EU", "EU_adoption", None),
|
||
("DE", "DE_adoption", None)):
|
||
res = score(jur, path, SERIES[adop], SERIES[sent] if sent else None, level=cd.targets(jur)["Ko"])
|
||
b = best(res)
|
||
lines += [f"## {jur}", "", f"Series: {SERIES[adop]['source']}" + (f"; {SERIES[sent]['source']}" if sent else ""),
|
||
"", "| Candidate | log-likelihood | weight | t0 | K_Ω at window end |", "|---|---|---|---|---|"]
|
||
for c in sorted(res, key=lambda k: -res[k]["w"]):
|
||
v = res[c]
|
||
lines.append(f"| {fmt(c)} | {v['ll']:.1f} | {v['w']:.3f} | {v['t0']:.2f} | "
|
||
f"{interp(path(jur, *c), v['t0'] + SERIES[adop]['points'][-1][0], 'Ko'):.3f} |")
|
||
anchors = {}
|
||
for (a, _lp), v in res.items():
|
||
anchors[a] = anchors.get(a, 0) + v["w"]
|
||
loops = {lp: sum(v["w"] for (a, l2), v in res.items() if l2 == lp) for lp in ("off", "on")}
|
||
lines += ["", f"**Best:** {fmt(b)} (weight {res[b]['w']:.2f}). Anchor weights: "
|
||
+ ", ".join(f"{a} {w:.2f}" for a, w in anchors.items())
|
||
+ f". Loop weights: off {loops['off']:.2f}, on {loops['on']:.2f}.", ""]
|
||
(OUT / "matching.md").write_text("\n".join(lines))
|
||
print("\n".join(lines))
|
||
|
||
|
||
def synthetic_rows(jur: str, first: tuple, second: tuple | None, t_switch: float) -> list[dict]:
|
||
"""Truth: candidate `first` until t_switch, then `second` parameters continuing from that state."""
|
||
p1, s1 = cd.start_state(jur, first[0])
|
||
r1 = g.simulate({**p1, "sent": 1.0 if first[1] == "on" else 0.0}, s1, T=t_switch, record_every=1)
|
||
if second is None:
|
||
return g.simulate({**p1, "sent": 1.0 if first[1] == "on" else 0.0}, s1, T=T_PATH, record_every=1)
|
||
p2, _s2 = cd.start_state(jur, second[0])
|
||
state = {k: r1[-1][k] for k in g.STOCKS}
|
||
r2 = g.simulate({**p2, "sent": 1.0 if second[1] == "on" else 0.0}, state, T=T_PATH - t_switch, record_every=1)
|
||
for r in r2[1:]:
|
||
r["t"] = round(r["t"] + t_switch, 6)
|
||
return r1 + r2[1:]
|
||
|
||
|
||
_PARAMS: dict = {}
|
||
|
||
|
||
def regime_params(jur: str, cand: tuple) -> dict:
|
||
"""Parameters of a regime (anchor dynamics, loop switch) for jurisdiction `jur`."""
|
||
if (jur, cand) not in _PARAMS:
|
||
p, _s = cd.start_state(jur, cand[0])
|
||
_PARAMS[(jur, cand)] = {**p, "sent": 1.0 if cand[1] == "on" else 0.0}
|
||
return _PARAMS[(jur, cand)]
|
||
|
||
|
||
def filter_step(jur: str, state: dict, pts: list, scale_ls: float, log_sigma: float,
|
||
level: dict | None, window: float, sent_pts: list | None = None, sent_sd: float = 0.02) -> dict:
|
||
"""Score every regime by running its dynamics forward from `state` over the window.
|
||
|
||
sent_pts: optional sentiment observations (tau, Sigma) within the window."""
|
||
out = {}
|
||
for c in CANDIDATES:
|
||
rows = g.simulate(regime_params(jur, c), state, T=window, record_every=1)
|
||
sse = sum((math.log(y) - math.log(max(1e-9, interp(rows, tau, "Ko"))) - scale_ls) ** 2 for tau, y in pts)
|
||
ll = -0.5 * sse / log_sigma ** 2
|
||
for tau, y in sent_pts or []:
|
||
ll += -0.5 * ((interp(rows, tau, "Sig") - y) / sent_sd) ** 2
|
||
if level is not None:
|
||
ll += -0.5 * cd.z(rows[-1]["Ko"], level) ** 2
|
||
out[c] = {"ll": ll, "end": {k: rows[-1][k] for k in g.STOCKS}}
|
||
mx = max(v["ll"] for v in out.values())
|
||
zs = sum(math.exp(v["ll"] - mx) for v in out.values())
|
||
for v in out.values():
|
||
v["w"] = math.exp(v["ll"] - mx) / zs
|
||
return out
|
||
|
||
|
||
def synthetic_run(truth: list[dict], start: float, years: float, window: float, seed: int,
|
||
log_sigma: float = 0.03, sentiment: bool = False) -> list[tuple[float, dict]]:
|
||
"""Quarterly regime filter on noisy biweekly observations (research/ModelMatching.md § 3).
|
||
|
||
First window: whole-path alignment fit (t0, scale) gives the initial state and the scale. Then,
|
||
every quarter, each regime runs forward from the current state estimate over the last window;
|
||
the state estimate advances with the best regime."""
|
||
rnd = random.Random(seed)
|
||
obs, t = [], start
|
||
while t <= start + years:
|
||
obs.append((t, interp(truth, t, "Ko") * 25.0 * math.exp(rnd.gauss(0, log_sigma))))
|
||
t += 2 / 52
|
||
sobs = [(start + k, interp(truth, start + k, "Sig") + rnd.gauss(0, 0.02)) for k in range(int(years) + 1)] \
|
||
if sentiment else []
|
||
te0 = start + window
|
||
pts0 = [(tt - start, y) for tt, y in obs if start <= tt <= te0]
|
||
first = score("US", path, {"points": pts0, "log_sigma": log_sigma},
|
||
level={"lo": interp(truth, te0, "Ko") - 0.075, "hi": interp(truth, te0, "Ko") + 0.075})
|
||
b = best(first)
|
||
scale_ls = first[b]["ls"]
|
||
t_est = first[b]["t0"] + window
|
||
rows_b = path("US", *b)
|
||
i_est = int(round(t_est / STEP))
|
||
est = {k: rows_b[i_est][k] for k in g.STOCKS}
|
||
evals = [(te0, first)]
|
||
step = 0.25
|
||
te = te0 + step
|
||
while te <= start + years + 1e-9:
|
||
pts = [(tt - (te - step), y) for tt, y in obs if te - step - 1e-9 <= tt <= te]
|
||
ko = interp(truth, te, "Ko")
|
||
spts = [(tt - (te - step), y) for tt, y in sobs if te - step < tt <= te + 1e-9]
|
||
res = filter_step("US", est, pts, scale_ls, log_sigma, {"lo": ko - 0.075, "hi": ko + 0.075}, step, spts)
|
||
prev = evals[-1][1]
|
||
n = len(CANDIDATES)
|
||
prior = {c: (1 - SWITCH) * prev[c]["w"] + SWITCH * (1 - prev[c]["w"]) / (n - 1) for c in CANDIDATES}
|
||
mx = max(res[c]["ll"] for c in CANDIDATES)
|
||
post = {c: prior[c] * math.exp(res[c]["ll"] - mx) for c in CANDIDATES}
|
||
zs = sum(post.values())
|
||
for c in CANDIDATES:
|
||
res[c]["w"] = post[c] / zs # evidence accumulates across quarters
|
||
est = res[best(res)]["end"]
|
||
evals.append((te, res))
|
||
te += step
|
||
return evals
|
||
|
||
|
||
def main_synthetic() -> None:
|
||
lines = [f"# {g.MODEL_VERSION}: matching and flip detection on synthetic data", "",
|
||
"Generated by `simulations/gadm02/matching.py --synthetic` (ASSEMBLY-WP-0012-T06).",
|
||
"US calibrated state. Truth simulated from a known candidate (optionally switching regime at",
|
||
"t_switch), observed biweekly as 25 × K_Ω with log-noise sd 0.03, from t 2 to t 8. Quarterly",
|
||
"regime filter: the first year fixes the initial state and the scale (whole-path alignment); then",
|
||
"every quarter each regime runs forward from the current state estimate over the last quarter, and",
|
||
"the estimate advances with the best regime. Weights accumulate across quarters (forward filter,",
|
||
f"switch probability {SWITCH} per quarter). Level evidence: true K_Ω ± 0.075 (as wide as the real",
|
||
f"US target). Flip: new best with weight ≥ {FLIP_WEIGHT} twice in a row.", "",
|
||
"| Indicators | Truth | Seed | Best at first evaluation | Flips detected (t, from → to) | Verdict |",
|
||
"|---|---|---|---|---|---|"]
|
||
cases = [(("B3", "off"), None, None), (("B2", "off"), None, None),
|
||
(("B2", "off"), ("B4", "off"), 4.0), (("B3", "off"), ("B3", "on"), 4.0)]
|
||
for first, second, ts in cases:
|
||
truth = synthetic_rows("US", first, second, ts or T_PATH)
|
||
for seed, sentiment in ((1, False), (2, False), (3, False), (1, True), (2, True), (3, True)):
|
||
ev = synthetic_run(truth, 2.0, 6.0, 1.0, seed, sentiment=sentiment)
|
||
fl = flips(ev)
|
||
b0 = best(ev[0][1])
|
||
if second is None:
|
||
verdict = "ok (no false flip)" if not fl else "FALSE FLIP"
|
||
else:
|
||
hit = [e for e in fl if e[0] >= ts and e[2][0] == second[0] and (second[1] == first[1] or e[2][1] == second[1])]
|
||
verdict = f"detected, delay {hit[0][0] - ts:.2f} y" if hit else "missed"
|
||
label = fmt(first) + (f" → {fmt(second)} at t {ts}" if second else "")
|
||
lines.append(f"| {'adoption + annual sentiment' if sentiment else 'adoption'} | {label} | {seed} | {fmt(b0)} | "
|
||
+ ("; ".join(f"{t:.2f}: {fmt(a)} → {fmt(b)}" for t, a, b in fl) or "none")
|
||
+ f" | {verdict} |")
|
||
lines.append("")
|
||
(OUT / "matching_synthetic.md").write_text("\n".join(lines))
|
||
print("\n".join(lines))
|
||
|
||
|
||
|
||
|
||
# --- WP-0017-T02/T03: interacting-multiple-model filter with several indicators -----------------
|
||
|
||
def observe(rows: list[dict], t: float, kind: str, dO: float = None) -> float:
|
||
"""Model value of an indicator at time t (rows recorded every STEP)."""
|
||
if kind == "adopt":
|
||
return interp(rows, t, "Ko")
|
||
if kind == "invest": # gross organism investment: dK_Omega/dt + dO K_Omega
|
||
dO = g.PARAMS["dO"] if dO is None else dO
|
||
a, b = interp(rows, max(0.0, t - STEP), "Ko"), interp(rows, t, "Ko")
|
||
return max(1e-9, (b - a) / STEP + dO * b)
|
||
if kind == "pd":
|
||
return interp(rows, t, "Pd")
|
||
if kind == "sig":
|
||
return interp(rows, t, "Sig")
|
||
raise KeyError(kind)
|
||
|
||
|
||
NOISE = {"adopt": ("log", 0.03), "invest": ("log", 0.05), "pd": ("abs", 0.01), "sig": ("abs", 0.02)}
|
||
CADENCE = {"adopt": 2 / 52, "invest": 0.25, "pd": 1 / 12, "sig": 0.5}
|
||
|
||
|
||
def obs_ll(rows: list[dict], obs: list, scales: dict) -> float:
|
||
ll = 0.0
|
||
for tau, kind, y in obs:
|
||
mode, sd = NOISE[kind]
|
||
m_ = observe(rows, tau, kind)
|
||
if mode == "log":
|
||
ll += -0.5 * ((math.log(y) - math.log(max(1e-9, m_)) - scales[kind]) / sd) ** 2
|
||
else:
|
||
ll += -0.5 * ((y - m_) / sd) ** 2
|
||
return ll
|
||
|
||
|
||
IMM_SWITCH = 0.10 # switching prior per quarter for the IMM filter (tested in main_imm)
|
||
|
||
|
||
def imm_run(truth: list[dict], start: float, years: float, seed: int, kinds: tuple,
|
||
switch: float = IMM_SWITCH, step: float = 0.25) -> list[tuple[float, dict]]:
|
||
"""IMM regime filter (research/ModelMatching.md § 4.1).
|
||
|
||
Initial states and log-scales from a whole-path alignment of the first year of adoption data;
|
||
then every quarter: mix regime states by switching probabilities, run each regime from its
|
||
mixed state, update weights by the likelihood of all indicators in the quarter."""
|
||
rnd = random.Random(seed)
|
||
obs = []
|
||
for kind in kinds:
|
||
t = start
|
||
while t <= start + years + 1e-9:
|
||
mode, sd = NOISE[kind]
|
||
v = observe(truth, t, kind)
|
||
y = v * 25.0 * math.exp(rnd.gauss(0, sd)) if kind == "adopt" else \
|
||
v * 3.0 * math.exp(rnd.gauss(0, sd)) if kind == "invest" else v + rnd.gauss(0, sd)
|
||
obs.append((t, kind, y))
|
||
t += CADENCE[kind]
|
||
te0 = start + 1.0
|
||
first_pts = [(t - start, y) for t, k, y in obs if k == "adopt" and t <= te0]
|
||
ko = interp(truth, te0, "Ko")
|
||
first = score("US", path, {"points": first_pts, "log_sigma": NOISE["adopt"][1]},
|
||
level={"lo": ko - 0.075, "hi": ko + 0.075})
|
||
b = best(first)
|
||
t_b = first[b]["t0"]
|
||
scales = {"adopt": first[b]["ls"]}
|
||
inv = [(t - start, y) for t, k, y in obs if k == "invest" and t <= te0]
|
||
if inv:
|
||
rb = path("US", *b)
|
||
scales["invest"] = sum(math.log(y) - math.log(observe(rb, t_b + tau, "invest")) for tau, y in inv) / len(inv)
|
||
states = {}
|
||
for c in CANDIDATES:
|
||
rc = path("US", *c)
|
||
i = min(len(rc) - 1, int(round((first[c]["t0"] + 1.0) / STEP)))
|
||
states[c] = {k: rc[i][k] for k in g.STOCKS}
|
||
w = {c: first[c]["w"] for c in CANDIDATES}
|
||
evals = [(te0, {c: {"w": w[c]} for c in CANDIDATES})]
|
||
n = len(CANDIDATES)
|
||
pi = {(i, c): (1 - switch) if i == c else switch / (n - 1) for i in CANDIDATES for c in CANDIDATES}
|
||
te = te0 + step
|
||
while te <= start + years + 1e-9:
|
||
window = [(t - (te - step), k, y) for t, k, y in obs if te - step < t <= te + 1e-9]
|
||
cbar = {c: sum(pi[(i, c)] * w[i] for i in CANDIDATES) for c in CANDIDATES}
|
||
lls, ends = {}, {}
|
||
for c in CANDIDATES:
|
||
mu = {i: pi[(i, c)] * w[i] / cbar[c] for i in CANDIDATES}
|
||
x0 = {k: sum(mu[i] * states[i][k] for i in CANDIDATES) for k in g.STOCKS}
|
||
rows = g.simulate(regime_params("US", c), x0, T=step, record_every=1)
|
||
lls[c] = obs_ll(rows, window, scales)
|
||
ends[c] = {k: rows[-1][k] for k in g.STOCKS}
|
||
mx = max(lls.values())
|
||
post = {c: cbar[c] * math.exp(lls[c] - mx) for c in CANDIDATES}
|
||
zs = sum(post.values())
|
||
w = {c: post[c] / zs for c in CANDIDATES}
|
||
states = ends
|
||
evals.append((te, {c: {"w": w[c]} for c in CANDIDATES}))
|
||
te += step
|
||
return evals
|
||
|
||
|
||
def main_imm() -> None:
|
||
kinds = ("adopt", "invest", "pd", "sig")
|
||
lines = [f"# {g.MODEL_VERSION}: IMM regime filter on synthetic data", "",
|
||
"Generated by `simulations/gadm02/matching.py --imm` (ASSEMBLY-WP-0017-T02/T03/T04).",
|
||
"US calibrated state, sentiment anchored (BR-009 fix). Truth from a known candidate, optionally",
|
||
"switching at t 4; observed from t 2 to t 8. Indicators and noise: adoption biweekly (log sd 0.03),",
|
||
"investment quarterly (log sd 0.05), discovery channel Pd monthly (sd 0.01), sentiment twice a",
|
||
f"year (sd 0.02). IMM filter, switch probability {IMM_SWITCH} per quarter. Flip: new best with weight",
|
||
f"≥ {FLIP_WEIGHT} twice in a row, judged on anchor weights (summed over the loop state) and, for the",
|
||
"loop case, on loop weights (summed over anchors). False flips count both.", "",
|
||
"**Acceptance (anchor regimes):** anchor correct at the first evaluation in ≥ 5 of 6 no-switch runs,",
|
||
"no false flips in 6 no-switch runs, and ≥ 4 of 6 anchor switches detected within 1.5 years. The",
|
||
"loop switch is reported separately: with sentiment anchored, switching the loop on is invisible",
|
||
"until perceived benefit moves.", "",
|
||
"| Truth | Seed | Best at start | Flips (t, from → to) | Verdict |", "|---|---|---|---|---|"]
|
||
correct = false = det = 0
|
||
loop_det = 0
|
||
cases = [(("B3", "off"), None, None, "none"), (("B2", "on"), None, None, "none"),
|
||
(("B2", "off"), ("B4", "off"), 4.0, "anchor"), (("B3", "off"), ("B2", "off"), 4.0, "anchor"),
|
||
(("B3", "off"), ("B3", "on"), 4.0, "loop")]
|
||
for first_c, second, ts, kind in cases:
|
||
truth = synthetic_rows("US", first_c, second, ts or T_PATH)
|
||
for seed in (1, 2, 3):
|
||
ev = imm_run(truth, 2.0, 6.0, seed, kinds)
|
||
fl = flips(ev, by=0) if kind != "loop" else flips(ev, by=1)
|
||
b0 = best(ev[0][1])
|
||
if kind == "none":
|
||
correct += b0[0] == first_c[0]
|
||
fl = flips(ev, by=0) + flips(ev, by=1)
|
||
false += bool(fl)
|
||
verdict = ("anchor correct" if b0[0] == first_c[0] else "anchor wrong") + ("; FALSE FLIP" if fl else "")
|
||
else:
|
||
if kind == "anchor":
|
||
hit = [e for e in fl if ts <= e[0] <= ts + 1.5 and e[2][0] == second[0]]
|
||
det += bool(hit)
|
||
else:
|
||
hit = [e for e in fl if ts <= e[0] <= ts + 1.5 and e[2][0] == second[1]]
|
||
loop_det += bool(hit)
|
||
verdict = f"detected, delay {hit[0][0] - ts:.2f} y" if hit else "missed"
|
||
label = fmt(first_c) + (f" → {fmt(second)} at t {ts}" if second else "")
|
||
lines.append(f"| {label} | {seed} | {fmt(b0)} | "
|
||
+ ("; ".join(f"{t:.2f}: {a[0]} → {b[0]}" for t, a, b in fl) or "none")
|
||
+ f" | {verdict} |")
|
||
ok = correct >= 5 and false == 0 and det >= 4
|
||
lines += ["", f"Anchor correct at start {correct}/6, false flips {false}/6, anchor switches detected {det}/6, "
|
||
f"loop switches detected {loop_det}/3. Acceptance: {'**passed**' if ok else '**not passed**'}.", ""]
|
||
(OUT / "matching_imm.md").write_text("\n".join(lines))
|
||
print("\n".join(lines))
|
||
|
||
|
||
if __name__ == "__main__":
|
||
if "--imm" in sys.argv:
|
||
main_imm()
|
||
elif "--synthetic" in sys.argv:
|
||
main_synthetic()
|
||
else:
|
||
main_observed()
|