Regulator wiring: set-pool modes, inert-drug and conserved-participant switches (all default-neutral), specs/032-038

CLAUDE.md

@@ -332,6 +332,19 @@ env vars only override for benchmark sweeps:
   `containment.csv`; the solve reports `self_inhibitor_rule` (on / off /
   inert) and how many inhibitor slots were flagged. Measured numbers are in
   specs/023. specs/022 has the rule; specs/023 has the adoption.
+- `DS_DRUG_MODE=propagate` (default) | `inert` (specs/032): drug-derived
+  nodes (a Reactome `Drug`, a complex containing one, a set that is all drugs),
+  listed by the generator's `drugs.csv`, held at baseline like cofactors, for a
+  cell without the drug. The solve reports `drug_rule` and `drugs_held`.
+- `DS_SET_POOL_MODE=product` (default) | `extreme` | `geomean` | `mean` | `max`
+  (specs/033): how a set POOL node (the generator's `LNG_SET_POOL`: one node
+  per set-valued catalyst or regulator, fed by `set_member` edges) combines its
+  members. It has no effect on networks without `set_member` edges.
+- `DS_CONSERVED_MODE=off` (default) | `inert` (specs/035): a loop input released
+  UNCHANGED by its own loop (every producer consumes a complex containing it:
+  RAS GAPs, scaffold partners, F-actin) is held at baseline like a cofactor.
+  Signal carriers transformed inside the loop are not touched. The solve reports
+  `conserved_rule` and `conserved_held`.
 - Export aggregation default: `stoichiometry_weighted`.
 
 ### Where design decisions live
@@ -413,6 +426,27 @@ for the last feature that touched it rather than re-deriving from the code.
   reaction forms (a curation-review list). 96 entity → reaction breaks are copy
   multiplicity; the targeted join is untried, near a family that failed three
   times. About 41 cases for an any-copy pool. Analysis only.
+- `specs/032-inert-drugs/` — `DS_DRUG_MODE`: drug-derived nodes held at
+  baseline. Not adopted: held-out +30, but experimental −5, RAF worse.
+- `specs/033-set-pool-nodes/`, `specs/038-pool-max/` — one OR node per
+  set-valued catalyst or regulator (`LNG_SET_POOL`), feeding the reaction once.
+  Adoptable on faithfulness under `DS_SET_POOL_MODE=product`: held-out +6 / −2,
+  experimental +3 / +8 (the +8 concentrated in 4 perturbations) on two
+  measurements. `max` (any member suffices) is refuted (held-out −120,
+  experimental −57); `mean` does not gain.
+- `specs/034-regulator-wiring/` — inventory of how catalysts and regulators are
+  wired; the RAF collapse traced by iteration; MP-BioPath's network and its
+  1,800 hand edits compared with ours. UUIDs already resolve 86% of their
+  loop artefacts. Our loops are mostly curated interconversion cycles, and the
+  219 inhibitions they added are mostly the reversal enzyme of such cycles
+  (187 absent from ours). The design input for what comes next.
+- `specs/035-conserved-participants/` — `DS_CONSERVED_MODE`: loop inputs
+  released unchanged held at baseline. Not adopted: it held PIP3 itself, and
+  experimental was −122.
+- `specs/036-capped-bundles-as-pools/`, `specs/037-combined-raf-fixes/` — the
+  variant cap bundles a set's alternatives into one all-required node
+  (`LNG_CAP_POOLS` fixes the structure, default off). Combined with de-duplication
+  it was held-out −162; de-duplication is the damage. Not adopted.
 - `specs/009-solver-defaults/` — the one-variable-at-a-time re-measurement
   behind the `DS_*` defaults above (cited in that section too).
 - `.specify/memory/constitution.md` — project principles the specs are
@@ -426,19 +460,22 @@ behaviour, not the whole directory; `ls specs/` is the complete list.
 Per-feature numbers belong in that feature's `research.md`, not here.
 
 ### Testing Strategy
-**Most of the test suite cannot fail.** Thirteen files contain assertions; the
+**Most of the test suite cannot fail.** Sixteen files contain assertions; the
 other seven execute code and print output.
 
 | file | assertions | note |
 |---|---|---|
