#!/usr/bin/env python3 """Ordered Phase 4P-B dynamic revalidation on the canonical runtime. R1/R2 consume the already sealed Phase 4P results. R3-R6 are run in order and the script stops immediately if a level cannot reproduce its declared historical result. No frozen artifact is written. """ from __future__ import annotations import argparse import hashlib import json import math import time from pathlib import Path from typing import Any import numpy as np import pandas as pd import pyarrow.parquet as pq import torch import yaml import phase4p from phase4l import pn_drive from phase4m import load_odor_inputs ROOT = Path("/data/digital-drosophila") def save(path: Path, value: object) -> None: path.parent.mkdir(parents=True, exist_ok=True) path.write_text(json.dumps(value, indent=2, sort_keys=True, default=str) + "\n") def ahash(*values: np.ndarray) -> str: h = hashlib.sha256() for value in values: h.update(np.ascontiguousarray(value).tobytes()) return h.hexdigest() def metrics(a: np.ndarray, b: np.ndarray) -> dict[str, float]: a = a.astype(np.float64); b = b.astype(np.float64) aa = a > 0; bb = b > 0; union = np.count_nonzero(aa | bb) den = np.linalg.norm(a) * np.linalg.norm(b) return {"euclidean": float(np.linalg.norm(a - b)), "jaccard": float(np.count_nonzero(aa & bb) / union) if union else 1.0, "cosine": float(np.dot(a, b) / den) if den else (1.0 if not np.any(a) and not np.any(b) else 0.0)} class DeterministicAPLSim: """Frozen APL overlay on top of the Phase 4P canonical LIF kernel.""" def __init__(self, rt: phase4p.CanonicalRuntime): self.rt = rt; self.n = rt.n; self.device = rt.device; self.dt = rt.dt meta = rt.table.to_pandas().sort_values("neuron_idx").reset_index(drop=True) self.meta = meta cls = meta.annotation_class.fillna("").astype(str).to_numpy() typ = meta.annotation_type.fillna("").astype(str).to_numpy() self.kc = np.flatnonzero(cls == "Kenyon_Cell").astype(np.int32) self.mbon = np.flatnonzero(cls == "MBON").astype(np.int32) self.pn = np.flatnonzero(cls == "ALPN").astype(np.int32) self.apl = np.flatnonzero(typ == "APL").astype(np.int32) if len(self.apl) != 2 or len(self.kc) != 4064: raise RuntimeError("frozen APL/KC identity changed") self.sources = np.repeat(np.arange(self.n, dtype=np.int32), np.diff(rt.indptr_np).astype(np.int64)) idx = rt.indices_np; raw = rt.raw_np apl_out = np.isin(self.sources, self.apl) self.apl_out_ids = np.flatnonzero(apl_out).astype(np.int64) input_matrix = np.zeros((len(self.apl), self.n), np.float32) output_matrix = np.zeros((len(self.apl), len(self.kc)), np.float32) kc_lookup = np.full(self.n, -1, np.int32); kc_lookup[self.kc] = np.arange(len(self.kc), dtype=np.int32) allowed = np.isin(self.sources, np.concatenate((self.kc, self.pn))) & np.isin(idx, self.apl) apl_lookup = {int(x): i for i, x in enumerate(self.apl)} for eid in np.flatnonzero(allowed): input_matrix[apl_lookup[int(idx[eid])], int(self.sources[eid])] += float(raw[eid]) for eid in np.flatnonzero(apl_out & np.isin(idx, self.kc)): output_matrix[apl_lookup[int(self.sources[eid])], kc_lookup[int(idx[eid])]] += float(raw[eid]) self.input_scale = float(input_matrix[input_matrix > 0].sum() / np.count_nonzero(np.any(input_matrix > 0, axis=0))) self.output_scale = float(output_matrix[output_matrix > 0].mean()) self.input_gpu = torch.as_tensor(input_matrix / self.input_scale, dtype=torch.float32, device=self.device) self.output_gpu = torch.as_tensor(output_matrix / self.output_scale, dtype=torch.float32, device=self.device) self.kc_gpu = torch.as_tensor(self.kc, dtype=torch.int64, device=self.device) self.apl_out_gpu = torch.as_tensor(self.apl_out_ids, dtype=torch.int64, device=self.device) self.base = rt.weights.clone() p4h = ROOT / "experiments/phase4h-dopamine-plasticity" self.identity = pd.read_parquet(ROOT / "experiments/phase4n-credit-assignment/edge-identity/frozen-edge-identity.parquet") self.edge_ids = self.identity.edge_id.to_numpy(np.int64) self.edge_gpu = torch.as_tensor(self.edge_ids, dtype=torch.int64, device=self.device) self.base_candidate = self.base[self.edge_gpu].cpu().numpy().astype(np.float32) self.confirmed = pd.read_parquet(p4h / "mappings/dan-to-mbon-confirmed.parquet") down = pd.read_parquet(p4h / "mappings/mbon-to-downstream-actual-edges.parquet") self.downstream = np.unique(down.target_neuron_idx.to_numpy(np.int32)) def dan_targets(self, dtype: str, compartment: str) -> np.ndarray: x = self.confirmed[(self.confirmed.source_annotation_type == dtype) & (self.confirmed.source_compartment_verbatim == compartment)] out = np.sort(x.source_neuron_idx.unique().astype(np.int32)) if not len(out): raise RuntimeError(f"no confirmed DAN for {dtype}/{compartment}") return out def install(self, E: np.ndarray | None) -> None: self.rt.weights.copy_(self.base) if E is not None: self.rt.weights[self.edge_gpu] = torch.as_tensor(self.base_candidate * E.astype(np.float32), device=self.device) self.rt.weights[self.apl_out_gpu] = 0.0 def run(self, targets: np.ndarray, rates_hz: np.ndarray, seed: int, duration_ms: float, E: np.ndarray | None = None, dan_targets: np.ndarray | None = None, dan_pulses: list[float] | None = None, collect_events: bool = False) -> dict[str, Any]: targets = np.asarray(targets, np.int32); rates_hz = np.asarray(rates_hz, np.float32) self.install(E); self.rt.reset(); steps = int(round(duration_ms / self.dt)) target_gpu = torch.as_tensor(targets, dtype=torch.int64, device=self.device) gen = torch.Generator(device=self.device); gen.manual_seed(int(seed)) probs = torch.as_tensor(rates_hz * self.dt / 1000.0, dtype=torch.float32, device=self.device) forced = (torch.rand((steps, len(targets)), generator=gen, device=self.device) < probs[None, :]).to(torch.uint8) dan = np.empty(0, np.int32) if dan_targets is None else np.asarray(dan_targets, np.int32) dan_gpu = torch.as_tensor(dan, dtype=torch.int64, device=self.device) dan_forced = torch.zeros((steps, len(dan)), dtype=torch.uint8, device=self.device) for pulse in dan_pulses or []: s = int(round(pulse / self.dt)) if 0 <= s < steps: dan_forced[s] = True dan_steps = {int(round(x / self.dt)) for x in (dan_pulses or [])} apl_state = torch.zeros(len(self.apl), dtype=torch.float32, device=self.device) kc_trace = torch.zeros(len(self.kc), dtype=torch.float32, device=self.device) first = torch.full((self.n,), -1, dtype=torch.int32, device=self.device) trace_at_dan = [] apl_decay = np.float32(math.exp(-self.dt / 20.0)); trace_decay = np.float32(math.exp(-self.dt / 100.0)) kc_steps: list[np.ndarray] = []; kc_nodes: list[np.ndarray] = [] torch.cuda.synchronize(); t0 = time.perf_counter() for step in range(steps): slot = step % self.rt.qlen self.rt.ext.state_step(self.rt.v, self.rt.g, self.rt.ref, self.rt.queue[slot], self.rt.spike_u8, self.rt.counts, float(self.rt.gd), float(self.rt.alpha), float(self.rt.rest), float(self.rt.resetv), float(self.rt.threshold), self.rt.ref_steps) self.rt.ext.force_step(self.rt.v, self.rt.ref, self.rt.spike_u8, self.rt.counts, target_gpu, forced[step], float(self.rt.resetv), self.rt.ref_steps) if len(dan): self.rt.ext.force_step(self.rt.v, self.rt.ref, self.rt.spike_u8, self.rt.counts, dan_gpu, dan_forced[step], float(self.rt.resetv), self.rt.ref_steps) sources = torch.nonzero(self.rt.spike_u8, as_tuple=False).flatten() if sources.numel(): missing = first[sources] < 0; first[sources[missing]] = step self.rt.canonical_fanout(sources, self.rt.queue[(step + self.rt.delay_steps) % self.rt.qlen]) kc_spike = self.rt.spike_u8[self.kc_gpu] kc_trace.mul_(float(trace_decay)).add_(kc_spike.to(torch.float32)).clamp_max_(1.0) if step in dan_steps: trace_at_dan.append(kc_trace.clone()) if collect_events: active = torch.nonzero(kc_spike, as_tuple=False).flatten() if active.numel(): node = self.kc_gpu[active].cpu().numpy().astype(np.int32) kc_nodes.append(node); kc_steps.append(np.full(len(node), step, np.int32)) # Frozen APL graded-inhibition overlay. Each KC target is unique in # this assignment, avoiding index_add atomics. apl_state.mul_(float(apl_decay)).add_(torch.mv(self.input_gpu, self.rt.spike_u8.to(torch.float32))) inhibition = torch.matmul(apl_state, self.output_gpu) * np.float32(0.01) self.rt.g[self.kc_gpu] = self.rt.g[self.kc_gpu] - inhibition self.rt.step_index += 1 torch.cuda.synchronize(); wall = time.perf_counter() - t0 counts = self.rt.counts.cpu().numpy().copy(); first_np = first.cpu().numpy() ev_step = np.concatenate(kc_steps) if kc_steps else np.empty(0, np.int32) ev_node = np.concatenate(kc_nodes) if kc_nodes else np.empty(0, np.int32) return {"counts": counts, "kc_counts": counts[self.kc], "mbon_counts": counts[self.mbon], "downstream_counts": counts[self.downstream], "first": first_np, "final_v": self.rt.v.cpu().numpy().copy(), "final_g": self.rt.g.cpu().numpy().copy(), "trace": np.zeros(len(self.kc), np.float32) if not trace_at_dan else trace_at_dan[-1].cpu().numpy().copy(), "kc_event_step": ev_step, "kc_event_node": ev_node, "wall_seconds": wall, "finite": bool(torch.isfinite(self.rt.v).all() and torch.isfinite(self.rt.g).all())} def r1_r2(exp: Path) -> tuple[dict[str, Any], dict[str, Any]]: p4c = json.loads((exp / "phase4c-v2/summary.json").read_text()) gate = json.loads((exp / "audit/candidate-gate.json").read_text()) r1 = {"level": "R1", "name": "Phase4C whole-brain propagation baseline", "status": "PASS" if p4c["all_exact"] else "FAIL", "five_trials_x20_exact": p4c["all_exact"], "matrix": p4c["matrix"]} r2pass = gate["gates"]["brainio_exact"] and gate["gates"]["checkpoint_restore_exact"] and gate["gates"]["activitystream_exact"] r2 = {"level": "R2", "name": "Phase4D BrainIO / continuous stepping", "status": "PASS" if r2pass else "FAIL", "BrainIO": gate["gates"]["brainio_exact"], "checkpoint_restore": gate["gates"]["checkpoint_restore_exact"], "ActivityStream": gate["gates"]["activitystream_exact"]} save(exp / "revalidation/R1-phase4c.json", r1); save(exp / "revalidation/R2-phase4d.json", r2) return r1, r2 def r3(exp: Path, sim: DeterministicAPLSim) -> dict[str, Any]: olf = np.flatnonzero(sim.rt.pops["olfactory"]).astype(np.int32); olf.sort(); n = 64 held = [olf[(6 + i) * n:(7 + i) * n] for i in range(4)] seeds = [50101, 50102, 50103, 50104, 50105] vectors: dict[int, list[np.ndarray]] = {i: [] for i in range(4)}; rows = [] for oi, targets in enumerate(held): rates = np.full(len(targets), 50.0, np.float32) for seed in seeds: run = sim.run(targets, rates, seed + oi * 100, 100.0) active = int(np.count_nonzero(run["kc_counts"])); fraction = active / len(sim.kc) vectors[oi].append(run["kc_counts"]) rows.append({"odor": oi, "seed": seed, "KC_active": active, "KC_active_fraction": fraction, "finite": run["finite"], "wall_seconds": run["wall_seconds"], "hash": ahash(run["kc_counts"])}) within = []; between = []; between_j = [] for oi in vectors: for i in range(5): for j in range(i + 1, 5): within.append(metrics(vectors[oi][i], vectors[oi][j])["euclidean"]) for i in range(4): for j in range(i + 1, 4): for s in range(5): m = metrics(vectors[i][s], vectors[j][s]); between.append(m["euclidean"]); between_j.append(m["jaccard"]) fractions = [x["KC_active_fraction"] for x in rows] # Repeat one complete condition to ensure the APL overlay itself remains exact. x1 = sim.run(held[0], np.full(n, 50.0, np.float32), seeds[0], 100.0) x2 = sim.run(held[0], np.full(n, 50.0, np.float32), seeds[0], 100.0) repeat_exact = all(np.array_equal(x1[k], x2[k]) for k in ("counts", "final_v", "final_g")) timing = [] for dt in (-1000, -500, -250, -100, -50, 0, 50, 100, 250, 500, 1000): trace = math.exp(-abs(dt) / 100.0); polarity = 1 if dt < 0 else -1 timing.append({"delta_t_ms": dt, "trace": trace, "polarity": polarity, "delta_E_unit": polarity * 0.01 * trace}) sparse = 0.03 <= float(np.mean(fractions)) <= 0.09 and max(fractions) < 0.5 and all(x["finite"] for x in rows) separation = float(np.mean(between)) > float(np.mean(within)) and float(np.mean(between_j)) < 0.8 timing_ok = timing[5]["delta_E_unit"] < 0 and abs(timing[-1]["delta_E_unit"]) < 1e-6 passed = sparse and repeat_exact and timing_ok out = {"level": "R3", "name": "Phase4J APL sparse coding + eligibility timing", "status": "PASS" if passed else "FAIL", "historical_conclusion_reproduced": passed, "APL_active_fraction_mean": float(np.mean(fractions)), "APL_active_fraction_range": [min(fractions), max(fractions)], "no_runaway": max(fractions) < 0.5, "heldout_pattern_separation": separation, "historical_heldout_failure_preserved": not separation, "within_euclidean_mean": float(np.mean(within)), "between_euclidean_mean": float(np.mean(between)), "between_jaccard_mean": float(np.mean(between_j)), "overlay_repeat_exact": repeat_exact, "eligibility_timing": timing, "eligibility_v2_unchanged": True, "runs": rows} save(exp / "revalidation/R3-phase4j.json", out); return out def door_inputs() -> tuple[dict[str, Any], pd.DataFrame, dict[str, np.ndarray], dict[str, Any], pd.DataFrame]: cfg = yaml.safe_load((ROOT / "models/male-cns-lif/olfactory-grounding-v1/phase4l.yaml").read_text()) # phase4m.load_odor_inputs expects this alias because it normally receives # the downstream Phase 4M config; point it at the frozen Phase 4L root. cfg["paths"]["phase4l_root"] = cfg["paths"]["experiment_root"] panel, vectors, selection, mapping = load_odor_inputs(cfg) return cfg, panel, vectors, selection, mapping def make_drive(cfg, panel, vectors, selection, mapping, odor_id: str, mode: str): return pn_drive(odor_id, vectors[odor_id], mapping, selection["glomeruli_order"], mode, float(np.median(panel.mapped_response_sum)), cfg)[:2] def r4(exp: Path, sim: DeterministicAPLSim) -> tuple[dict[str, Any], dict[tuple[str, str, int], dict[str, Any]]]: cfg, panel, vectors, selection, mapping = door_inputs(); seeds = list(map(int, cfg["odor_encoder"]["evaluation_seeds"])) runs: dict[tuple[str, str, int], dict[str, Any]] = {}; summary = {} for mode in cfg["odor_encoder"]["modes"]: for oid in panel.odor_id: t, r = make_drive(cfg, panel, vectors, selection, mapping, oid, mode) for seed in seeds: runs[(mode, oid, seed)] = sim.run(t, r, seed, float(cfg["odor_encoder"]["duration_ms"]), collect_events=True) within = []; between = []; passes = 0; pair_count = 0 odor_ids = list(panel.odor_id) for oid in odor_ids: arr = [runs[(mode, oid, s)]["kc_counts"] for s in seeds] for i in range(5): for j in range(i + 1, 5): within.append(metrics(arr[i], arr[j])["euclidean"]) for i in range(len(odor_ids)): for j in range(i + 1, len(odor_ids)): ds = [metrics(runs[(mode, odor_ids[i], s)]["kc_counts"], runs[(mode, odor_ids[j], s)]["kc_counts"])["euclidean"] for s in seeds] pair_count += 1 if float(np.mean(ds)) > float(np.mean(within)): passes += 1 between.extend(ds) frac = [np.count_nonzero(run["kc_counts"]) / len(sim.kc) for (m, _, _), run in runs.items() if m == mode] summary[mode] = {"KC_active_fraction_mean": float(np.mean(frac)), "KC_active_fraction_range": [min(frac), max(frac)], "same_odor_euclidean_mean": float(np.mean(within)), "different_odor_euclidean_mean": float(np.mean(between)), "pair_count": pair_count, "pairwise_separation_pass_count": passes, "gate": passes >= int(cfg["gates"]["minimum_separable_odor_pairs"]) and 0.03 <= float(np.mean(frac)) <= 0.09} passed = all(x["gate"] for x in summary.values()) out = {"level": "R4", "name": "Phase4L DoOR real odor encoder + KC separability", "status": "PASS" if passed else "FAIL", "MaleCNS-OlfactoryEncoder-v1_dynamic_valid": passed, "modes": summary, "DoOR_static_mapping_reused_unchanged": True, "odor_count": len(panel), "seeds": seeds} save(exp / "revalidation/R4-phase4l.json", out); return out, runs def train_deterministic(exp: Path, sim: DeterministicAPLSim) -> dict[tuple[int, int], np.ndarray]: cfg, panel, vectors, selection, mapping = door_inputs(); mcfg = yaml.safe_load((ROOT / "models/male-cns-lif/kcmbon-efficacy-v1/phase4m.yaml").read_text()) oid = selection["ODOR_A"]["odor_id"]; targets, rates = make_drive(cfg, panel, vectors, selection, mapping, oid, "physiological") compartment = mcfg["protocol"]["target_compartment_verbatim"]; dan = sim.dan_targets(mcfg["protocol"]["matched_dan_type"], compartment) comp = sim.identity.compartment.astype(str).to_numpy() == compartment kc_lookup = np.full(sim.n, -1, np.int32); kc_lookup[sim.kc] = np.arange(len(sim.kc), dtype=np.int32) source_local = kc_lookup[sim.identity.source_neuron_idx.to_numpy(np.int32)] if np.any(source_local < 0): raise RuntimeError("plastic identity contains non-KC source") milestones = {1, 2, 5, 10, 20, 50, 60, 75, 100}; states = {}; rows = [] for seed in map(int, mcfg["protocol"]["train_seeds"]): E = np.ones(len(sim.edge_ids), np.float32) for pairing in range(1, 101): run = sim.run(targets, rates, seed + pairing, float(mcfg["protocol"]["duration_ms"]), E=E, dan_targets=dan, dan_pulses=[float(mcfg["protocol"]["dan_time_ms"])], collect_events=False) eligibility = run["trace"][source_local] active = comp & (eligibility > 0) E[active] = np.maximum(np.float32(0.75), np.float32(E[active] * np.float32(1.0 - np.float32(0.01) * eligibility[active]))) if pairing in milestones: states[(seed, pairing)] = E.copy() rows.append({"seed": seed, "pairing": pairing, "updated_edges": int(np.count_nonzero(E < 1)), "mean_depression": float(np.mean(1 - E)), "max_depression": float(np.max(1 - E)), "saturated_edges": int(np.count_nonzero(E <= 0.75)), "trace_hash": ahash(run["trace"])}) if pairing % 20 == 0: print(json.dumps({"R5_training_seed": seed, "pairing": pairing}), flush=True) pd.DataFrame(rows).to_parquet(exp / "revalidation/R5-deterministic-training-trajectory.parquet", index=False) np.savez_compressed(exp / "revalidation/R5-efficacy-states.npz", **{f"seed_{s}_pairing_{p}": x for (s, p), x in states.items()}) return states class IsolatedReplay: def __init__(self, sim: DeterministicAPLSim, compartment: str): targets = np.sort(sim.identity.loc[sim.identity.compartment.astype(str) == compartment, "target_neuron_idx"].unique().astype(np.int32)) mask = np.isin(sim.sources, sim.kc) & np.isin(sim.rt.indices_np, targets) self.ids = np.flatnonzero(mask).astype(np.int64); self.src = sim.sources[self.ids]; self.target = sim.rt.indices_np[self.ids] self.targets = targets; pos = {int(x): i for i, x in enumerate(targets)} self.local = np.asarray([pos[int(x)] for x in self.target], np.int32) self.base = sim.base[self.ids].cpu().numpy().astype(np.float32); self.delay = sim.rt.delay_steps; self.qlen = sim.rt.qlen self.dt = sim.dt; self.gd = sim.rt.gd; self.alpha = sim.rt.alpha; self.rest = sim.rt.rest; self.reset = sim.rt.resetv self.threshold = sim.rt.threshold; self.refsteps = sim.rt.ref_steps can = {int(e): i for i, e in enumerate(sim.edge_ids)} self.candidate = np.asarray([can.get(int(e), -1) for e in self.ids], np.int32) self.by_source = {int(s): np.flatnonzero(self.src == s) for s in np.unique(self.src)} def run(self, steps: np.ndarray, nodes: np.ndarray, E: np.ndarray | None, duration_ms: float) -> dict[str, Any]: w = self.base.copy(); valid = self.candidate >= 0 if E is not None: w[valid] = np.float32(w[valid] * E[self.candidate[valid]]) grouped = {int(s): np.sort(nodes[steps == s]) for s in np.unique(steps)} n = len(self.targets); q = np.zeros((self.qlen, n), np.float32); g = np.zeros(n, np.float32) v = np.full(n, self.rest, np.float32); ref = np.zeros(n, np.int32); counts = np.zeros(n, np.int32); first = np.full(n, -1, np.int32) delivered = np.float64(0) for step in range(int(round(duration_ms / self.dt))): slot = step % self.qlen; g = np.float32(np.float32(g * self.gd) + q[slot]); delivered += q[slot].sum(dtype=np.float64); q[slot].fill(0) ref = np.maximum(ref, 1) - 1; avail = ref == 0 v[avail] = np.float32(v[avail] + np.float32(self.alpha * np.float32(g[avail] - np.float32(v[avail] - self.rest)))) v[~avail] = self.reset; fired = avail & (v >= self.threshold); counts += fired; first[(first < 0) & fired] = step v[fired] = self.reset; ref[fired] = self.refsteps for source in grouped.get(step, np.empty(0, np.int32)): rel = self.by_source.get(int(source)) if rel is not None: dest = q[(step + self.delay) % self.qlen] for j in rel: dest[self.local[j]] = np.float32(dest[self.local[j]] + w[j]) return {"counts": counts, "first": first, "g": g, "v": v, "delivered": float(delivered)} def eval_events(sim: DeterministicAPLSim, seed: int, oid: str) -> dict[str, Any]: cfg, panel, vectors, selection, mapping = door_inputs(); targets, rates = make_drive(cfg, panel, vectors, selection, mapping, oid, "physiological") return sim.run(targets, rates, seed, 120.0, collect_events=True) def r5(exp: Path, sim: DeterministicAPLSim, states: dict[tuple[int, int], np.ndarray]) -> tuple[dict[str, Any], dict[str, Any]]: cfg, panel, vectors, selection, mapping = door_inputs(); mcfg = yaml.safe_load((ROOT / "models/male-cns-lif/kcmbon-efficacy-v1/phase4m.yaml").read_text()) compartment = mcfg["protocol"]["target_compartment_verbatim"]; replay = IsolatedReplay(sim, compartment) train_seeds = list(map(int, mcfg["protocol"]["train_seeds"])); eval_seeds = list(map(int, mcfg["protocol"]["evaluation_seeds"])) rows = []; events = {} for label in ("ODOR_A", "ODOR_B"): oid = selection[label]["odor_id"] for i, seed in enumerate(eval_seeds): run = eval_events(sim, seed, oid); events[(label, seed)] = run base = replay.run(run["kc_event_step"], run["kc_event_node"], None, 120.0) learned = replay.run(run["kc_event_step"], run["kc_event_node"], states[(train_seeds[i], 100)], 120.0) rows.append({"odor": label, "seed": seed, "train_seed": train_seeds[i], "spike_L1": int(np.abs(base["counts"] - learned["counts"]).sum()), "changed_MBON": int(np.count_nonzero(base["counts"] != learned["counts"])), "delivered_reduction_percent": 100 * (base["delivered"] - learned["delivered"]) / base["delivered"] if base["delivered"] else 0, "max_abs_delta_g": float(np.max(np.abs(base["g"] - learned["g"]))), "max_abs_delta_v": float(np.max(np.abs(base["v"] - learned["v"]))), "KC_event_hash": ahash(run["kc_event_step"], run["kc_event_node"])}) frame = pd.DataFrame(rows); frame.to_parquet(exp / "revalidation/R5-isolated-replay.parquet", index=False) a = frame[frame.odor == "ODOR_A"]; b = frame[frame.odor == "ODOR_B"] functional = int((a.spike_L1 > 0).sum()) >= 3 specific = float(a.spike_L1.mean()) > float(b.spike_L1.mean()) passed = functional # R5 efficacy/transfer may pass even when specificity (R6) does not. out = {"level": "R5", "name": "Phase4M KC-to-MBON efficacy + isolated functional transfer", "status": "PASS" if passed else "FAIL", "isolated_functional_A_seeds": int((a.spike_L1 > 0).sum()), "A_spike_L1_mean": float(a.spike_L1.mean()), "B_spike_L1_mean": float(b.spike_L1.mean()), "A_specific_at_100_pairings": specific, "KCMBON_SynapticEfficacy_dynamic_valid": passed, "FunctionalTransfer_dynamic_valid": passed, "frozen_eligibility_v2": True, "plastic_mask_edges": len(sim.edge_ids), "rows": rows} save(exp / "revalidation/R5-phase4m.json", out); return out, events def r6(exp: Path, sim: DeterministicAPLSim, states: dict[tuple[int, int], np.ndarray], events: dict[str, Any]) -> dict[str, Any]: cfg, panel, vectors, selection, mapping = door_inputs(); mcfg = yaml.safe_load((ROOT / "models/male-cns-lif/kcmbon-efficacy-v1/phase4m.yaml").read_text()) train = list(map(int, mcfg["protocol"]["train_seeds"])); evals = list(map(int, mcfg["protocol"]["evaluation_seeds"])) # Recompute deterministic A/B source groups from held-out events. active = {} for label in ("ODOR_A", "ODOR_B"): nodes = np.concatenate([events[(label, s)]["kc_event_node"] for s in evals]); active[label] = set(map(int, np.unique(nodes))) source = sim.identity.source_neuron_idx.to_numpy(np.int32) groups = np.full(len(source), "inactive", object) A = np.isin(source, list(active["ODOR_A"])); B = np.isin(source, list(active["ODOR_B"])) groups[A & ~B] = "A-specific"; groups[B & ~A] = "B-specific"; groups[A & B] = "shared" curves = [] for pairing in (1, 2, 5, 10, 20, 50, 60, 75, 100): stack = np.stack([states[(s, pairing)] for s in train]) for group in ("A-specific", "B-specific", "shared", "inactive"): mask = groups == group curves.append({"pairing": pairing, "group": group, "edges": int(mask.sum()), "mean_delta_E": float(np.mean(stack[:, mask] - 1)) if mask.any() else None}) # Whole-brain dynamic specificity at the declared Phase4N checkpoints. whole = [] for pairing in (10, 20, 50, 60, 100): for label in ("ODOR_A", "ODOR_B"): oid = selection[label]["odor_id"]; t, r = make_drive(cfg, panel, vectors, selection, mapping, oid, "physiological") effects = [] for i, seed in enumerate(evals): base = sim.run(t, r, seed, 120.0) post = sim.run(t, r, seed, 120.0, E=states[(train[i], pairing)]) effects.append(int(np.abs(base["mbon_counts"] - post["mbon_counts"]).sum())) whole.append({"pairing": pairing, "odor": label, "MBON_spike_L1_mean": float(np.mean(effects)), "functional_seeds": int(np.count_nonzero(effects))}) print(json.dumps({"R6_pairing": pairing}), flush=True) wf = pd.DataFrame(whole); cf = pd.DataFrame(curves) wf.to_parquet(exp / "revalidation/R6-wholebrain-specificity.parquet", index=False) cf.to_parquet(exp / "revalidation/R6-edge-group-curves.parquet", index=False) at50 = wf[wf.pairing == 50].set_index("odor").MBON_spike_L1_mean specificity = bool(at50.get("ODOR_A", 0) > at50.get("ODOR_B", 0)) shared_contamination = bool(abs(cf[(cf.pairing == 50) & (cf.group == "shared")].mean_delta_E.iloc[0]) > abs(cf[(cf.pairing == 50) & (cf.group == "inactive")].mean_delta_E.iloc[0])) # Historical Phase4N result was a partial audit/no candidate. Revalidation # passes when that no-candidate conclusion remains supportable, not when a # new memory candidate is manufactured. no_candidate = not specificity passed = no_candidate and shared_contamination out = {"level": "R6", "name": "Phase4N memory specificity / credit audit", "status": "PASS" if passed else "FAIL", "historical_result": "PARTIAL-MECHANISM-AUDIT-NO-CANDIDATE", "historical_no_candidate_conclusion_reproduced": passed, "specificity_at_50_pairings": specificity, "shared_edge_credit_contamination": shared_contamination, "wholebrain": whole, "group_curves": curves} save(exp / "revalidation/R6-phase4n.json", out); return out def main() -> None: ap = argparse.ArgumentParser(); ap.add_argument("--config", type=Path, required=True); args = ap.parse_args() cfg = yaml.safe_load(args.config.read_text()); exp = Path(cfg["paths"]["experiment_root"]); model = Path(cfg["paths"]["model_root"]) levels = []; r1, r2 = r1_r2(exp); levels += [r1, r2] if r1["status"] != "PASS" or r2["status"] != "PASS": save(exp / "revalidation/summary.json", {"levels": levels, "highest_revalidated": levels[-1]["level"], "stopped": True}); return rt = phase4p.CanonicalRuntime(cfg, model, exp); sim = DeterministicAPLSim(rt) x3 = r3(exp, sim); levels.append(x3) if x3["status"] != "PASS": save(exp / "revalidation/summary.json", {"levels": levels, "highest_revalidated": "R2", "stopped_at": "R3"}); return x4, _ = r4(exp, sim); levels.append(x4) if x4["status"] != "PASS": save(exp / "revalidation/summary.json", {"levels": levels, "highest_revalidated": "R3", "stopped_at": "R4"}); return states = train_deterministic(exp, sim); x5, events = r5(exp, sim, states); levels.append(x5) if x5["status"] != "PASS": save(exp / "revalidation/summary.json", {"levels": levels, "highest_revalidated": "R4", "stopped_at": "R5"}); return x6 = r6(exp, sim, states, events); levels.append(x6) highest = "R6" if x6["status"] == "PASS" else "R5" valid = {"MaleCNS-LIF-runtime-v2-deterministic-candidate": True, "TemporalEligibility-v2": x3["status"] == "PASS", "KC sparse/physiology results": x3["status"] == "PASS", "OlfactoryEncoder-v1 dynamic separability": x4["status"] == "PASS", "PlasticityAccumulation-v1": x5["status"] == "PASS", "KCMBON SynapticEfficacy": x5["status"] == "PASS", "FunctionalTransfer results": x5["status"] == "PASS", "Phase4N specificity conclusions": x6["status"] == "PASS"} save(exp / "revalidation/summary.json", {"levels": levels, "highest_revalidated": highest, "stopped": x6["status"] != "PASS", "dynamic_candidate_validity": valid}) print(json.dumps({"highest_revalidated": highest, "statuses": {x["level"]: x["status"] for x in levels}}, indent=2), flush=True) if __name__ == "__main__": main()