# NOTES — seahorse / AD network-pharmacology null-model control

## What actually happened, in order

1. **Environment check.** `/tmp` already had a warm cache from a prior probe of this task
   (`ctd.tsv.gz`, `string.gz`, `kegg.gmt`, `dis.gmt`, `gwas.gmt` — see timestamps ~11:xx today).
   Verified each file's header/format by hand (`zcat | head`) before writing any code, since the
   plan's column indices (CTD col 8 = OrganismID) and GMT formats needed to match reality, not
   assumption.
2. One file was missing from cache: `9606.protein.info.v12.0.txt.gz` (ENSP→gene-symbol map,
   needed to turn STRING edges into gene symbols). Downloaded fresh in 2s — network egress works
   fine from this sandbox, contrary to the plan's implicit assumption that everything was already
   cached.
3. **Format gotcha in DisGeNET's GMT** (`dis.gmt`): unlike KEGG/GWAS Enrichr GMTs (`term \t \t
   gene1 \t gene2 ...`), DisGeNET's has an extra category token at position 2 (`disease` or
   `group`) before the gene list starts. Handled by filtering those two literal tokens out of the
   parsed gene list (`GMT_CATEGORY_TOKENS` in `run_experiment.py`) rather than hard-coding a
   different column offset — simpler and safe against future format drift.
4. **Disease-swap term matching.** The plan named diseases in plain English ("Hypertension",
   "Osteoarthritis", ...) but DisGeNET doesn't have those exact strings. Grepped `dis.gmt` for
   the closest UMLS-style term for each (`Diabetes Mellitus, Non-Insulin-Dependent` for T2D,
   `Hypertensive disease`, `Generalized osteoarthritis`, `Colorectal Carcinoma`, `Asthma` — the
   last two matched verbatim). Logged the mapping in `DISEASE_SWAP_TERMS` and printed it at
   runtime so it's auditable.
5. **Design decision: which arm builds the PPI subnetwork from which gene set.** The plan's
   baseline description ("CTD targets → intersect AD geneset → STRING subnetwork → top-10 hub ...
   → KEGG+GO hypergeometric") is ambiguous about whether the *subnetwork* is built from the full
   herb target set or the herb∩disease intersection. Went with the standard TCM-network-pharm
   convention (subnetwork on the intersection, enrichment on the full target set) since that's
   what published seahorse/AD papers actually do — this is exactly the pipeline being tested for
   artifacts, so it needs to match real practice, not a convenient simplification.
6. **Dropped GO_BP_2025** entirely (never downloaded). None of the 6 pre-registered metrics
   reference it — they only need the KEGG hsa05010 term, hub gene names, overlap counts, and
   disease ranks. Keeping it out saved a nontrivial download+compute cost (GO BP libraries
   typically carry several thousand terms) for zero effect on any reported number. Recorded as a
   deviation.
7. **Speed surprise.** Assumed 1000 null draws × 2 null models × 3 seeds would need to be cut
   down to fit 40 minutes (per-iteration cost looked worrying on paper: STRING subnetwork +
   hypergeometric test). In practice each draw is a handful of set operations on already-loaded
   in-memory dicts — a full null arm of 1000 draws takes ~0.4s (size-matched) to ~3.7s
   (degree-matched, slower because of the log2-bin lookup/fallback logic). Ran the **full planned
   scale** (1000 draws × 3 seeds × 2 null models = 6000 draws) in under 20 seconds, no reduction
   needed. This freed up budget to also run a supplementary 3×3 sweep over
   `string_confidence_cutoff` × `ad_geneset_source` (reduced to n_null=300, seed=0 only, since the
   full 1000×3-seed grid across 9 combinations was not attempted — see `deviations` in
   results.json and `results_raw/sweep_results.csv`).
8. **Sanity check before trusting the "confirmed" outcome:** looked at the actual gene-level hub
   lists, not just the aggregate rate. Both null arms recover AKT1/TNF/IL6 at similarly high rates
   as any real herb would (these three are extremely high-degree "hub of everything" genes in
   STRING — a well-known promiscuity artifact), while JUN/PPARG/TP53 are recovered ~0% by either
   real seahorse data or the null (they're just not in the intersection set often, real or random,
   at this confidence cutoff) and APP is recovered only slightly more often by chance (~3-9%) than
   its low base rate would suggest. This is consistent with the aggregate numbers, not an
   artifact of how `hub_recovery_rate` is defined.
9. Total wall-clock for `bash run.sh` end-to-end (downloads/cache-copies + prefilter + primary run
   + sweep + figures): **49 seconds**, comfortably inside the 40-minute budget the plan allowed —
   the actual bottleneck the plan anticipated (STRING's ~13.7M raw edge rows) is handled by a
   4-second `awk` prefilter to `combined_score>=400` before anything touches pandas/Python dicts.

## Result summary (see results.json for full numbers)

- **Baseline (real seahorse compound list, n=20, 19 with CTD human targets after excluding
  Proline which has zero CTD human records):** 1270 deduped CTD target genes, 71 overlap with
  KEGG hsa05010 (Alzheimer disease), hypergeometric p=3.2e-7 (FDR=1.2e-6) — looks like a strong,
  "publishable" hit in isolation.
- **Null-degree-matched (primary/strong null, 1000 draws × 3 seeds):** random compound sets of
  matched size and matched per-compound CTD target-degree reproduce a top-10 hub that intersects
  {APP,JUN,PPARG,AKT1,TNF,IL6,TP53} **97.3%** of the time (vs. the pre-registered ≥50% confirm
  threshold). The seahorse overlap (71) and KEGG p-value are **not** significant against this null
  (emp_p_ad_overlap=0.478, emp_p_ad_kegg=0.403 — both far above 0.05).
- **Disease-swap control:** with the seahorse target set held fixed, AD ranks **4th of 6** by
  enrichment strength — Type 2 Diabetes, Colorectal Carcinoma and Hypertension all show *stronger*
  enrichment than AD in the exact same target set. The pipeline has essentially no
  disease-specificity for this compound list.
- **Sweep (string cutoff ∈ {400,700,900} × AD source ∈ {KEGG, DisGeNET, GWAS}, reduced N=300,
  seed 0):** hub_recovery_rate ranges 64%–100% across all 9 combinations, always well above the
  20% refute threshold — the "confirmed" (artifact) call is not an artifact of one particular
  cutoff or gene-set choice.

**Falsifiable prediction outcome: CONFIRMED.** The pre-registered refute condition (hub
recovery<20% AND both emp_p<0.05) did not occur in any seed or sweep cell. This means: under this
pipeline, using these databases, the published seahorse↔AD core-target/pathway narrative
(APP/JUN/PPARG-style "core module") is statistically indistinguishable from what any random
size/degree-matched compound set from CTD would produce. That is a negative result for the
pipeline's discriminative power — not a claim about whether seahorse has real anti-AD
pharmacology (no cell/animal/binding data were used here, per `out_of_scope`).

## What would change this conclusion

- A different, non-CTD/STRING/Enrichr target-prediction/PPI/pathway toolchain might behave
  differently — this result is a statement about *this specific standard pipeline's* database
  biases, not about seahorse biology.
- The peptide/protein-class actives of seahorse (not covered by CTD, which is small-molecule/gene
  interactions only) are entirely out of scope here and could carry a real, undetected signal.
- If the diagnosis is "AKT1/TNF/IL6 are generic high-degree STRING hubs contaminating every
  herb's hub list," a plausible next step (not run here, budget exhausted) would be to re-run hub
  ranking with STRING degree normalized against genome-wide background degree (e.g. z-score vs.
  random gene sets of the same size) instead of raw within-subnetwork degree — this might recover
  specificity that raw degree ranking washes out. Flagged as a follow-up, not attempted.

## Files

- `run_experiment.py` — primary pipeline (baseline + both null arms + disease-swap), argparse CLI.
- `sweep.py` — supplementary cutoff × AD-source robustness sweep, imports `run_experiment` as a
  module to reuse the CTD/STRING loading and pipeline functions without re-implementing them.
- `make_figures.py` — 4 figures from `results.json` + `results_raw/*.csv`.
- `run.sh` — one command, downloads-or-reuses-cache → prefilters STRING → runs primary → runs
  sweep → makes figures → writes `logs/run.log` (via `exec > >(tee ...)` at the top of the script,
  so the log is captured regardless of how run.sh is invoked).
- `results_raw/` — every intermediate CSV/JSON referenced by results.json's numbers (raw null
  draws per seed, KEGG enrichment table, disease-swap table, baseline hub list).