-| `test/test_config_validation.jl` | 196 | |
+| `test/test_config_validation.jl` | 201 | |
 | `test/test_loop_elasticity.jl` | 173 | + 1 `@test_broken` |
 | `test/test_propagator_invariants.jl` | 121 | |
 | `test/test_loop_pool.jl` | 81 | |
 | `test/test_solver_determinism.jl` | 80 | |
 | `test/test_and_curves.jl` | 71 | |
 | `test/test_scc_break_roles.jl` | 56 | |
 | `test/test_self_inhibition.jl` | 80 | specs/022, added 2026-09-25 |
+| `test/test_drugs.jl` | 59 | specs/032, added 2026-09-26 |
+| `test/test_set_pool.jl` | 43 | specs/033, added 2026-09-26 |
+| `test/test_conserved.jl` | 15 | specs/035, added 2026-09-26 |
 | `test/test_cli_observations.jl` | 39 | |
 | `test/test_api_errors.jl` | 44 | |
 | `test/test_cycle_handling.jl` | 34 | + 2 `@test_broken` |

bench/analysis/regulator_inventory.py

@@ -0,0 +1,195 @@
+#!/usr/bin/env python3
+"""Inventory of how catalysts and regulators are wired into a catalog build
+(specs/034). Every catalyst/regulator EDGE into a reaction node is classified on
+four axes:
+
+  role    catalyst | positive | negative | requirement | pos_gene_expr | neg_gene_expr
+          (from the Regulation subclass in Neo4j; the generator only knows
+          catalyst / PositiveRegulation / NegativeRegulation)
+  shape   of the curated entity: protein | complex | complex_with_set | set |
+          set_of_complexes | small_molecule | drug-derived (drugs.csv when present)
+  copy    of the SOURCE node that is wired in:
+            root          no incoming edge (what the benchmark pins)
+            upstream      produced, and the reaction does not reach it
+            downstream    produced, and the reaction reaches it (a loop closes)
+            pool          an LNG_SET_POOL node
+  combine the edge's and_or flag
+
+With --cases, each curator/experimental case is attributed to the regulator
+edges its perturbed gene's pinned nodes can reach through a regulator edge, and
+accuracy is tallied per (role, shape, copy) cell.
+
+Usage: regulator_inventory.py --catalog BUILD_DIR [--cases ARM_DIR]
+"""
+from __future__ import annotations
+
+import argparse
+import collections
+import csv
+import sys
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+ROLE_Q = """
+UNWIND $rids AS rid MATCH (r:ReactionLikeEvent {stId: rid})
+CALL {
+  WITH r MATCH (r)-[:catalystActivity]->(:CatalystActivity)-[:physicalEntity]->(e)
+  RETURN 'catalyst' AS role, e.stId AS e
+  UNION WITH r MATCH (r)-[:regulatedBy]->(rg)-[:regulator]->(e)
+  RETURN CASE WHEN rg:Requirement THEN 'requirement'
+              WHEN rg:PositiveGeneExpressionRegulation THEN 'pos_gene_expr'
+              WHEN rg:NegativeGeneExpressionRegulation THEN 'neg_gene_expr'
+              WHEN rg:PositiveRegulation THEN 'positive'
+              ELSE 'negative' END AS role, e.stId AS e
+}
+RETURN rid AS r, role, e
+"""
+
+SHAPE_Q = """
+UNWIND $ids AS s MATCH (e:PhysicalEntity {stId: s})
+RETURN s,
+  CASE WHEN e:SimpleEntity THEN 'small_molecule'
+       WHEN e:EntitySet AND size([(e)-[:hasMember|hasCandidate]->(m:Complex) | m]) > 0 THEN 'set_of_complexes'
+       WHEN e:EntitySet THEN 'set'
+       WHEN e:Complex AND size([(e)-[:hasComponent*1..6]->(x:EntitySet) | x]) > 0 THEN 'complex_with_set'
+       WHEN e:Complex THEN 'complex'
+       ELSE 'protein' END AS shape,
+  [(e)-[:hasComponent|hasMember|hasCandidate*0..8]->(l:PhysicalEntity) | l.stId] AS leaves
+"""
+
+
+def main() -> int:
+    ap = argparse.ArgumentParser()
+    ap.add_argument("--catalog", type=Path, required=True)
+    ap.add_argument("--cases", type=Path, help="arm dir with {curator,experimental}_cases.tsv")
+    ap.add_argument("--cell", help="role/shape/copy: list its edges and the FAILED cases routed through them")
+    a = ap.parse_args()
+    import benchmark_vs_mpbiopath as BM
+    from py2neo import Graph
+    g = Graph(BM.NEO4J_URL, auth=(BM.NEO4J_USER, BM.NEO4J_PASSWORD))
+
+    names = {}
+    for line in open(Path(__file__).resolve().parents[1] / "catalog_pathways.tsv"):
+        if line.startswith(("#", "id\t")):
+            continue
+        sid, name = line.rstrip("\n").split("\t")
+        names[sid] = name
+
+    cells = collections.Counter()          # (role, shape, copy, and_or) -> edges
+    want = tuple(a.cell.split("/")) if a.cell else None
+    cell_edges = collections.Counter()     # (pathway, regulator stid, reaction stid) -> failed cases
+    case_cells = collections.defaultdict(lambda: [0, 0])   # (role, shape, copy) -> [n, correct]
+    cases = {}
+    if a.cases:
+        for ax in ("curator", "experimental"):
+            p = a.cases / f"{ax}_cases.tsv"
+            if p.exists():
+                for r in csv.DictReader(open(p), delimiter="\t"):
+                    if r["valid"] == "1" and not r.get("exclusion"):
+                        cases.setdefault(r["pathway"], []).append((ax, r))
+
+    for d in sorted(a.catalog.glob("R-HSA-*")):
+        pid = d.name
+        lab = {x["uuid"]: x["stable_id"] for x in csv.DictReader(open(d / "stid_to_uuid_mapping.csv"))}
+        kinds = {x["uuid"]: x["node_kind"] for x in csv.DictReader(open(d / "nodes.csv"))}
+        drugs = set()
+        if (d / "drugs.csv").exists():
+            drugs = {x["stable_id"] for x in csv.DictReader(open(d / "drugs.csv"))}
+        E = list(csv.DictReader(open(d / "logic_network.csv")))
+        fwd = collections.defaultdict(set)
+        indeg = collections.Counter()
+        for e in E:
+            fwd[e["source_id"]].add(e["target_id"])
+            indeg[e["target_id"]] += 1
+        reach_cache: dict = {}
+
+        def reach(u):
+            if u not in reach_cache:
+                seen, st = {u}, [u]
+                while st:
+                    for t in fwd[st.pop()]:
+                        if t not in seen:
+                            seen.add(t)
+                            st.append(t)
+                reach_cache[u] = seen
+            return reach_cache[u]
+
+        reg_edges = [e for e in E if e["edge_type"] in ("catalyst", "regulator")]
+        rx_stids = sorted({lab.get(e["target_id"], "") for e in reg_edges} - {""})
+        roles = collections.defaultdict(list)       # (reaction stid) -> [(role, entity stid)]
+        for r in g.run(ROLE_Q, rids=rx_stids).data():
+            roles[r["r"]].append((r["role"], r["e"]))
+        ent = sorted({e for v in roles.values() for _, e in v})
+        shape, leaves = {}, {}
+        for r in g.run(SHAPE_Q, ids=ent).data():
+            shape[r["s"]], leaves[r["s"]] = r["shape"], set(r["leaves"])
+        edge_cell = {}
+        for e in reg_edges:
+            src, tgt = e["source_id"], e["target_id"]
+            rstid, sstid = lab.get(tgt, ""), lab.get(src, "")
+            want_neg = e["pos_neg"] == "neg"
+            # the curated annotation this edge comes from: the entity equal to,
+            # or containing, the source node's entity, with a matching sign
+            cands = [(ro, en) for ro, en in roles.get(rstid, [])
+                     if (ro in ("negative", "neg_gene_expr")) == want_neg
+                     and (en == sstid or sstid in leaves.get(en, set()))]
+            if e["edge_type"] == "catalyst":
+                cands = [c for c in cands if c[0] == "catalyst"] or cands
+            else:
+                cands = [c for c in cands if c[0] != "catalyst"] or cands
+            role, en = cands[0] if cands else ("unmatched", sstid)
+            shp = "drug" if (en in drugs or sstid in drugs) else shape.get(en, "?")
+            if kinds.get(src) == "set_pool":
+                copy = "pool"
+            elif indeg[src] == 0:
+                copy = "root"
+            elif src in reach(tgt):
+                copy = "downstream"
+            else:
+                copy = "upstream"
+            cells[(role, shp, copy, e["and_or"])] += 1
+            edge_cell[(src, tgt)] = (role, shp, copy)
+        # attribute cases: a case uses a regulator edge if its pins reach the source
+        for ax, r in cases.get(names.get(pid, ""), []):
+            pins = [u for u in r["gene_uuids"].split("|") if u]
+            if not pins:
+                continue
+            got = set()
+            seen, st = set(pins), list(pins)
+            while st:
+                for t in fwd[st.pop()]:
+                    if t not in seen:
+                        seen.add(t)
+                        st.append(t)
+            outs = set(u for u in r["output_uuids"].split("|") if u)
+            ok = r["predicted"] == r["expected"]
+            for (s, t), c in edge_cell.items():
+                if s in seen and outs & reach(t):
+                    got.add(c)
+                    if want and c == want and not ok:
+                        cell_edges[(ax, names.get(pid, pid)[:30], r["gene"], r["direction"],
+                                    lab.get(s, s), lab.get(t, t))] += 1
+            for c in got:
+                key = (ax,) + c
+                case_cells[key][0] += 1
+                case_cells[key][1] += ok
+
+    print("role / shape / copy / and_or -> edges")
+    for k, v in sorted(cells.items(), key=lambda kv: -kv[1]):
+        print(f"  {v:7d}  {' / '.join(k)}")
+    if want:
+        print(f"\nFAILED cases through {a.cell} (axis, pathway, gene, dir, regulator node, reaction):")
+        for k, v in cell_edges.most_common(25):
+            print(f"  {v:4d}  {k}")
+        return 0
+    if case_cells:
+        print("\ncases whose route uses a regulator edge of the cell (a case may use many cells):")
+        for k, (n, ok) in sorted(case_cells.items(), key=lambda kv: -kv[1][0]):
+            if n >= 20:
+                print(f"  {k[0]:12} {' / '.join(k[1:]):45} n={n:6d} acc={ok / n:.1%}")
+    return 0
+
+
+if __name__ == "__main__":
+    sys.exit(main())

bench/benchmark_mpbiopath_cases.py

@@ -260,6 +260,10 @@ def load_dbid_to_uuids(pathway_dir: Path) -> dict[str, list[str]]:
     if nodes_file.exists():
         with nodes_file.open(newline="") as handle:
             for row in csv.DictReader(handle):
+                # specs/033: a set pool is not a readout or pin target (see
+                # benchmark_vs_mpbiopath.load_stid_to_uuids).
+                if normalize_cell(row.get("node_kind")) == "set_pool":
+                    continue
                 stable_ids = {normalize_cell(row.get("diagram_entity_id"))}
                 stable_ids.update(normalize_cell(item) for item in normalize_cell(row.get("member_leaves")).split("|"))
                 for stable_id in stable_ids - {""}:

bench/benchmark_vs_mpbiopath.py

@@ -330,6 +330,14 @@ def load_stid_to_uuids(pathway_dir: Path):
         seen = set()
         with open(nodes_csv) as f:
             for row in csv.DictReader(f):
+                # specs/033: a set POOL node is the set as one catalyst/regulator,
+                # not a readout or a pin target. Indexed, it would capture set-
+                # valued readouts (p-ERK dimers, p-AKT) and join the mean of every
+                # member readout, changing what is MEASURED rather than the model
+                # (review of the set-pool change). Readouts and pins stay as in
+                # an unpooled build.
+                if (row.get("node_kind") or "") == "set_pool":
+                    continue
                 uuid = str(row["uuid"])
                 keys = set()
                 de = (row.get("diagram_entity_id") or "").strip()
@@ -889,10 +897,19 @@ def run_pathway(pathway_id: str, pathway_name: str, gene_to_stids_cache=None,
     # Record what was actually pinned, so a run's protocol is on the record
     # rather than inferred from its flags (specs/023: 89% of pins had silently
     # become mid-pathway complexes for two months).
+    # specs/032: a pin on a drug node would override its inert hold and carry
+    # the perturbation through the drug, so count them (expected 0).
+    drug_uuids = set()
+    drugs_csv = pathway_dir / "drugs.csv"
+    if drugs_csv.exists():
+        with open(drugs_csv) as f:
+            for r in csv.DictReader(f):
+                drug_uuids.update(stid_to_uuids.get(r["stable_id"], []))
     for us in gene_to_uuids.values():
         PIN_TALLY["perturbations"] += bool(us)
         PIN_TALLY["pinned"] += len(set(us))
         PIN_TALLY["pinned_roots"] += sum(1 for u in set(us) if indeg[u] == 0)
+        PIN_TALLY["pinned_drugs"] += sum(1 for u in set(us) if u in drug_uuids)
 
     # Resolve key_output dbIds to their REAL stIds (may be R-ALL-, R-NUL-, not
     # just R-HSA-) so species that are genuinely in the network are found.
@@ -1002,7 +1019,10 @@ def run_pathway(pathway_id: str, pathway_name: str, gene_to_stids_cache=None,
     # bundle's cofactor list (the solve would report self_inhibitor_rule inert).
     network_payload = {"nodes": parsed["nodes"], "edges": edges, "pathways": parsed["pathways"],
                        "containment": parsed.get("containment"),
-                       "cofactor_stids": parsed.get("cofactor_stids")}
+                       "cofactor_stids": parsed.get("cofactor_stids"),
+                       # specs/032: without it an edge-drop arm under DS_DRUG_MODE=inert
+                       # reached the solver with no drug table (review of the branch).
+                       "drug_stids": parsed.get("drug_stids")}
     # Use the server-cached network by id (fast: send only observations per
     # solve). But if SKIP_EDGE_TYPES modified the edges above, the cached
     # network is stale, so send the full modified payload instead.
@@ -1041,9 +1061,16 @@ def run_pathway(pathway_id: str, pathway_name: str, gene_to_stids_cache=None,
             raise SystemExit(f"{pathway_name}: the solve reports self_inhibitor_rule="
                              f"{ds_result['self_inhibitor_rule']!r} -- the network reached the solver "
                              "without its containment table, so this is not the default model.")
+        if ds_result and ds_result.get("drug_rule") == "inert: no drug table":
+            raise SystemExit(f"{pathway_name}: DS_DRUG_MODE=inert but the network has no drugs.csv "
+                             "(a build that predates it), so no drug is held (specs/032).")
         if ds_result:
             SOLVE_TALLY["solves"] += 1
             SOLVE_TALLY["converged"] += bool(ds_result.get("converged"))
+            SOLVE_TALLY[f"drug_rule={ds_result.get('drug_rule', 'unknown')}"] += 1
+            SOLVE_TALLY["drugs_held"] += int(ds_result.get("drugs_held") or 0)
+            SOLVE_TALLY[f"conserved_rule={ds_result.get('conserved_rule', 'unknown')}"] += 1
+            SOLVE_TALLY["conserved_held"] += int(ds_result.get("conserved_held") or 0)
         activities = ds_result.get("node_activities", {}) if ds_result else {}
 
         for r in rows:
@@ -1330,8 +1357,18 @@ def main():
                     f"{r['total']}\t{r['correct']}\t{r['accuracy']:.6f}\t"
                     f"{r['valid_total']}\t{r['valid_correct']}\t{r['valid_accuracy']:.6f}\n")
     print(f"\nConverged: {SOLVE_TALLY['converged']} of {SOLVE_TALLY['solves']} solves")
+    rules = sorted(k.split("=", 1)[1] for k in SOLVE_TALLY if k.startswith("drug_rule="))
+    print(f"\nDrugs: rule {'/'.join(rules) or 'unknown'}, {SOLVE_TALLY['drugs_held']} node-solves held (specs/032)")
+    if "inert" in rules and SOLVE_TALLY["drugs_held"] == 0:
+        # An inert arm that held nothing measured the default model.
+        raise SystemExit("DS_DRUG_MODE=inert but no solve held a drug node; the arm is the control.")
+    crules = sorted(k.split("=", 1)[1] for k in SOLVE_TALLY if k.startswith("conserved_rule="))
+    print(f"\nConserved: rule {'/'.join(crules) or 'unknown'}, {SOLVE_TALLY['conserved_held']} node-solves held (specs/035)")
+    if "inert" in crules and SOLVE_TALLY["conserved_held"] == 0:
+        raise SystemExit("DS_CONSERVED_MODE=inert but no solve held a conserved node; the arm is the control.")
     print(f"\nPinned: {PIN_TALLY['pinned']} nodes over {PIN_TALLY['perturbations']} "
-          f"perturbations, {PIN_TALLY['pinned_roots']} of them roots (DS_PIN_SCOPE={PIN_SCOPE})")
+          f"perturbations, {PIN_TALLY['pinned_roots']} of them roots, "
+          f"{PIN_TALLY['pinned_drugs']} drug nodes (DS_PIN_SCOPE={PIN_SCOPE})")
     print(f"\nPer-pathway report: {args.report}")
 
     if args.dump_cases: