Metagenome-Prioritized Phage Cocktails for Crohn's Disease and IBD

1. Four reproducible IBD ecotypes with clear disease stratification

Training on 8,489 MetaPhlAn3 samples (fact_taxon_abundance, CMD_HEALTHY + CMD_IBD cohorts) with two independent methods — LDA on pseudo-counts and GMM on CLR + PCA-20 — across K ∈ {2..8}. Per-method fit measures (LDA held-out perplexity, GMM BIC) monotonically decrease with K, as expected for flexible latent-factor models. The discriminating signal is cross-method ARI between LDA and GMM, which has a local maximum at K = 4 (ARI = 0.131) and a second peak at K = 7 (0.140). A parsimony rule — smallest K within 0.02 ARI of the maximum — selects K = 4. Per-sample method agreement at K = 4 is 48.9 %.

The four consensus ecotypes are biologically coherent:

This is H1a directionally supported: ≥ 3 reproducible ecotypes. E0 / E1 / E2 / E3 map recognizably onto the original Bacteroides / Prevotella / Ruminococcus enterotype framework (Arumugam 2011, Costea 2018), with E1 / E3 reflecting the Bacteroides2 (Bact2) low-cell-count dysbiosis signature documented in CD by Vandeputte et al. 2017.

(Notebooks: NB01_ecotype_training.ipynb, NB01b_ecotype_refit.ipynb)

2. UC Davis CD patients span three ecotypes, none in E2

All 26 Kuehl_WGS samples (23 unique patients) projected onto the K = 4 reference via the synonymy layer. 262 unique Kaiju-classified species normalized to 97 canonical species in the training feature space. UC Davis distributes:

• E0 — diverse commensal: 7 samples (27 %)
• E1 — Bacteroides2 transitional: 11 samples (42 %)
• E2 — Prevotella copri enterotype: 0 samples
• E3 — severe Bacteroides-expanded: 8 samples (31 %)

χ²(3) vs uniform = 10.0, p = 0.019. The distribution is non-random. UC Davis looks Western (no E2 Prevotella-dominant patients), with active disease dominating (73 % E1 or E3). H1b directionally supported — patients distribute across multiple ecotypes rather than concentrating in one, validating the stratified-targeting premise of the project.

Longitudinal patients: 1112 → E3 at both timepoints; p1/p1reseq → E3 both; p2 → E1 both; 1460 (calprotectin 7,280 μg/g) → E1; patient 6967 flips E1 ↔ E3 between two samples. The 6967 finding is the first direct observation of intra-patient ecosystem instability — relevant to Pillar 5 H5d (dosing-schedule implications).

(Notebook: NB02_ecotype_projection.ipynb)

Methodological aside: Kaiju ↔ MetaPhlAn3 classifier mismatch asymmetry

Projecting Kuehl (Kaiju) onto a MetaPhlAn3-trained embedding exposes an asymmetric robustness between the two ecotype methods. LDA on pseudo-counts is robust: 54 % of Kuehl feature rows outside the training feature space is handled by treating absence as not-detected. GMM on CLR + PCA is fragile: the same sparsity forces all 26 Kuehl samples into a single Gaussian (E3) at confidence > 0.97 — an artifact, not biology. Documented as a project-level finding and committed to docs/discoveries.md. LDA is the primary Kuehl projection call; GMM is advisory.

3. Clinical covariates alone are insufficient for within-IBD ecotype assignment

Two classifiers trained on the pooled CMD cohort (LightGBM) to predict K = 4 consensus ecotype from clinical covariates:

• Minimal — {is_ibd, sex, age}, n = 8,489: macro OvR AUC = 0.799.
• Extended — adds {hbi_max, sccai_max, calp_max}, n = 1,675 subset: macro AUC = 0.810.

Both exceed the H1c threshold of 0.70. On paper, H1c passes. But applied to UC Davis patients:

• Minimal classifier vs NB02 metagenomic call: 41 % agreement (9/22).
• Extended classifier vs NB02: 36 % agreement (8/22).
• 12 / 22 patients disagree under both classifiers.

The minimal classifier predicts E1 for 19/22 UC Davis patients. In the training cohort, IBD samples split ~58 % E1 / ~40 % E3 / ~2 % E0 / ~0 % E2, so the classifier's dominant learned rule is "is_ibd = 1 → E1." When applied to UC Davis (all-CD, is_ibd constant), this rule collapses to the marginal mode. The extended classifier's training subset is 702 E1 / 959 E3 / 3 E0 / 11 E2 — effectively an E1-vs-E3 binary problem — and severity markers do not separate the two reliably.

H1c revised interpretation: clinical covariates distinguish HC vs IBD trivially (dominated by is_ibd) but do not separate IBD ecotypes (E1 transitional vs E3 severe). For UC Davis-type cohorts, metagenomics remains required for ecotype assignment. The "AUC 0.80 on paper / 41 % patient agreement in practice" gap is itself a methodologically important finding — OvR-AUC on a pooled cohort with a strong cohort-axis feature overstates per-patient classifier usefulness.

(Notebook: NB03_clinical_ecotype_classifier.ipynb)

4. Compositional correction partially, but not fully, resolves the C. scindens paradox

Starting observation: the preliminary project's pooled Mann-Whitney differential-abundance analysis called Clostridium scindens CD-enriched at log₂FC +2.67 — contradicting its established role as a bile-acid-producing protective species (~79 % prevalence in healthy individuals). Three explanations were possible: (1) compositional artifact, (2) strain heterogeneity, (3) ecotype mixing in the pooled analysis.

NB00 tests explanation (1) by running the same curated battery of 8 protective + 6 reported-pathobiont species under both raw Mann-Whitney on relative abundance and CLR-based Mann-Whitney (Gloor 2017, Lin & Peddada 2020). Findings:

• Compositional correction recovers depletion signal for 4+ protective species that raw Mann-Whitney misses: F. prausnitzii, A. muciniphila, R. hominis, L. eligens, A. rectalis flip from "n.s." (raw) to "CD↓" (CLR).
• Roseburia intestinalis shows a sign flip: raw CD↑ → CLR CD↓. The C. scindens artifact pattern reproduced on a second species.
• Top reported pathobionts agree across methods: M. gnavus (R. gnavus), E. bolteae, E. lenta, H. hathewayi all CD-enriched in both.
• C. scindens remains CD-enriched under both methods (raw log₂FC +5.66; CLR Δ +1.25). Compositional correction alone is insufficient; explanations (2) and (3) remain live.

Implication: the norm-N1 decision to re-analyze is justified. Pooled Mann-Whitney on relative abundance is systematically under-sensitive for protective-species depletion, an artifact of compositional bias plus pooled-cohort heterogeneity. C. scindens is not resolvable at this level — the resolution required a design that eliminates study confounding (see §5, NB04c's within-IBD-substudy CD-vs-nonIBD meta). Under that confound-free design C. scindens is CD↑, which is not a paradox but rather the expected behavior of a species that happens to be more prevalent in the IBD-cohort source studies than the healthy-cohort source studies in cMD.

(Notebook: NB00_data_audit.ipynb)

5. Within-ecotype × within-substudy meta-analysis defines ecotype-specific Tier-A (rigor-controlled)

Retraction and rigor repair — what NB04 claimed vs what this section presents

The original NB04 analysis (within-ecotype CD-vs-HC CLR Mann-Whitney, committed 2026-04-24 early, now superseded) made three headline claims, of which two are retracted here:

1. ~~H2c — the C. scindens paradox is resolved by within-ecotype stratification.~~ Retracted. Under the confound-free within-IBD-substudy CD-vs-nonIBD contrast, C. scindens is genuinely CD↑ (+1.18 CLR-Δ, FDR 1e-8, 4/4 sign-concordance across sub-studies); under leave-one-species-out refit, C. scindens is CD↑ within both E1 and E3 once it is not part of the clustering input. The NB04 within-ecotype n.s. call was a feature-leakage artifact (clustering samples on taxon abundance then testing the same taxon within cluster is selection-on-outcome confounding). The paradox was not a paradox.
2. H2b — target sets differ between ecotypes (Jaccard = 0.14). ~~Interpretation retained; statistic replaced.~~ Jaccard 0.14 was near the random-overlap baseline (~0.10 for top-30 of ~300 filtered species). Under a permutation null (200 random-label permutations, null mean 0.785 ± 0.054), the observed Jaccard of 0.104 has empirical p = 0.000 — H2b survives strongly and the divergence is real, but the claim in NB04 rested on an effect size without a reference distribution, which is the wrong argument structure.
3. ~~NB04 Tier-A list (33 species: 18 E1, 15 E3).~~ Retracted. The within-ecotype DA that produced this list was substantially driven by feature leakage (held-out-species sensitivity Jaccard: E1 = 0.230, E3 = 0.064, vs > 0.5 leakage-bounded threshold). When tested against an independent confound-free CD-vs-control contrast, 14 of 18 E1 candidates had negative within-substudy effects — they were ecotype-markers, not CD drivers.

The rigor-repair pipeline that produced the replacement analysis in this section is NB04b (bootstrap CIs + leakage-bound sensitivity + LOO refit + Jaccard permutation null + ecotype stability) → NB04c (confound-free within-IBD-substudy CD-vs-nonIBD meta + LinDA in pure Python) → NB04d (rigor-controlled stopping rule) → NB04e (within-ecotype × within-substudy CD-vs-nonIBD meta). See FAILURE_ANALYSIS.md for the full arc. Two generalizable pitfalls from this repair are in docs/pitfalls.md: cMD substudy-nesting unidentifiability and feature leakage in cluster-stratified DA.

5a. Confound-free design and why it works

curatedMetagenomicData pools samples across ≈ 51 source studies. In the ecotype-assigned slice (8,489 samples), 45 sub-studies have ≥ 10 HC and 5 sub-studies have ≥ 10 CD, but zero sub-studies contain both HC and CD. CMD's healthy-cohort samples come from HC-only studies (LifeLinesDeep_2016, AsnicarF_2021, YachidaS_2019, …); CMD-IBD CD samples come from IBD-cohort studies (HallAB_2017, VilaAV_2018, LiJ_2014, IjazUZ_2017, NielsenHB_2014). A pooled CD-vs-HC LME with substudy random effect is therefore structurally unidentifiable (the random effect perfectly predicts diagnosis) — we verified empirically that statsmodels.mixedlm silently fails to converge on every battery species under this design.

The confound-free contrast that is identifiable is CD-vs-nonIBD within IBD sub-studies — four cMD studies have ≥ 10 CD and ≥ 10 nonIBD (HallAB_2017: 89/73, LiJ_2014: 76/10, IjazUZ_2017: 56/38, NielsenHB_2014: 21/248). Within each sub-study, CLR-Δ is computed on the same samples processed the same way by the same group, with no study-level confound. Inverse-variance weighted meta-analysis across sub-studies combines into a cohort-level effect; stratifying the contrast by ecotype before meta-analysis produces the ecotype-specific effect. Because ecotype and CD-effect are computed on disjoint axes (samples partitioned by ecotype; CLR-Δ computed within each partition across sub-studies), there is no feature-leakage self-reference.

5b. NB04c cohort-level within-substudy CD-vs-nonIBD meta recovers the canonical CD signature

Across the four IBD sub-studies (242 CD / 369 nonIBD pooled), the 14-species curated battery produces the expected canonical CD signature — pathobionts up, protective commensals down — with strong sign concordance:

This flatly contradicts NB04's within-ecotype calls for several species. NB04 called F. prausnitzii / R. hominis / L. eligens CD↑ within both E1 and E3 (the "Simpson's paradox" of the original section 5) — the confound-free analysis shows they are CD↓, consistent with their classical protective-commensal role. NB04's C. scindens "n.s." within both ecotypes is similarly contradicted. The within-ecotype DA in NB04 was producing direction reversals as a compound artifact of feature leakage plus the pooled-cohort substudy × diagnosis confound — both compositional-bias-aware DA methods we tried on the within-ecotype subsets (CLR-MW and LinDA) share the bias, so n_evidence from within-ecotype methods alone does not resolve it.

5c. NB04e ecotype-specific Tier-A under within-ecotype × within-substudy meta

Stratifying the within-substudy CD-vs-nonIBD contrast by ecotype before meta-analysis tests whether the canonical CD signature differs by ecotype and produces ecotype-specific Tier-A lists that are structurally free of leakage. Three (substudy × ecotype) cells meet the ≥ 10 CD AND ≥ 10 nonIBD eligibility bar:

E0 and E2 are not viable — these are the healthy-cohort ecotypes; no IBD-study nonIBD samples live in them. E1 is meta-viable across two sub-studies; E3 is single-study-only (HallAB_2017).

E1 Tier-A (meta-analysis, 51 candidates, all 100 % sign-concordant across sub-studies, FDR < 0.10, pooled CLR-Δ > 0.5):

Full list: data/nb04e_within_ecotype_meta.tsv.

Biological coherence of the E1 Tier-A. The rigor-controlled E1 list organizes into three biologically coherent groups that the retracted NB04 list did not:

• Oral-derived streptococci as ectopic colonizers (ranks 2, 3, 5) — S. salivarius, S. thermophilus, S. parasanguinis are canonical oral-cavity species; their enrichment in CD stool is consistent with the "oral-gut axis" ectopic-colonization literature for IBD (Xiang 2024, PMID 39188957; Guo 2024, PMID 38545880; Tanwar 2023, PMID 37645044) and with the specific finding that S. salivarius is a salivary biomarker for orofacial granulomatosis co-occurring with CD (Goel 2019, PMID 30796823). Caveat: S. thermophilus appears in anti-inflammatory multi-strain probiotic formulations (Biagioli 2020, PMID 32629887) — its CD↑ signal at the species level does not automatically imply pathobiont status; strain-level evidence (A3 literature, A5 engraftment) is the NB05 disambiguation step.
• Vancomycin-resistant pathobionts (rank 4) — Erysipelatoclostridium innocuum is specifically documented as a vancomycin-resistant pathobiome in IBD with clinically consequential phenotypes (creeping-fat formation and intestinal strictures in CD; reduced UC remission rates); FMT is under evaluation as an intervention (Le 2025, PMID 40074633). Of the E1 Tier-A, E. innocuum is the candidate with the strongest stand-alone clinical-association evidence independent of this project's data.
• Clostridiales-expansion pathobionts (ranks 6–10) — Enterocloster asparagiformis, Intestinibacter bartlettii, Hungatella symbiosa, Erysipelatoclostridium ramosum, and the related clostridial reclassifications (E. bolteae, E. clostridioformis, E. citroniae in the lower-ranked list) are species in the Lachnospiraceae / Erysipelotrichaceae expansion characteristic of the Bacteroides-2 dysbiosis subtype (Vandeputte 2017) that defines E1. Their appearance in E1 Tier-A reflects the underlying ecology rather than being an ecotype-marker artifact — the within-substudy meta-analysis controls for the selection-on-outcome pattern that produced the original NB04 E1 list.
• Polyphenol-metabolism taxa as ambiguous CD-associated (rank 9) — Gordonibacter pamelaeae produces urolithins from dietary ellagitannins (Selma 2014, PMID 24744017) and increases during microbiome recovery from dysbiosis-inducing insult (Tierney 2023, PMID 36840551), which makes its CD↑ signal in E1 difficult to interpret as pathobiont activity. Flag for Tier-A-A4 protective-analog exclusion in NB05.

E3 Tier-A (provisional, 40 candidates, single-study HallAB_2017): top candidates H. symbiosa (+4.64), M. gnavus (+4.46), B. coccoides (+4.22), R. faecis (+4.14), C. spiroforme (+4.11), S. salivarius (+4.04), E. innocuum (+3.68). F. plautii replicates at +2.26 (FDR 0.02). Blautia wexlerae does not replicate (+0.25, FDR 0.80) — the NB04d "rock-solid E3 triad" that included B. wexlerae relied on NB04c's cohort-level within-substudy evidence rather than the E3-restricted within-substudy evidence; under the stricter E3 × HallAB_2017 test, B. wexlerae is removed from the rock-solid set.

The E3 list should be treated as provisional until a second cMD-IBD sub-study that populates E3 with ≥ 10 CD and ≥ 10 nonIBD samples becomes available (candidate: HMP2 once PENDING_HMP2_RAW is resolved).

5d. Classical pathobiont reality check (revised)

Under the confound-free design, the classical engraftment pathobionts from donor 2708 → P1 → P2 are unambiguously CD-enriched — the opposite of what NB04's within-ecotype analysis suggested:

5 of 6 engraftment pathobionts pass the confound-free CD↑ test. The NB04 within-ecotype "these pathobionts are ecotype-markers, not CD drivers" narrative is retracted — they are CD drivers under any analysis that controls for study confounding; the within-ecotype DA was systematically reversing their sign because the HC samples in each ecotype came from entirely different source studies than the CD samples in the same ecotype.

5e. Stopping rule and NB05 input

NB04d formalized a rigor-controlled stopping rule for NB05:

Criterion 1 fails in both ecotypes and documents that NB04's original within-ecotype Tier-A cannot be trusted directly — but NB04e's rigor-controlled within-ecotype × within-substudy meta is structurally immune to the leakage (clustering axis and DA axis are disjoint). Criterion 4 fails and documents that the ecotype framework is internally marginally stable (ARI 0.13–0.17 across 5 × 80 % subsamples); the framework is usable for downstream stratification (NB02 projection is deterministic once fit) but is not externally replicated and must be flagged.

Per-ecotype NB05 verdict:

• E1 — PROCEED with the 51-candidate meta-viable Tier-A. Confound-free, multi-study support, canonical-pathobiont enrichment at the top.
• E3 — PROCEED WITH CAVEAT using the 40-candidate single-study Tier-A. Replication is the first follow-up (HMP2 ingestion is the unblock).
• Cross-ecotype — the 5 engraftment-confirmed pathobionts are cross-ecotype candidates (NB04c §3).

(Notebooks: NB04_within_ecotype_DA.ipynb superseded; NB04b_analytical_rigor_repair.ipynb, NB04c_rigor_repair_completion.ipynb, NB04d_stopping_rule.ipynb, NB04e_option_A_viability.ipynb are the rigor-controlled pipeline.)

5f. Pillar 2 strengthening — LOSO stability, pathway-feature refit, and HMP2 external replication

Three additional analyses (NB04f, NB04g, NB04h) tested the ecotype framework and Tier-A claims against three distinct failure modes: (i) cross-study generalization, (ii) feature-leakage residual, (iii) external-cohort replication. The results are honest: the ecotype framework has real cross-study variance, but the operational Tier-A replicates strongly on an external cohort.

NB04f — Leave-one-substudy-out (LOSO) ecotype stability. For each of the top 8 cMD sub-studies by sample count (LifeLinesDeep_2016, AsnicarF_2021, NielsenHB_2014, VilaAV_2018, LiJ_2014, HallAB_2017, YachidaS_2019, HansenLBS_2018), held out the substudy, refit K=4 LDA on the remaining samples, projected held-out back, Hungarian-aligned, and computed ARI against the NB01b consensus_ecotype on the held-out samples. Mean LOSO ARI = 0.113 (range 0.000–0.282); mean per-sample agreement 55.0 % (range 16.9 %–85.5 %). This fails the > 0.30 "stable" threshold and is more honest than the bootstrap ARI 0.13–0.17 NB04b reported. Interpretation: some sub-studies (LifeLinesDeep ARI 0.21 / agreement 85.5 %, HansenLBS 65 %) align well with the consensus framework; others (AsnicarF 38 %, VilaAV 17 %) do not. The NB01b consensus was LDA+GMM at 48.9 % cross-method agreement, so part of the LOSO gap is the LDA↔consensus disagreement baseline. The within-cMD ecotype framework is cross-study variable; the "four reproducible ecotypes" framing must be qualified accordingly.

NB04g — Pathway-feature ecotype refit (Option B structural test). Refit K=4 LDA on fact_pathway_abundance (3,145 CMD_IBD samples, 2,000 top-variance HUMAnN3 pathways after prevalence + informative filtering) and compared to the taxon-based consensus_ecotype on the same samples. ARI = 0.113; per-sample agreement 50.6 %. Per-ecotype agreement: E1 = 65.3 % (n=1,839), E2 = 47.4 % (n=19), E3 = 30.7 % (n=1,244), E0 = 0 % (n=43; E0 is the healthy-dominant ecotype, rare in the CMD_IBD pathway cohort). E1 taxon-ecotype is substantially recoverable from a disjoint (pathway) feature basis — that's the meaningful result. E3 is less recoverable (30.7 %), which combined with NB04e's single-study E3 evidence supports the provisional E3 Tier-A framing. The ecotype structure is mixed ecological + taxonomic, not purely one or the other. Scope limitation: fact_pathway_abundance is CMD_IBD only; no HC pathway coverage, so the refit is within-disease only.

NB04h — HMP_2019_ibdmdb (HMP2) external replication. Pulled HMP2 MetaPhlAn3 profiles directly from curatedMetagenomicData v3.18 (1,627 samples, 130 subjects, 582 species; 255/335 training-feature overlap after synonymy). HMP2 is explicitly NOT in our CMD_IBD training set (HMP_2019_ibdmdb was absent from the CMD_IBD substudy list), so it is a genuinely held-out cohort with the same MetaPhlAn3 classifier namespace as the training data.

• Ecotype distribution replicates directionally: HMP2 subject-level (n=130) concentrates in E1 (106 subjects, 82 %), E2 (11), E3 (10), E0 (3). Subject-level χ² for ecotype × {CD, UC, nonIBD} = 15.61, p = 0.016 — ecotype stratifies disease in HMP2 at statistical significance. HMP2 skews more toward E1 than cMD does (cMD is more balanced E1+E3); this likely reflects HMP2's recruitment of newly-diagnosed / milder IBD rather than flare-dominated severe disease.
• Projection confidence is high: median max LDA posterior = 0.861; 80.4 % of samples have max posterior > 0.70. No Kaiju↔MetaPhlAn3 fragility (unlike the UC Davis GMM projection) because HMP2 uses the same MetaPhlAn3 pipeline as the training data.
• E1 Tier-A replicates strongly: per-species CD-vs-nonIBD CLR-Δ within HMP2-projected E1 samples (593 CD / 337 nonIBD), cross-referenced against the 51-candidate NB04e E1 Tier-A list. 45 / 51 (88.2 %) are sign-concordant (both CD↑). Top replicators include M. gnavus (HMP2 effect +1.08, FDR 3e-13), E. asparagiformis (+0.89, FDR 1e-21), H. symbiosa (+1.18, FDR 6e-22), E. innocuum (+0.28, FDR 4e-16), E. bolteae (+1.27, FDR 2e-18), E. clostridioformis (+1.04, FDR 2e-21). Only 2 of the top 20 fail: S. thermophilus (sign-discordant — HMP2 E1 effect slightly negative; potentially reflects differential dairy exposure in HMP2 vs HallAB/NielsenHB cohorts) and Bacteroides stercoris (sign-discordant, n.s.).

Synthesis. The three tests collectively upgrade Pillar 2 from "rigor-controlled on a single cohort with marginally-stable ecotype framework" to "rigor-controlled on cMD + externally replicated on HMP2 with honest documentation of cross-study ecotype variance." The operational Tier-A for NB05 is validated; the ecotype-framework-reproducibility caveat is honestly stated but bounded (the framework is cross-study variable but externally usable because projected ecotypes stratify disease and Tier-A replicates at 88 %).

(Notebooks: NB04f_loso_ecotype_stability.ipynb, NB04g_pathway_ecotype_refit.ipynb, NB04h_hmp2_external_replication.ipynb. HMP2 pull via pull_hmp2_metaphlan3.R against curatedMetagenomicData v3.18.)

5g. NB05 Tier-A scoring — prioritized target list for Pillar 4/5

Four criteria (A3–A6 from RESEARCH_PLAN.md §Criteria) applied to the 71 unique rigor-controlled candidates (51 E1 + 40 E3 provisional + 5 cross-ecotype engraftment, dedup'd). Scoring:

• A3 Literature + cohort CD-association (0–5): five independent signals per candidate — NB04c confound-free meta, HMP2 external replication concordance, ref_cd_vs_hc_differential (Kumbhari reference; log₂FC > 0.5 + FDR < 0.10), ref_species_ibd_associations (UHGG-indexed dxIBD mixed-effects), ref_phage_biology (curated top-tier targets). Distribution: 13 candidates score 0; 22 score 1; 25 score 2; 5 score 3; 5 score 4; 1 scores 5.
• A4 Protective-analog exclusion (0 / 1): fails if within-IBD-substudy CD-vs-nonIBD effect (NB04c §3) is negative — protective-analog risk — or if candidate is on the curated-protective-species list (F. prausnitzii, A. muciniphila, R. intestinalis, R. hominis, L. eligens, A. rectalis, C. scindens, C. eutactus, B. adolescentis, B. longum). Three candidates fail: Anaerostipes hadrus (−0.32 confound-free effect), Clostridium scindens (curated protective list), Roseburia faecis (−2.74 effect).
• A5 Engraftment / strain adaptation (0 / 0.5 / 1): 1.0 for the 5 donor-2708-engraftment pathobionts (M. gnavus, E. lenta, E. coli, E. bolteae, H. hathewayi); 0.5 for Kumbhari strain-competition disease-dominance or IBD-adapted-strain gene signal. 3 candidates hit the 0.5 tier (A. hadrus, B. cellulosilyticus, F. plautii).
• A6 BGC inflammatory mediator (0 / 0.5 / 1): BGC count per candidate from ref_bgc_catalog (synonymy-inverted matching so pre-GTDB names like "Ruminococcus gnavus" match canonical "Mediterraneibacter gnavus"). 1.0 if ≥ 1 BGC contains CD-enriched CB-ORFs (effect > 0.5 + FDR < 0.05 per ref_cborf_enrichment). 14 candidates score 1.0, including M. gnavus (39 BGCs, 26 CD-enriched CB-ORFs — the largest count), E. coli (93 BGCs with MIBiG matches to Colibactin / Yersiniabactin / Enterobactin), S. salivarius (98 BGCs with Salivaricin 9 / A / Cochonodin I), Streptococcus parasanguinis (51 BGCs, 11 CD-enriched), Hungatella hathewayi (6 BGCs, 4 CD-enriched).

Total Tier-A score = A3/5 + A4 + A5 + A6 (range 0–4); actionable threshold = 2.5.

6 actionable candidates (top 6 of 71 scored):

Tier-B candidates (score 2.2–2.4, sub-threshold): Enterocloster asparagiformis, Streptococcus salivarius, E. citroniae, E. clostridioformis, Blautia coccoides, Veillonella atypica, S. parasanguinis, Actinomyces oris, V. parvula. These have BGC + A4-pass + A3 = 1–2 signals but lack direct engraftment or strain-adaptation evidence; Pillar 4 phage-targetability scoring may promote any of these based on B-tier phage-availability evidence.

(Notebooks: NB05_tier_a_scoring.ipynb + run_nb05.py. The scored TSV data/nb05_tier_a_scored.tsv is the authoritative hand-off to NB06 co-occurrence networks and Pillar 4 phage-target scoring. Note: this notebook was executed via run_nb05.py rather than nbconvert due to an environment-specific numpy.bool serialization issue in the nbconvert notebook-save path; outputs are authoritative and pre-populated in the committed .ipynb.)

5h. NB06 co-occurrence networks per ecotype (H2d test)

Four per-subnet correlation networks built via CLR transform + rank-based Pearson (= Spearman rho), per-edge BH-FDR, thresholded at |rho| > 0.3 AND FDR < 0.05, with Louvain community detection (networkx.community.louvain_communities, edge-weighted by |rho|). Networkx 3.5's built-in Louvain was sufficient; FastSpar / SpiecEasi installation was held back as unnecessary for the H2d question given clear module structure at CLR-Spearman.

H2d verdict — nominally PARTIAL, biologically SUPPORTED for the pathobiont module:

The raw mean-actionable-per-module is 1.38 on E1_all + E3_all (below the ≥ 2 bar stated in the plan), but the signal is not uniformly distributed across modules. In every subnet, a single module contains 4-5 of the 6 actionable Tier-A candidates:

The remaining modules per subnet are commensal / Prevotella / diverse-healthy communities that naturally contain 0 Tier-A hits by construction. The mean-per-module statistic is diluted by these biologically-irrelevant-to-the-question modules.

Biological interpretation: Tier-A pathobionts form a single ecologically-linked co-occurrence module within CD ecotypes. Multi-target phage cocktails are therefore appropriate for the pathobiont-module members — they co-favour similar conditions (likely bile-acid dysregulation + low-oxygen inflammation) and the ecological coupling suggests a cocktail hitting 3+ of {M. gnavus, E. lenta, E. bolteae, H. hathewayi, E. coli} will have compounding effects.

Ecotype-specific module membership:

• F. plautii is in the main pathobiont module in E1 but in the generalist module in E3. Relevant for Pillar 5 per-patient cocktails: for E1 patients, F. plautii + main-pathobiont co-targeting is ecologically coherent; for E3 patients, F. plautii is less linked and may need a separate phage.
• E. coli is in the pathobiont module in E3 only, not E1. Consistent with AIEC being more characteristic of severe-Bacteroides-expanded E3 than transitional E1.

Module-anchor commensals (top-degree non-Tier-A hubs in the pathobiont modules, useful for Pillar 3 functional-driver anchoring):
- E1_all module 1: Firmicutes bacterium CAG 110, Collinsella massiliensis, Phascolarctobacterium sp CAG 266
- E3_all module 1: Butyricicoccus pullicaecorum, Anaerostipes caccae, Lactococcus lactis

Literature grounding — butyrate producers anchor the pathobiont module despite their anti-inflammatory biology. Butyricicoccus pullicaecorum is extensively studied as a butyrate-producing Clostridial-cluster-IV IBD-probiotic candidate (Geirnaert 2015a, Steppe 2014, Jeraldo 2016), with published safety data and anti-inflammatory short-chain-fatty-acid profile. Anaerostipes caccae is another canonical butyrate producer. Both being top-degree hubs in the E3 pathobiont module — not in a separate healthy-commensal module — is a biologically interesting finding: the ecological niche the pathobionts occupy is shared with butyrate-producing commensals that are CD-depleted in most pooled analyses but co-vary with pathobionts under within-ecotype co-occurrence. This is consistent with a metabolic-partner / cross-feeding interpretation (pathobiont-produced substrates support the butyrate-producing commensal; the commensal's butyrate doesn't suppress the pathobiont in this context) and suggests Pillar 3 should look specifically at cross-feeding metabolite exchange in this module. It also cautions against "preserve butyrate-producers" as a naive phage-targeting goal — these species may actually track with the pathobionts, not against them, in the CD ecological context.

(Notebook: NB06_cooccurrence_networks.ipynb; executed via run_nb06.py with pre-populated outputs in the committed .ipynb — same workaround as NB05 for the nbconvert numpy.bool issue.)

Pillar 2 close-out

With NB06 complete, Pillar 2 is fully closed: rigor-controlled Tier-A (NB04b-e) → externally replicated on HMP2 (NB04h) → scored + prioritized (NB05, 6 actionable of 71) → co-occurrence structure mapped (NB06, single-pathobiont-module finding). The set of scored + module-assigned + hub-ranked Tier-A is the complete input package for Pillar 4 (phage-availability × target) and Pillar 5 (UC Davis per-patient cocktail drafts).

6. Taxonomy synonymy layer is the project's reusable foundation

data/species_synonymy.tsv — 2,417 alias → 1,848 canonical species, grounded in ref_taxonomy_crosswalk NCBI taxid matching with GTDB r214+ genus renames supplemented. This was motivated by a failure mode discovered in NB00: fact_taxon_abundance contains three divergent taxon-name formats between cohorts (CMD_IBD short names, CMD_HEALTHY full MetaPhlAn3 lineage, KUEHL_WGS Kaiju), and a naive pivot splits the same species into multiple zero-overlap rows, producing log₂FC ≈ 28 artifacts.

The synonymy layer was built once in NB01b and is joined against by NB00, NB01, NB02, NB03, and every downstream notebook. Documented as a project-level pitfall (docs/pitfalls.md) and a candidate BERIL convention — any project integrating multi-cohort microbiome data needs this layer, and the pattern (NCBI taxid + GTDB-version-aware rename table) generalizes.

7. NB07a — Pathway DA + H3a v1.7 three-clause falsifiability (Pillar 3 opener)

First Pillar 3 notebook, executed under RESEARCH_PLAN.md v1.7 norms (post-adversarial-review). Per norm N12, primary contrast is within-IBD-substudy CD-vs-nonIBD meta on fact_pathway_abundance (HUMAnN3 MetaCyc, CMD_IBD only). Per norm N15, substudy meta-viability re-verified for the pathway modality: 3 robust (HallAB_2017, IjazUZ_2017, NielsenHB_2014) + 1 boundary (LiJ_2014, nonIBD = 10) — not the "4 meta-viable" framing v1.6 had inherited from NB04e's taxonomic-modality counts. VilaAV_2018 excluded (CD = 216, nonIBD = 0). 575 unstratified MetaCyc pathways → 409 after 10%-prevalence filter.

H3a v1.7 verdict: PARTIALLY SUPPORTED — 2 of 3 clauses pass; clause (b) is structurally degenerate, not a fundamental refutation.

Top pathway-pathobiont attribution recapitulates known AIEC biology. The top 25 pairs are all Escherichia coli pathways with biological coherence:

Tier-A pair-count distribution: E. coli 105 pairs > 0.4 (76% of total signal); H. hathewayi 16; M. gnavus 8; F. plautii 7; E. lenta 1; E. bolteae 0. E. coli's domination reflects three factors: highly specialized AIEC functional repertoire well-matched to MetaCyc, higher relative-abundance variance, and E. coli genome content being especially well-represented in HUMAnN3 (vs more obligate-anaerobe Tier-A members like E. lenta / E. bolteae). E. bolteae's zero pathway-attribution despite being NB05-actionable signals that its CD-up signal may be at the BGC-level (testable via NB08, executed) or Kumbhari strain-frequency level (NB10a fact_strain_competition) rather than HUMAnN3-pathway level.

Honest interpretation of clause (b) failure. The test was structurally underpowered: with only 3 of 52 CD-up pathways landing in the 7 a-priori categories, the top-3-concentration statistic is trivially 100% under both observed and random-allocation null. The cMD pathway DA at the unstratified MetaCyc level does not preferentially load on the classical IBD-themed categories (bile-acid 7α-dehydroxylation, mucin degradation, sulfidogenesis, TMA/TMAO, eut/pdu, polyamine, AA-decarb). Instead, the CD signal captures broader bacterial-fitness-in-inflamed-gut themes (heme/iron, glyoxylate, fat metabolism, allantoin/purines) that are not in the prior-literature category set as constructed. Two non-mutually-exclusive interpretations:

1. Category-set choice. The 7 a-priori categories were drawn from prior pathobiont-mechanism literature; broadening to include "alternative-electron-acceptor metabolism," "iron acquisition," "fat-utilization-in-inflammation" categories may recover the (b) signal.
2. Stratified-pathway resolution. The unstratified-pathway level may be dominated by housekeeping pathways shared across many genera. Per-species stratification (NB07b) — does M. gnavus gain bile-acid-deconjugation pathways CD-vs-nonIBD? — is the appropriate next test of H3a (b).

NB07b is therefore well-positioned: stratified-pathway (PWY-XXX|g__species) gives ~ 42K pathway-species combinations to query directly. Combined with the (c) attribution signal already strong here, H3a (b) is best resolved as a follow-up question, not as a current refutation.

Output artifacts:
- data/nb07a_pathway_meta.tsv — 409 pathways × meta statistics
- data/nb07a_pathway_pathobiont_pairs.tsv — 2,454 pathway × Tier-A-core species pairs with per-substudy + meta ρ
- data/nb07a_h3a_verdict.json — formal H3a v1.7 verdict with all permutation-null statistics

(Notebook: NB07a_pathway_DA_H3a_falsifiability.ipynb; executed via run_nb07a.py per the established nbconvert numpy.bool workaround.)

8. NB07b — Stratified-pathway DA per Tier-A-core species (H3a (b) species-resolved re-test)

NB07a clause-(b) failed because only 3 of 52 CD-up unstratified MetaCyc pathways landed in the 7 a-priori IBD categories — test was structurally underpowered. NB07b tests the alternative: at the species-resolved (HUMAnN3 stratified-pathway form PWY-XXX|g__species), do per-species CD-up pathways concentrate in those categories?

H3a (b) verdict at species-resolved level: NOT SUPPORTED for all 6 Tier-A core species — but the verdict is structural (the 7 a-priori categories are sparse in per-species pathway repertoires), not biological.

Combined NB07a + NB07b H3a (b) conclusion: the 7 a-priori IBD-mechanism categories (bile-acid, mucin, sulfide, TMAO, eut/pdu, polyamine, AA-decarb) are too narrow to capture HUMAnN3 MetaCyc CD signal at either unstratified or species-resolved level. CD-associated pathways are dominated by biosynthesis / niche-shift signals not the prior-literature themes. H3a (b) refutation is real for the v1.7-stated category set but does NOT mean "no compositional themes exist" — it means a different category schema (broader, programmatically derived from MetaCyc taxonomy) is needed.

But NB07b reveals biology NB07a couldn't see at the per-species level.

Hungatella hathewayi (33 passing pathways) — coherent within-carrier biosynthesis-vs-degradation shift CD-vs-nonIBD:

Interpretation: under CD, H. hathewayi shifts to a biosynthetic / catabolic-stress state (pentose phosphate, anabolic AA biosynthesis, purine salvage) and away from sugar-utilization (lactose, galactose, propanediol, inositol). The propanediol-degradation CD-down within H. hathewayi is consistent with niche partitioning — E. coli may dominate this niche under CD (consistent with NB07a's E. coli propanediol-biosynthesis CD-up at cohort level).

Escherichia coli (20 passing) shows a within-carrier CD-DOWN per-pathway abundance — opposite of the cohort-level CD-up direction:

This is not contradictory to NB07a's cohort-level E. coli CD-up:
- At cohort level: E. coli relative abundance is CD-up, so total pathway flux scales up
- Within carriers: per-cell pathway repertoire is CD-down — E. coli in CD samples show less metabolic versatility per cell

Two non-mutually-exclusive interpretations:
1. CD's E. coli are an AIEC-specialized subset that has shed peripheral metabolic capabilities (consistent with NB05 §5g Yersiniabactin / Colibactin / Dubinsky 2022 IBD-specific lineage finding)
2. CD's E. coli face metabolic competition from co-abundant Klebsiella / other Enterobacteriaceae, so per-cell read-mapping share is lower across pathways

Both align with the AIEC narrative. Deeper strain-level analysis on cMD raw reads is dropped per plan v1.9 (no-raw-reads scope); the no-raw-reads alternative is genome-content survey on kbase_genomes + Kumbhari fact_strain_competition (Future Direction #9).

Other 4 Tier-A core species (M. gnavus, E. lenta, F. plautii, E. bolteae) show small within-carrier shifts (≤ 4 passing pathways each). Their CD signal is dominated by carriage prevalence, not within-carrier metabolic shift. Implication: pathway-level mechanism is not the right resolution for these species; the no-raw-reads alternatives are BGC-level (NB08, executed) and Kumbhari strain-frequency (NB10a) for the species in fact_strain_competition.

Output artifacts:
- data/nb07b_stratified_pathway_da.tsv — 621 species-pathway × meta rows
- data/nb07b_h3a_b_species_verdict.json — per-species formal verdict
- figures/NB07b_stratified_H3a_b.png — 4-panel diagnostic

(Notebook: NB07b_stratified_pathway_DA.ipynb; executed via run_nb07b.py.)

9. NB07 v1.8 H3a (b) retest — MetaCyc class hierarchy + Fisher per-theme enrichment

The v1.7 H3a (b) FAIL verdict (NB07a §6-8 + NB07b §4-5) was driven by regex-on-pathway-names limitations: only 44/409 background pathways matched the 7 a-priori category patterns; pathways like PWY-5920 (superpathway of heme biosynthesis from glycine) were silently categorized as "0_other" because "heme" / "iron" were not in the regex set. v1.8 replaces the regex approach with structured MetaCyc class assignments from /global_share/KBaseUtilities/ModelSEEDDatabase/Biochemistry/Aliases/Provenance/MetaCyc_Pathways.tbl and expands to 12 IBD-relevant themes including iron/heme acquisition, fat metabolism / glyoxylate, anaerobic respiration, purine/pyrimidine recycling, and aromatic AA / chorismate / indole.

Per plan norm N17 (added in v1.8): prefer ontology / class hierarchy over name-pattern regex for pathway / gene / metabolite categorization wherever feasible. Regex remains a sensitivity check.

v1.8 cohort-level H3a (b) verdict: SUPPORTED. Per-theme Fisher's exact enrichment (CD-up × in-theme) on the NB07a passing pathways across 12 IBD themes; BH-FDR across themes; theme supported if FDR < 0.10 AND odds ratio > 1.5:

Iron/heme acquisition is the dominant CD-up biochemical theme (8.1× over background; 29 % of CD-up pathways vs 8 % expected). The 15 CD-up pathways in this theme include heme biosynthesis (PWY-5920, ρ = 0.64 with E. coli in NB07a §c), heme degradation, and siderophore-related pathways.

This completely reverses the v1.7 H3a (b) "FAIL" verdict — driven entirely by category-schema choice, not by any change in the underlying data. The v1.7 regex on pathway descriptive names matched "iron" / "heme" only via "PWY-5920: superpathway of heme biosynthesis from glycine" → "0_other"; the v1.8 ModelSEED class hierarchy correctly assigns PWY-5920 to HEME-SYN, Heme-b-Biosynthesis, Cofactor-Biosynthesis, Tetrapyrrole-Biosynthesis.

v1.8 species-level H3a (b): H. hathewayi has two themes supported (16 CD-up pathways):
- 11_purine_pyrimidine_recycling: OR=4.86, FDR=0.048 (7/16 CD-up pathways)
- 04_TMA_choline: OR=9.33, FDR=0.048 (4/16 CD-up pathways)

Other Tier-A core species have insufficient species-level CD-up pathway counts for per-theme power (M. gnavus, E. lenta, F. plautii have 0 CD-up; E. coli, E. bolteae have 1-2).

Four-way convergence on iron biology as the dominant CD pathobiont specialization in this dataset:

1. NB05 §5g — E. coli MIBiG matches: Yersiniabactin + Enterobactin (both iron siderophores) + Colibactin (same pks pathogenicity island)
2. NB07a §c — top pathway-pathobiont attribution: heme biosynthesis ↔ E. coli (ρ = 0.640)
3. NB07 v1.8 H3a (b) — iron/heme is the dominant CD-up theme (OR = 8.11, FDR = 7e-6)
4. AIEC literature — Dalmasso 2021 (yersiniabactin), Prudent 2021 (LF82 IBC formation via yersiniabactin), Dogan 2014 (AIEC iron-pathway enrichment) all flag iron acquisition as central AIEC fitness mechanism

This is a robust, multi-line-of-evidence-supported claim about CD pathobiont biology that emerges only after the v1.8 schema fix. The v1.7 "FAIL — degenerate" verdict was masking a 8.1-fold iron-theme enrichment.

v1.8 also reveals H. hathewayi-specific themes: purine/pyrimidine recycling and TMA/choline metabolism. H. hathewayi's NB07b CD-up pathways included pentose phosphate (precursor for nucleotide biosynthesis), purine nucleobases degradation, and PWY-6803 phosphatidylcholine acyl editing — exactly the building-block pathways of the two enriched themes. H. hathewayi is a known TMA producer (CutC/D activity) and the choline-metabolism enrichment is mechanistically coherent.

Methodological lesson (added to docs/discoveries.md): v1.7 → v1.8 is a major scientific reversal driven entirely by category-schema choice. v1.7 "no compositional themes" (FAIL) → v1.8 "iron/heme is the dominant theme" (SUPPORTED, OR 8.1, FDR 7e-6). Same data, same DA. The lesson: regex-on-pathway-names is a poor substitute for curator-validated ontology / class hierarchy when one is available; ModelSEEDDatabase ships a usable MetaCyc class hierarchy with 90 % coverage of HUMAnN3 outputs and should be the default for any pathway category-enrichment test in BERIL projects.

Output artifacts:
- data/nb07_h3a_v18_pathway_classes.tsv — pathway × MetaCyc-classes × IBD-themes (audit trail)
- data/nb07_h3a_v18_cohort_enrichment.tsv — per-theme Fisher enrichment cohort-level
- data/nb07_h3a_v18_species_enrichment.tsv — per-species per-theme Fisher enrichment
- data/nb07_h3a_v18_verdict.json — formal v1.8 H3a (b) verdict
- figures/NB07_H3a_v18_class_enrichment.png — visualization

(Script: run_nb07_h3a_v18.py. Builds on NB07a + NB07b pathway-DA outputs; adds class-based theme enrichment per plan v1.8.)

10. NB07c — Module-anchor commensal × pathobiont metabolic coupling (H3a-new)

NB06 surfaced the H2d co-occurrence finding that 4–5 of 6 actionable Tier-A pathobionts co-cluster in a single CD-specific module per ecotype. The follow-up question — and the H3a-new test specified in plan v1.7 X4 fix — is whether the non-pathobiont anchors of these CD-specific modules show metabolic coupling with the pathobionts they sit beside. Two CD-specific modules from data/nb06_module_hubs.tsv provide the substrate:

• E1_CD module 0 (75 nodes): anchor commensals Clostridiales bacterium 1_7_47FAA, Anaerostipes caccae (the only genuine butyrate-producer in module-anchor commensals), Bacteroides nordii; with the 5 actionable Tier-A pathobionts of the E1_CD subnet (H. hathewayi, F. plautii, E. bolteae, E. lenta, M. gnavus).
• E3_CD module 1 (57 nodes): anchor commensals Actinomyces sp. oral-taxon-181, Actinomyces sp. HMSC035G02 (both oral cavity ectopic colonizers), Lactonifactor longoviformis (lactate utilizer); with 4 module pathobionts (E. lenta, H. hathewayi, E. coli, M. gnavus).

Per (anchor, pathobiont) pair: within-IBD-substudy Spearman ρ across CMD_IBD samples; Fisher z-meta across the 3 robust substudies (ZellerG_2014, NielsenHB_2014, IjazUZ_2017); sign concordance. Plus an iron-context layer: triple correlation (anchor × pathobiont × iron-pathway) over the 15 v1.8 iron/heme pathways, to test whether the v1.8 iron-theme is a community-wide signature or pathobiont-specific.

E1_CD coupling — A. caccae shows clean strong-positive coupling with all 5 module pathobionts (sign concordance 1.0 across all 3 substudies):

The A. caccae × pathobiont pattern is consistent with butyrate-producer cross-feeding embedded in the CD pathobiont module — pathobiont-released substrates (M. gnavus glucorhamnan / mucin sugars; F. plautii bile-acid metabolites; H. hathewayi lactate) feed A. caccae's butyrogenic fermentation. The §2 iron-pathway layer rules against an iron-cross-feeding mechanism: ρ(A. caccae × iron-pwy) = +0.13 (mean over 15 iron pathways), much weaker than ρ(E. coli × iron-pwy) = +0.45. So the A. caccae coupling is not iron-mediated — most likely substrate / sugar / lactate mediated. Cross-feeding vs shared-environment disambiguation is deferred to NB09c (metabolite corroboration).

The B. nordii negative coupling with M. gnavus / F. plautii is consistent with niche competition for similar polysaccharide substrates — B. nordii is a generalist Bacteroides that competes with the same mucin / glycan-degrading pathobionts, and CD selects one over the other.

E3_CD coupling — only Lactonifactor × E. lenta is strong (ρ_meta = +0.27). The two oral Actinomyces anchors couple weakly (~ρ=0.17–0.19 with M. gnavus / E. lenta). Oral Actinomyces in the gut are co-trafficked ectopic colonizers under CD inflammation rather than metabolic partners; their NB06 module membership reflects shared inflammation-driven colonization, not metabolic coupling.

Iron-pathway co-variation concentrates on E. coli — narrowing the v1.8 iron-theme interpretation:

The 15 iron pathways include ENTBACSYN-PWY (Enterobactin biosynthesis, E. coli-canonical), HEMESYN2-PWY (heme biosynthesis II), and 8 menaquinol-biosynthesis pathways. E. coli's ρ = +0.45 across 2,674 CMD_IBD samples means iron-pathway abundance scales proportionally with E. coli abundance — i.e., the iron-pathway signal is dominantly carried by E. coli. Other Tier-A core species show weak iron-pathway coupling, consistent with them having other CD specialization mechanisms (TMA / choline for H. hathewayi per v1.8 §9; bile-acid 7α-dehydroxylation for F. plautii; glucorhamnan / mucin for M. gnavus).

This narrows the v1.8 iron-theme interpretation: rather than "all CD pathobionts have iron specialization," the more accurate framing is "CD's E. coli (AIEC subset) drives the iron-acquisition theme; other Tier-A pathobionts have non-iron specializations." This is mechanistically coherent with NB05 §5g (only E. coli of the actionable Tier-A had iron-siderophore MIBiG matches: Yersiniabactin + Enterobactin) and the NB07b within-carrier E. coli CD-DOWN per-pathway pattern (AIEC strain-level specialization at the cost of generalist metabolic capabilities).

H3a-new verdict: PARTIALLY SUPPORTED. A. caccae × pathobiont coupling is clean in E1_CD (4 strong-positive pairs, all sign-concordant); E3_CD anchors lack metabolic-coupling structure (oral-gut co-trafficking dominates). Cross-feeding vs shared-environment disambiguation deferred to NB09c (metabolite-level test).

Pillar 4 cocktail-design implication — metabolic-coupling cost. The A. caccae × pathobiont coupling has a direct clinical implication for cocktail design: a phage cocktail that depletes M. gnavus / F. plautii / H. hathewayi / E. bolteae may incidentally reduce A. caccae abundance through loss of substrate. A. caccae is the only genuine butyrate-producer in the module-anchor commensals, and butyrate is anti-inflammatory; its incidental depletion could partially offset the cocktail's therapeutic benefit. NB05 actionable Tier-A targets need a "metabolic-coupling cost" annotation before cocktail finalization — for each candidate target, what beneficial commensals depend on it as a substrate source, and what is the predicted net inflammatory balance? This is the species-pair-level extension of the H2d concern surfaced in NB06 (single-pathobiont-module → cocktail-design implication).

Output artifacts:
- data/nb07c_anchor_pathobiont_species_rho.tsv — 27 (anchor × pathobiont) pairs with ρ_meta, per-substudy values, sign concordance
- data/nb07c_anchor_pathobiont_iron_triple.tsv — 405 (anchor × pathobiont × iron-pathway) triples
- data/nb07c_h3a_new_verdict.json — formal verdict
- figures/NB07c_anchor_pathobiont_coupling.png — 2-panel heatmap

(Script: run_nb07c.py. Builds on NB06 module hubs + NB07_v1.8 iron-pathway list.)

11. NB08a — BGC × pathobiont enrichment (H3c) — genomic mechanism layer

NB05 §5g qualitatively flagged E. coli MIBiG matches Yersiniabactin + Enterobactin + Colibactin. NB07 v1.8 §9 found iron/heme acquisition is the dominant CD-up MetaCyc pathway-class theme (OR=8.1, FDR 7e-6). NB07c §10 found iron-pathway co-variation concentrates on E. coli at the sample-correlation level. NB08a is the genomic-content level test: do Tier-A pathobiont genomes (per ref_bgc_catalog, Elmassry 2025; 10,060 BGCs across 6,221 species-annotated entries) carry an over-represented iron-siderophore / genotoxin biosynthetic gene-cluster signature, and is this signature uniformly distributed across actionable Tier-A or concentrated on E. coli?

Test 1 — BGC-theme enrichment (Fisher's exact, Tier-A core BGCs in theme vs background BGCs in theme; BH-FDR across themes; 4 IBD-relevant themes covering iron-siderophore MIBiG matches, genotoxin/microcin MIBiG matches, RiPP-bacteriocin classes, and NRPS-PKS-hybrid classes):

The iron_siderophore Fisher OR of 44.4 is one of the largest enrichments in the project. It complements v1.8 §9's pathway-class iron OR of 8.1 (different evidence stream — pathway-cohort co-occurrence vs genomic gene-cluster content) and converges on the same biology: iron acquisition is a CD pathobiont-defining genomic capability, not just a pathway-level cohort signal.

Test 2 — Per-Tier-A-core species iron + genotoxin MIBiG breakdown:

E. coli alone carries the iron+genotoxin BGC signature within actionable Tier-A core. Its 54 iron BGCs comprise 19 Yersiniabactin, 16 Enterobactin, plus 19 BGCs of class=siderophore (some redundant per-strain assemblies); the 25 genotoxin BGCs comprise 8 Colibactin, 15 Microcin B17, and 2 Microcin J25. The other 5 Tier-A core species sit in MIBiG dark matter — they carry substantial BGC content (E. lenta 41, M. gnavus 58, E. bolteae 18) but no MIBiG-annotated iron or genotoxin clusters. This is the genomic-level confirmation of the v1.8 + NB07c narrowing: iron biology in CD pathobionts is essentially an E. coli (AIEC) phenomenon. Other Tier-A pathobionts have non-iron CD-association mechanisms (consistent with v1.8 §9: H. hathewayi purine + TMA/choline themes).

Test 3 — CB-ORF CD-vs-HC enrichment per Tier-A core (read-level; from ref_cborf_enrichment 5,157 CB-ORFs):

5 of 6 Tier-A core species have CD-up CB-ORF rates above the 2.5 % catalog background (range 14–82 %); E. lenta is the exception, with CB-ORFs CD-DOWN at 54 %. This complements the BGC-MIBiG analysis: even where MIBiG annotations are absent (the 5 dark-matter Tier-A core), per-CB-ORF CD-vs-HC RPKM is independently elevated in CD samples. The E. lenta CD-DOWN pattern is consistent with NB07b's species-resolved finding that E. lenta per-pathway abundance is mostly carriage-prevalence-driven (not within-carrier abundance-shifted) and aligns with the canonical Eggerthella CD-association mechanism being drug-metabolism (cardiac glycoside inactivation, Koppel et al. 2018) rather than BGC-encoded inflammatory mediators.

Test 4 — ebf/ecf cohort meta CD-vs-HC (per ref_ebf_ecf_prevalence, 1,349 samples × 4 cohorts; Mann-Whitney CD-vs-HC per cohort + Stouffer's z-meta):

All 4 cohorts (HMP2-IBDMDB, MetaHIT, LLDEEP-NLIBD, PRISM) show CD > HC for both ebf and ecf, with cliff-deltas 0.17–0.73. The Elmassry 2025 immunoactive fatty acid amide BGC family CD-up finding replicates cleanly in our cohort-meta design at p < 1e-31 — the largest single effect in the project so far.

H3c verdict — PARTIALLY SUPPORTED. The hypothesis ("BGC-encoded inflammatory mediators localize to a minority of Tier-A pathobionts and show CD-enrichment beyond what species-level abundance captures"):

1. Localization to minority: confirmed — iron+genotoxin BGCs are uniquely E. coli within actionable Tier-A core (1 of 6).
2. Beyond species-abundance: partially confirmed — CB-ORF CD-vs-HC is independently CD-up enriched (read-level, not species-derived) for 5 of 6 Tier-A core; ebf/ecf is independently CD-up at p<1e-31 across 4 cohorts. However, the strict species × BGC interaction-term test specified in the original H3c is untested (would require species-stratified per-sample BGC abundance, not in the current pre-computed mart slice).

Five-line iron-acquisition convergence narrative: NB05 §5g (per-actionable MIBiG lookup) → NB07a §c (pathway × pathobiont attribution, ρ=0.640 heme↔E.coli) → NB07 v1.8 §9 (cohort pathway-class enrichment, OR=8.1) → NB07c §2 (sample-level co-variation, ρ=0.45 E.coli×iron-pwy) → NB08a §2 (genomic BGC content, OR=44.4 driven by E.coli's 54 iron BGCs) — five independent evidence streams converging on AIEC iron-acquisition as a central CD-pathobiont specialization mechanism.

Pillar 4 cocktail-design implication (sharpened): phage cocktail design should distinguish:
- E. coli component: target AIEC subset specifically — Yersiniabactin/Enterobactin/Colibactin-positive strains. Per plan v1.9, AIEC strain-resolution from cMD raw reads is dropped; the no-raw-reads alternative is kbase_genomes pks + iron-BGC genome-content query (Future Direction #9);
- M. gnavus, E. lenta, E. bolteae, H. hathewayi, F. plautii: BGC-mechanism dark — design relies on Tier-A scoring + NB06 module membership + NB07b within-carrier metabolic signature, not on BGC presence;
- ebf/ecf RPKM as a sample-level CD biomarker (no per-species attribution available) — usable for treatment-response monitoring, not for cocktail target selection.

Output artifacts:
- data/nb08a_tier_a_bgc_repertoire.tsv — per-species BGC repertoire summary
- data/nb08a_bgc_theme_enrichment.tsv — 4 themes × Fisher OR + FDR (3 supported)
- data/nb08a_tier_a_iron_genotoxin_per_species.tsv — per-species iron + genotoxin MIBiG breakdown
- data/nb08a_cborf_enrichment_per_tier_a.tsv — per-species CB-ORF CD-up rate vs background
- data/nb08a_ebf_ecf_cd_vs_hc.tsv — ebf/ecf cohort meta z-stats
- data/nb08a_h3c_verdict.json — formal H3c verdict
- figures/NB08a_bgc_pathobiont_enrichment.png — 3-panel summary

(Script: run_nb08a.py. Builds on NB05 actionable Tier-A; ref_bgc_catalog + ref_cborf_enrichment + ref_ebf_ecf_prevalence per RESEARCH_PLAN.md NB08a spec.)

12. NB09a — HMP2 metabolomics CD-vs-nonIBD (H3d-DA)

NB09a is the first metabolomics analysis in the project. It tests whether the metabolite axes that distinguish CD from nonIBD in HMP2 are coherent with the iron / TMA / fat / bile-acid / fatty-acid-amide themes already established by the NB07–NB08a pathway and BGC analyses. Per plan v1.9 (raw-read deprecation), no raw-FASTQ reprocessing — uses precomputed HMP2 metabolomics from the mart + cMD R-package HMP2 sample metadata for diagnosis labeling. 468 of 546 metabolomics samples (86 %) match cMD HMP2 metagenomics samples directly by sample-ID code (CSM*); the matched set covers 106 subjects (50 CD + 30 UC + 26 nonIBD) with diagnosis. Subject-level analysis uses one sample per subject (first-occurrence visit) to avoid within-subject correlation; UC excluded from primary contrast.

Per-metabolite Mann-Whitney + cliff_delta + BH-FDR on 592 named (HMDB-annotated) HMP2 metabolites (out of 81,867 total, mostly unnamed peaks). 579 testable (≥5 samples in each group). Passing threshold: FDR < 0.10 + |cliff_delta| > 0.20 → 52 metabolites pass (50 CD-up, 2 CD-down).

Per-theme Fisher's exact across 11 IBD-relevant chemical-class themes (BH-FDR across themes; supported = FDR < 0.10 AND OR > 1.5):

Polyamines CD-up (5 of 9): putrescine (cliff=+0.45, FDR=0.08), N1-acetylspermine (+0.56, 0.06), N-acetylputrescine (+0.46, 0.06), anserine (+0.43, 0.08), diacetylspermine (+0.43, 0.08). Established IBD biomarker class (reviewed Pegg 2014; reported as CD biomarker by Wang 2018, Franzosa 2019). Mechanistically interesting: the v1.8 §9 pathway-level 06_polyamine_urea theme was CD-DOWN at pathway-level (OR=0.42) but CD-UP at metabolite-level here (OR=14.6) — this is mechanistically coherent: polyamine metabolite-pool accumulation can result from increased catabolism of dietary protein / mucin + reduced microbial polyamine clearance without requiring elevated biosynthesis flux. The metabolite-pool readout is the clinically-actionable observation; the pathway-level signal reflects production capacity.

Long-chain PUFAs CD-up (6 of 15): adrenate (C22:4, +0.55, FDR=0.027), arachidonate (C20:4, +0.48, 0.05), docosapentaenoate (C22:5, +0.45, 0.06), docosahexaenoate (C22:6, +0.41, 0.08), eicosapentaenoate (C20:5, +0.40, 0.09). Covers both n-6 (adrenate, arachidonate — eicosanoid precursors) and n-3 (DHA, DPA, EPA) classes. The CD-up signal can reflect (a) impaired host fatty-acid uptake / β-oxidation in inflamed mucosa, (b) increased dietary fat mobilization, (c) reduced microbial PUFA biohydrogenation by Lactobacillus / Roseburia spp. leading to free-PUFA pool accumulation. Mechanistically connected to v1.8 §9 10_fat_metabolism_glyoxylate theme (theme not Fisher-significant at pathway level — OR=0.88 — but the metabolite-pool elevation IS theme-significant here).

Bile acids — only 3 of 21 CD-up but the 3 are tauro-conjugated (free taurine, tauro-α-muricholate, tauro-β-muricholate). Free taurine (the conjugating amino acid) is CD-up at cliff=+0.47. Tauro-α/β-muricholate is CD-up at +0.40. Consistent with reduced microbial bile-acid 7α-dehydroxylation in CD — primary tauro-conjugated BAs accumulate when F. plautii / C. scindens / Eggerthellaceae dehydroxylation activity is impaired (canonical Franzosa 2019 finding). Corroborates NB05 actionable F. plautii's mechanistic role in 7α-dehydroxylation.

Acyl-carnitines C16 + C18:1 CD-up (2 of 22): C16 carnitine (+0.41, FDR=0.08), C18:1 carnitine (+0.39, 0.09). Long-chain fatty-acid β-oxidation intermediates. Mechanistically connected to v1.8 §9 10_fat_metabolism theme + H. hathewayi TMA/choline (carnitine sits in the same metabolic neighborhood as choline → TMA → fatty-acid β-oxidation). Theme-level enrichment OR=1.06 not significant because the 22-carnitine background is dominated by short/medium-chain forms that don't shift in CD.

Lipid classes (CE / SM / TAG / DAG / ceramide) — 19 CD-up but theme-level OR=1.64 (FDR=0.76, not significant). Pattern is informative: dominated by cholesteryl-esters of long-chain PUFAs (12 CEs: C20:4, C20:5, C18:2, C18:3, C16:0, C18:1, C16:1; 4 SMs: C14:0, C16:0, C16:1, C24:1; 1 ceramide: C24:1; 2 TAGs). The CE-PUFA elevation directly mirrors the free-PUFA elevation (CEs are the storage form of esterified PUFAs). Theme-level Fisher doesn't pass because the 163-metabolite background is too dilute (most lipid classes don't shift, just the long-chain-PUFA-conjugated ones).

SCFAs: 0 of 14 CD-DA. None of acetate, propionate, butyrate, valerate, hexanoate reach FDR<0.10. Either (a) HMP2 LC-MS untargeted methods undersample SCFAs (volatile, polar — typically need GC-MS), (b) the cohort-level subject-averaged contrast masks within-subject variation, or (c) gross SCFA pool differences between CD and nonIBD are smaller than within-subject variability. This is an important null for the NB07c cross-feeding hypothesis — butyrate isn't differentially abundant CD-vs-nonIBD at subject-level in HMP2, so the A. caccae × pathobiont cross-feeding hypothesis cannot be confirmed via cohort-level butyrate DA. The NB09c sample-level paired test (does butyrate co-vary with A. caccae × pathobiont co-occurrence within paired metabolomics+metagenomics samples) remains the right test and is the natural follow-up.

Urobilin CD-DOWN (cliff=-0.38, FDR=0.09). Urobilin is the gut-bacterial catabolic product of bilirubin (produced by Clostridium / Bacteroides species expressing bilirubin reductase). CD-DOWN urobilin = reduced gut-bacterial bilirubin reduction = consistent with dysbiosis / loss of urobilinoid-producing commensals (Hall 2024; Vital 2018).

Convergence summary table — NB09a metabolomics × NB07-pillar pathway findings:

H3d-DA verdict: SUPPORTED. 52 metabolites pass DA + 2 themes pass Fisher enrichment. Both falsifiability gates met. NB09a adds two new mechanism axes to the project narrative (polyamines + long-chain PUFAs) that complement the iron + bile-acid + TMA findings from NB07–NB08a.

Pillar 4/5 implications:
- Polyamine and PUFA elevations are sample-level CD biomarkers — usable for clinical follow-up (treatment response monitoring) but not per-species cocktail-design targets (no per-species attribution available).
- Tauro-muricholate elevation confirms the mechanistic premise for F. plautii targeting (the bile-acid 7α-dehydroxylation deficit is real at the metabolite-pool level, not just predicted from pathway DA).
- The SCFA null at cohort level sharpens the NB07c "metabolic-coupling-cost" annotation — depleting A. caccae may not produce a measurable cohort-level butyrate change, but the within-sample paired test (NB09c) remains the right test for the cross-feeding causal claim.

Output artifacts:
- data/nb09a_metab_da_cd_vs_nonibd.tsv — 579 named metabolites × Mann-Whitney + cliff + FDR + theme assignment
- data/nb09a_metab_theme_enrichment.tsv — 11 themes × CD-up/CD-down × Fisher OR + FDR
- data/nb09a_h3d_da_verdict.json — formal H3d-DA verdict (SUPPORTED)
- figures/NB09a_metabolomics_cd_vs_nonibd.png — 3-panel summary

(Script: run_nb09a.py. Subject-level analysis on cMD-fetched HMP2 sample-to-subject metadata + mart fact_metabolomics; per plan v1.9 no-raw-reads constraint.)

13. NB09c — Sample-level paired metabolomics × metagenomics: NB07c cross-feeding disambiguation + bile-acid 7α-dehydroxylation network

NB07c §10 left a key disambiguation deferred: the A. caccae × pathobiont species-level coupling at +0.39 (E. bolteae), +0.33 (H. hathewayi), +0.31 (M. gnavus), +0.29 (F. plautii) was consistent with either butyrogenic cross-feeding or shared-environment co-response to the same CD niche. The proposed disambiguation: paired sample-level metabolomic-metagenomic correlation could surface candidate intermediate metabolites (lactate, mucin-glycan products, bile-acid metabolites) shared between anchor and pathobiont. Per plan v1.9 (no-raw-reads scope), NB09c executes this paired analysis using 468 paired CSM* HMP2 samples that have both metabolomics and cMD MetaPhlAn3 metagenomics.

Cross-feeding-triangle test (strict criteria: anchor × pathobiont same sign + both |ρ|>0.20 + both FDR<0.10): only 7 triangles across 8 species × 583 named HMP2 metabolites. Top candidates: caffeine (B. nordii / F. plautii), linoleoylethanolamide (B. nordii / E. lenta; A. caccae / E. lenta), urobilin (A. caccae / E. lenta; B. nordii / E. lenta), cholate (A. caccae / F. plautii; B. nordii / F. plautii). None of the 7 candidates is a butyrogenic cross-feeding intermediate; the pattern is consistent with shared-niche health-direction co-occurrence — A. caccae, B. nordii, F. plautii, E. lenta are all (mostly) commensal species that together correlate negatively with primary BA cholate and positively with urobilin (NB09a CD-DOWN), i.e. co-occur in healthy / normobiotic samples vs CD-dysbiotic samples.

Curated cross-feeding panel — direction-of-association profile:

Three findings emerge:

1. Butyrate × A. caccae = +0.10 is significant but very weak. Consistent with A. caccae as a butyrate producer at the population level, but the cohort-level signal is weak — partly methodological (HMP2 LC-MS undersamples SCFAs; only butyrate and propionate detected) and partly because gut butyrate pool is dominated by other butyrate producers (Faecalibacterium, Roseburia*).

2. Lactate × A. caccae = +0.18* vs lactate × F. plautii = −0.23* and × E. bolteae = −0.20* — opposite signs. If lactate were a cross-feeding intermediate (pathobiont produces, A. caccae consumes), we would expect same-sign positive (both increase together) or asymmetric same-sign (production-consumption coupling). The opposite-sign pattern is the strongest evidence against the cross-feeding hypothesis — A. caccae and the metabolic-coupling-candidate pathobionts occupy different metabolic niches at the cohort level.

3. E. coli dominates the cohort-level correlation signal: choline +0.25, carnitine +0.19, tryptophan +0.25, putrescine +0.19, cadaverine +0.45 (the strongest correlation in the entire panel). These match canonical E. coli / Enterobacteriaceae metabolism: lysine decarboxylase → cadaverine; arginine decarboxylase → putrescine; choline + carnitine substrates. The v1.8 §9 H. hathewayi TMA/choline finding from cMD pathway-level analysis does NOT strongly replicate at HMP2 sample level — H. hathewayi × choline ρ=−0.07 (NS); E. coli dominates. This narrows v1.8 §9 — the TMA/choline signal is a combined Enterobacteriaceae + Lachnospiraceae signal, not specifically H. hathewayi. Possible reasons: HMP2 has lower H. hathewayi prevalence (25 % of paired samples) than the cMD studies that drove the v1.8 finding; at the sample level, E. coli's high prevalence (50 %) and abundance variance dominate the cohort-correlation signal.

Bile-acid 7α-dehydroxylation network — independent strong finding:

F. plautii, E. lenta, and E. bolteae — the canonical bile-acid 7α-dehydroxylating bacteria — show the predicted substrate-product signature: negative correlation with primary tauro-conjugated bile acids (substrates) and positive correlation with secondary unconjugated bile acids (products: deoxycholate, lithocholate). This is the direct sample-level confirmation of the bile-acid 7α-dehydroxylation network operating in HMP2 samples. By contrast, M. gnavus and E. coli show the OPPOSITE pattern: positive with primary tauro-BAs, negative with secondary BAs — they are not in the 7α-dehydroxylation network.

NB07c verdict reframed: shared-environment co-occurrence

The cross-feeding hypothesis (a) is not supported by sample-level paired evidence. The shared-environment hypothesis (b) is the more parsimonious explanation for A. caccae × pathobiont species-level coupling. NB09c reframes NB07c: the cocktail-design implication originally stated as "depleting pathobionts may incidentally reduce A. caccae" (NB07c §10) is less load-bearing than originally framed — the substrate-product relationship is not detectable at sample level, so phage targeting of pathobionts is unlikely to substantially deplete A. caccae through metabolic coupling.

Bile-acid coupling cost — primary Pillar 4 cocktail-design annotation

Replace "metabolic-coupling cost" (NB07c §10) with bile-acid coupling cost as the primary per-target ecological-cost annotation for the NB05 actionable Tier-A set:

• F. plautii targeting: highest cost. Strong negative correlation with primary cholate (−0.26) AND positive correlation with secondary lithocholate (+0.15) → 7α-dehydroxylation activity. Phage depletion shifts BA pool toward primary inflammatory forms.
• E. bolteae targeting: moderate cost. Positive with deoxycholate (+0.17) and lithocholate (+0.18) → secondary 7α-dehydroxylation contributor.
• E. lenta targeting: moderate cost. Negative with primary cholate (−0.13) and ketodeoxycholate (−0.14) — partial 7α-dehydroxylation pattern, consistent with Eggerthellaceae bile-acid metabolism activity (Devlin & Fischbach 2015 PMID 26412091).
• H. hathewayi / M. gnavus / E. coli targeting: low BA-coupling cost — these species are not in the 7α-dehydroxylation network in this dataset.

Six-line cross-corroboration extended to bile-acid biology

The bile-acid 7α-dehydroxylation network finding is independently corroborated across project pillars:

1. NB07b — F. plautii within-carrier CD-up pathways include F420-dependent bile-acid metabolism pathways
2. NB09a — tauro-α/β-muricholate + free taurine CD-up at subject level (consistent with reduced microbial 7α-dehydroxylation in CD)
3. NB09c — F. plautii / E. lenta / E. bolteae show the direct substrate-product signature in paired sample-level data (negative with primary tauro-BAs, positive with secondary unconjugated BAs)
4. NB05 F. plautii and E. lenta are actionable Tier-A pathobionts
5. Literature (Franzosa 2019, Devlin & Fischbach 2015) — bile-acid 7α-dehydroxylation deficit is a canonical IBD finding

This is a second six-line cross-corroboration narrative alongside the iron-acquisition narrative (Novel Contribution #14): bile-acid biology is the second mechanism axis with multi-line evidence convergence.

Outputs

• data/nb09c_species_metabolite_corr.tsv — 8 species × 583 metabolites paired Spearman ρ + FDR
• data/nb09c_cross_feeding_triangles.tsv — 7 strict cross-feeding-triangle candidates
• data/nb09c_cross_feeding_panel.tsv — curated 7-theme panel direction-of-association profile
• data/nb09c_cross_feeding_verdict.json — formal verdict
• figures/NB09c_cross_feeding_disambiguation.png — heatmap + triangle scatter

(Script: run_nb09c.py. Builds on NB07c module-anchor × pathobiont species coupling + NB09a metabolite annotations. Per plan v1.9 no-raw-reads.)

14. NB10a — Kumbhari strain-adaptation gene predictor (H3b)

The H3b hypothesis (per plan v1.7): strain-adaptation gene content discriminates disease-dominant from health-dominant strains within Kumbhari, and the discriminating genes are biologically interpretable rather than housekeeping. Per plan v1.7 H3b restriction, the test is restricted to the Kumbhari/LSS-PRISM cohort itself (the only project data with strain-level resolution). Per plan v1.9 (no raw reads), uses precomputed ref_kumbhari_s7_gene_strain_inference (219,121 gene × species rows; 59 species) cross-referenced with fact_strain_competition (15,520 strain × sample rows; 100 UHGG genome IDs in 94 LSS-PRISM participants).

Tier-A core intersection with Kumbhari: only F. plautii is in the 59-species panel. The other 5 actionable Tier-A core (H. hathewayi, M. gnavus, E. lenta, E. bolteae, E. coli) are NOT in Kumbhari — Kumbhari focused on commensal strain heterogeneity, not pathobiont mechanism. The H3b test therefore evaluates the methodology of strain-adaptation gene inference on the 59 commensals; only one actionable Tier-A species is directly tested.

Falsifiability test result: among the 23,579 FDR<0.10 strain-adaptation genes (12,159 IBD-biased + 11,420 health-biased across 59 species), the IBD-biased gene set is significantly:

The IBD-biased gene set is 1.38× ENRICHED for adaptation and 0.62× DEPLETED for housekeeping — the canonical SUPPORTED pattern. The strict housekeeping-domination falsifiability bound is rejected at p=6.4e-20.

Adaptation subcategory breakdown (IBD-biased): mucin_glycan utilization (94 genes), two-component signaling (52), antibiotic resistance (35), membrane transport specialty (33), virulence/secretion (15), bile acid (15), oxidative stress response (15). The dominant IBD-adaptation theme is mucin/glycan utilization — consistent with the Kumbhari focus on commensal IBD-niche adaptation.

Cross-species shared IBD-adaptation KOs — top KEGG KOs IBD-biased in 12-21 of 59 species:

The cross-species pattern is dominated by transport, regulation, stress response, and DNA repair — classic niche-adaptation signatures. Top gene symbols include corC (Mg²⁺/Co²⁺ efflux), luxQ (autoinducer-2 quorum sensing), pglA / rfbC (capsule polysaccharide / LPS), bepC (efflux pump), abgB (aminobenzoyl-glutamate utilization) — consistent across 9-10 species each.

F. plautii informative null: zero FDR<0.10 strain-adaptation genes in the Kumbhari analysis (3,245 genes total tested in F. plautii). This is biologically meaningful given F. plautii has confirmed CD-association at species level (NB04e, NB05) and at the bile-acid 7α-dehydroxylation activity level (NB09c §13: ρ × cholate -0.26; ρ × lithocholate +0.15). The interpretation: F. plautii CD-association operates through species-level abundance, not strain-level genomic adaptation. The 7α-dehydroxylation activity is encoded by core bai-operon genes that are presumably present in essentially all F. plautii strains; the CD signal in NB04e + NB07b reflects how much F. plautii (any strain) is present in the sample, not which F. plautii strain is dominant. Independent corroboration of NB07b within-carrier finding (small per-pathway shifts within carriers; CD signal dominated by carriage prevalence).

H3b verdict: SUPPORTED. The Kumbhari strain-adaptation gene-content analysis carries biologically interpretable signal that is statistically distinct from housekeeping artifact. The cross-species shared IBD-adaptation gene set (sigma factor, ABC transport, quorum sensing, capsule LPS, efflux) is consistent with a biologically real "IBD-niche" gene signature operating across multiple commensal species.

Caveats / limitations:
- One cohort only (Kumbhari/LSS-PRISM); cross-cohort sensitivity test was acknowledged-weak per plan v1.7.
- 95 % of FDR-passing genes are unclassified by the gene-symbol regex — the curated 13-category classifier captures only annotated functional roles. The functional enrichment signal is at the category-classified subset, not the full gene set.
- Per plan v1.7, the cross-cohort sensitivity test (NB10b — species-level prevalence-weighted strain-feature on cMD) was specified as sensitivity-only, not load-bearing. Per plan v1.9, the cross-cohort cMD path is dropped (raw reads not available); NB10b is therefore not executed in this project.
- F. plautii CD-association is species-abundance-mediated — implies phage targeting of F. plautii would deplete BA 7α-dehydroxylation activity proportionally to abundance reduction (no strain-level escape route via gene-content selection within F. plautii).

Implications for project narrative:
1. The strain-level analysis methodology is validated on commensal species (H3b SUPPORTED). The actionable Tier-A coverage is single-species (F. plautii only); the analysis is informative for methodology but not for direct pathobiont-mechanism extension.
2. F. plautii CD-association is species-abundance-mediated, confirming NB07b within-carrier reading independently.
3. Cross-species shared IBD-adaptation gene set (sigma factor, ABC transport, quorum sensing, capsule LPS) is a real biological signature — bullet for the Pillar 4 / 5 hand-off discussion of bacterial fitness in the IBD gut.

Outputs

• data/nb10a_per_species_bias_counts.tsv — 59 species × {IBD-biased, health-biased} gene counts
• data/nb10a_sig_genes_classified.tsv — 23,579 FDR-passing genes × functional classification
• data/nb10a_cross_species_ibd_kos.tsv — multi-species shared KEGG KO signal
• data/nb10a_cross_species_ibd_symbols.tsv — multi-species shared gene-symbol signal
• data/nb10a_f_plautii_strain_adaptation.tsv — F. plautii deep-dive (empty, no FDR-passing genes)
• data/nb10a_h3b_verdict.json — formal verdict (SUPPORTED)
• figures/NB10a_kumbhari_strain_adaptation.png — 3-panel summary

(Script: run_nb10a.py. Per plan v1.7 H3b + v1.9 no-raw-reads.)

15. NB11 — HMP2 serology × Tier-A pathobiont (H3e)

The H3e hypothesis (per plan v1.7): anti-microbial antibody titers correlate with Tier-A pathobiont abundance within IBD patients, n=67 subjects across 3 sites (CCHMC/Harvard/Emory), site as covariate, effect threshold |r|>0.40 calibrated to n=67 power 0.80 at α=0.05. Single-cohort caveat structural per plan v1.7. Per plan v1.9: uses cMD-fetched HMP2 MetaPhlAn3 abundance + mart fact_serology (2,520 measurements; 12 assays). All 67 serology subjects (E-prefix Emory, H-prefix Harvard, M-prefix MGH/CCHMC) match cMD MetaPhlAn3 metagenomics.

Six EU-deduplicated serology axes: ANCA, ASCA, CBir1, IgA-ASCA, IgG-ASCA, OmpC. Cross with 8 species (6 actionable Tier-A core + 2 NB07c module-anchor commensals) → 48 tests; partial Pearson r residualized on site dummy variables; BH-FDR.

Cohort sanity check — canonical IBD-serology patterns hold in n=32 CD + 14 UC + 21 nonIBD subjects:

All 6 axes show canonical IBD direction → data validates standard immunology biology. The H3e test asks the harder question: do antibody titers track individual Tier-A pathobiont abundance within IBD patients?

Top moderate associations (none clear |r|>0.40 + FDR<0.10):

Cohort-aggregate signal is dominated by ANCA × M. gnavus (+0.31). Site-stratified breakdown shows most associations are concentrated in 1 of 2 productive sites (Harvard or CCHMC) and dilute when cohort-aggregated. With Harvard 33 + CCHMC 31 + Emory 3 subjects, single-site partial r reaches |r|≈0.40–0.46 individually but the cohort pull toward null prevents the strict |r|>0.40 + FDR<0.10 plan threshold from being met. Top observed: |r|=0.31, FDR=0.40 → strict H3e: NOT SUPPORTED.

Biologically plausible directions (all top pairs are positive r): anti-microbial antibody titers ↑ with target Tier-A species abundance. ANCA × M. gnavus + ANCA × H. hathewayi + ANCA × F. plautii — pANCA antibody (canonically UC-associated) co-elevated with multiple CD-pathobiont species, consistent with high-pathobiont-burden states triggering broader anti-microbial humoral response. CBir1 × E. bolteae mechanistically coherent (CBir1 is anti-bacterial-flagellin; E. bolteae is flagellated). IgA-ASCA × E. coli is the most mucosal-immunity-relevant axis (IgA isotype is the gut-immunity isotype); the +0.23 correlation is consistent with E. coli outer-antigen recognition.

H3e verdict — PARTIAL (NOT SUPPORTED at strict plan threshold). Per plan v1.7 falsifiability bound, no (assay × species) pair clears the |r|>0.40 + FDR<0.10 effect threshold. The PARTIAL framing reflects (a) cohort sanity check passing on canonical IBD-serology patterns, (b) biologically plausible direction of all top pairs, (c) productive sites reaching |r|≈0.40–0.46 individually, and (d) the structural single-cohort caveat acknowledged in plan v1.7.

Implications for Pillar 4–5 cocktail design: serology is not yet a quantitative target-prioritization signal at the per-target level for this cohort. The canonical CD-vs-UC ecotype-stratification axes are intact (ANCA-UC, ASCA-CD, CBir1-CD), so serology can be used as a CD-vs-UC stratifier but not as a per-target abundance predictor — the |r|≈0.31 ceiling at the cohort level means a single subject's anti-M. gnavus ANCA titer cannot be used to confidently predict that subject's M. gnavus abundance. Multi-cohort meta-analysis (out of project scope) is the structural next step.

Limitations:
- n=67 across 3 unevenly-distributed sites is structurally low-power for |r|>0.40 detection at α=0.05 + 48-test multiple-testing correction. Plan v1.7 H3e was acknowledged-weak by design.
- Subject-level aggregation loses longitudinal information (each subject contributed multiple serology + metagenomics visits). Mixed-effects regression would have more power but unbalanced visit counts complicate interpretation.
- OmpC and ANCA cross-react with related antigens beyond their named targets, so per-species attributions are correlative rather than exclusive.
- Single-cohort caveat structural — no project-internal replication path; would need multi-cohort meta to firm the signals up.

Output artifacts:
- data/nb11_serology_species_correlations.tsv — 48 (assay × species) tests with partial Pearson r + raw Spearman ρ + BH-FDR
- data/nb11_serology_site_stratified.tsv — site-stratified breakdown for top 10 pairs
- data/nb11_h3e_verdict.json — formal verdict
- figures/NB11_serology_pathobiont.png — heatmap + top 10 pairs

(Script: run_nb11.py. Per plan v1.7 H3e + v1.9 no-raw-reads.)

16. NB09b — Cross-cohort metabolomics bridge (HMP2 NB09a → FRANZOSA_2019)

NB09a §12 found that polyamines (OR=14.6) + long-chain PUFAs (OR=7.9) are CD-up theme-level signatures in HMP2; tauro-α/β-muricholate + free taurine CD-up corroborates F. plautii 7α-dehydroxylation deficit; urobilin CD-DOWN. NB09b tests whether these signatures replicate in FRANZOSA_2019 (220 participants × 88 CD + 76 UC + 56 Control), an independent multi-site IBD metabolomics cohort. Per plan v1.9 no raw reads — uses precomputed fact_metabolomics only. DAVE_SAMP_METABOLOMICS dropped (mouse ileal AKR/J tissue, not human IBD).

Method: Franzosa metabolite IDs follow [m/z]_[RT] format with no compound-name annotations in the mart. Bridge HMP2 named metabolites (with HMDB ID + m/z + RT) to Franzosa peaks within method by m/z (±0.005 Da tolerance). RT scales differ between HMP2 and Franzosa runs (different gradients/columns) so RT is not used as a matching constraint. 122 of 592 named HMP2 metabolites match a Franzosa peak; 118 of 122 pass Franzosa CD-vs-Control DA filters.

Per-theme cross-cohort sign-concordance:

3 themes ≥75 % sign-concordant (urobilin, acyl-carnitines, long-chain PUFA), overall 64 % concordance (76 of 118), 9 strict cross-cohort replications (both FDR<0.10 + |cliff|>0.20 + sign-match) — these represent the strongest cross-cohort signal in the project's metabolomics analyses.

Findings:

1. Urobilin CD-DOWN replicates at 100 % (3/3 matches; median Franzosa cliff −0.747 — the largest cross-cohort effect size in the panel). Reflects loss of bilirubin-reducing commensals (Hall 2024; Vital 2018) consistently in CD across HMP2 and Franzosa. Strongest replicated metabolite-level finding in the project.

2. Acyl-carnitines CD-up replicates at 80 % (4/5 matches; median Franzosa cliff +0.531). The C16 + C18:1 long-chain acyl-carnitines that NB09a §12 flagged as CD-up are independently elevated in Franzosa CD samples. Mechanistically connected to v1.8 §9 fat-metabolism/glyoxylate theme + impaired β-oxidation in inflamed mucosa.

3. Long-chain PUFAs CD-up replicates at 75 % (9/12 matches; median Franzosa cliff +0.273). Adrenate/arachidonate/DHA/DPA replicate cleanly; eicosapentaenoate (EPA) shows opposite sign in 1 of 3 non-concordant cases (possibly cohort-specific dietary fish-oil supplementation differences). Inflammatory eicosanoid-precursor pool elevation is a robust cross-cohort signature.

4. Bile-acid signals are weaker cross-cohort than within-cohort. Primary BAs replicate at 57 % (4/7), secondary BAs at 40 % (2/5). The NB09c §13 paired sample-level substrate-product signature for F. plautii / E. lenta / E. bolteae (ρ within-HMP2-paired) is the more robust evidence stream than cohort-aggregate cross-cohort DA. BA pool sizes are more susceptible to dietary and antibiotic confounding than within-cohort species × metabolite correlation. NB09c §13 finding does not require Franzosa replication; this NB09b result documents that BA-pool DA is more cohort-dependent than within-cohort correlation.

5. Polyamine bridge fails — m/z-only matching does not recover any of the 9 polyamines that drove the NB09a §12 OR=14.6 theme enrichment. Polyamines are low-mass HILIC-pos compounds (putrescine 89.11; spermidine 146.17; spermine 203.22) and the m/z-tolerance + Franzosa-method-library overlap does not yield matches. The HMP2 polyamine OR=14.6 stands as a single-cohort finding until cross-cohort replication via external Franzosa supplementary tables (with compound-name annotations) becomes available — outside the project's no-raw-reads scope.

Verdict — STRONG cross-cohort replication on 3 of 8 testable themes (urobilin + acyl-carnitines + long-chain PUFAs). Overall 64 % sign-concordance is moderately above the 50 % chance baseline; the 9 strict replications are firmly above-chance. Polyamines unable to bridge by m/z. The bile-acid DA signal is partly cohort-specific (which strengthens the NB09c §13 paired-sample interpretation, where the F. plautii / E. lenta / E. bolteae substrate-product signature is direct rather than DA-aggregated).

Bile-acid 7α-dehydroxylation narrative — sample-level evidence is stronger than cohort-aggregate cross-cohort

NB09b clarifies the relative robustness of two evidence streams for the bile-acid network finding:
- NB09c §13 paired sample-level direct substrate-product signature (ρ(F. plautii × cholate) = −0.26; × lithocholate = +0.15) operates within paired metabolomics+metagenomics samples in HMP2. It is a species × metabolite correlation, not a CD-vs-Control aggregate DA.
- NB09b cohort-aggregate cross-cohort DA on bile acids replicates at 57 % primary / 40 % secondary — moderately above chance for primary, sub-chance for secondary.

The interpretation: BA pool sizes vary across cohorts due to dietary and antibiotic confounding, but the mechanistic substrate-product signature for active 7α-dehydroxylation is preserved at the within-paired-sample level. The NB09c six-line bile-acid narrative is anchored on the paired sample-level evidence + literature mechanism, not on cohort-aggregate DA replication. This is a methodological note for any future BERIL multi-modal analysis: cohort-aggregate cross-cohort DA and within-cohort paired-sample correlation are different evidence streams; one can be weaker without invalidating the other.

Outputs

• data/nb09b_cross_cohort_concordance.tsv — 122 (HMP2-name × Franzosa-peak) matched pairs with cliff_delta + sign-match flag
• data/nb09b_theme_replication.tsv — per-theme sign-concordance summary
• data/nb09b_cross_cohort_verdict.json — formal verdict
• figures/NB09b_cross_cohort_metabolomics.png — cross-cohort cliff δ scatter + per-theme sign-concordance

(Script: run_nb09b.py. Per plan v1.7 H3d-DA cross-cohort follow-up + v1.9 no-raw-reads.)

17. NB09d — H3d-clust metabolite-feature ecotype LOSO stability

The H3d-clust hypothesis (per plan v1.7): does metabolite-derived clustering achieve higher cross-cohort stability than the taxonomic-derived ecotype framework? If yes, refitting ecotypes on metabolite features would give a more cross-cohort-portable framework for clinical translation. Per plan v1.9 no raw reads — uses NB09b m/z-bridge feature panel + pooled HMP2 + Franzosa metabolomics.

Method: pool 326 subjects (106 HMP2 × 220 Franzosa) on the 111 unique m/z-bridge metabolite features (NB09b); standardize per metabolite; PCA → 15 components; K-means K=4. Compare to NB04f taxonomic-ecotype LOSO ARI baseline (0.113) and NB04b bootstrap baseline (0.16).

Cluster structure is COHORT-driven, not diagnosis-driven:

PC1 explains 79 % of total variance and separates HMP2 (PC1 ≈ +12) from Franzosa (PC1 ≈ −5) cleanly. PCs 2–15 collectively explain ~10 % of variance and carry the diagnosis information. Within each cohort, clusters DO carry diagnosis information (cluster × diagnosis χ²(6)=88.3, p=7e-17), but cluster boundaries are mostly co-aligned with cohort boundaries.

Stability metrics:

The within-pooled bootstrap ARI of 0.937 is misleading as a stability metric: it measures how reproducibly the pooled K-means recovers the cohort-batch structure under subsampling, which is high. It does NOT measure whether the clustering is biologically informative; the cohort batch effect is consistently visible in any subsample.

H3d-clust verdict: NOT SUPPORTED. Cross-cohort LOSO ARI = 0.000 ≪ 0.113 taxonomic baseline. The metabolite-feature framework on the m/z-bridge feature panel does NOT achieve higher cross-cohort stability than the taxonomic framework.

Methodological lesson (now in docs/discoveries.md): m/z-bridge metabolomics has unaddressed batch effects that dominate cluster structure. The taxonomic-feature ecotype framework was naturally cross-cohort-portable because MetaPhlAn3 relative-abundance values are unitless and compositionally constrained — that constraint creates a natural normalization that handles cohort differences. Metabolite-feature clustering inherits absolute-intensity scale differences between LC-MS runs (different instruments, ionization tuning, solvent batches) and requires explicit batch correction (ComBat / SVA / RUV / quantile normalization) prior to clustering. Within-cohort metabolite-feature ecotypes might still be useful for single-laboratory clinical translation; the failure here is specifically the cross-cohort-portability question.

Implications for project narrative

• The taxonomic ecotype framework (NB01b consensus K=4) remains the project's primary ecotype basis for Pillar 4–5 cocktail design.
• Metabolite-feature ecotype refit is not a drop-in replacement for the taxonomic framework in this project's data scope; it would require external batch correction infrastructure.
• The NB09c §13 paired sample-level finding (within-HMP2 species × metabolite correlation) is the more robust mechanism evidence stream; cohort-aggregate cross-cohort metabolite DA is dominated by batch effects (NB09b primary BAs 57 % concord; secondary BAs 40 %).
• The H3 falsifiability framework is now fully closed: H3a (a/c) + H3a (b) v1.8 + H3a-new + H3b + H3c + H3d-DA + H3d-clust + H3e all tested with formal verdicts.

Limitations:
- m/z-bridge feature space is small (111 unique features). A larger feature space (500-1000 metabolites with name annotations across cohorts) might give different results — but requires external Franzosa supplementary annotations (out of scope per plan v1.9).
- No batch correction applied — ComBat/SVA/RUV would likely improve cross-cohort LOSO ARI substantially. The current verdict is for "raw m/z-bridge clustering without batch correction."
- K=4 is the plan-imposed K matching the taxonomic framework; natural cluster count on this feature space might be K=2 (HMP2 vs Franzosa).

Output artifacts:
- data/nb09d_metabolite_ecotype_assignments.tsv — 326 pooled subjects × {cohort, diagnosis, metabolite_cluster}
- data/nb09d_h3d_clust_verdict.json — formal H3d-clust verdict (NOT SUPPORTED)
- figures/NB09d_metabolite_ecotype_stability.png — 3-panel summary

(Script: run_nb09d.py. Per plan v1.7 H3d-clust + v1.9 no-raw-reads.)

18. NB07d — MOFA+-style HMP2 multi-omics joint factor pilot (taxonomy + metabolomics)

The plan v1.7 NB07d called for 3-modality MOFA+ on HMP2 taxonomy + pathways + metabolomics. Per plan v1.9 (no raw reads) and the project's data-scope reality — HMP2 pathway abundance is not in the mart (fact_pathway_abundance contains only CMD_IBD_PATHWAYS) — NB07d falls back to 2-modality joint factor analysis on (HMP2 taxonomy + HMP2 metabolomics) using sklearn CCA on per-modality PC scores (mofapy2 unavailable).

Data scope: 106 paired HMP2 CSM* subjects (intersection from NB09c) × 130 ≥10%-prevalence species (CLR-transformed) + 582 named metabolites (log10-intensity, ≥30% non-NaN coverage). Per-modality PCA → 30 components each → CCA → 4 canonical pairs.

All 4 canonical correlations are very strong: r = 0.964, 0.928, 0.911, 0.889 — taxonomy and metabolomics share substantial cross-modality structure on this paired-subject set.

CC1 is THE CD-vs-nonIBD diagnosis-discriminative joint factor:

CC1 separates CD (+0.235) from nonIBD (−0.593) by ~0.83 SD on a single joint axis; UC sits at +0.123, between CD and nonIBD (consistent with UC as a milder dysbiosis state on the same axis).

CC1 species loadings recapitulate the entire actionable Tier-A set: ALL 6 actionable Tier-A core species load POSITIVE (CD-direction): M. gnavus +0.195, E. coli +0.173, F. plautii +0.153, H. hathewayi +0.144, E. lenta +0.109, E. bolteae +0.103. Plus oral-gut Tier-B candidates V. parvula +0.194, A. intestini +0.161, V. atypica +0.151. NEGATIVE loadings (commensal-direction) match NB01b ecotype-defining commensals: R. bromii −0.173, R. bicirculans −0.170, A. putredinis −0.169, L. eligens, B. intestinihominis. NB06 H2d pathobiont module structure independently rediscovered as a single principal direction in joint species-metabolite space — two analytical approaches (CLR+Spearman+Louvain modules vs CCA joint factors) converge on the same biology.

CC1 metabolite loadings recapitulate ALL Pillar 3 metabolomics narratives in one axis:

Pilot verdict: SUCCESSFUL. CC1 is the unified Pillar 3 CD-vs-nonIBD axis in joint species-metabolite space. It is the cleanest single-factor representation of "what is CD biology" that the project has produced — a single principal direction that:
1. recapitulates the entire NB05 actionable Tier-A set (all 6 species CD-positive)
2. recapitulates the NB06 H2d pathobiont module structure
3. recapitulates the NB09a §12 polyamine + PUFA themes
4. recapitulates the NB09b §16 urobilin cross-cohort CD-DOWN finding
5. recapitulates the NB09c §13 BA-pool depletion + cadaverine × E. coli signature
6. recapitulates the NB08a §11 ebf/ecf fatty-acid-amide BGC-family signature

This is the cross-pillar unified narrative axis in a single multi-omics joint factor. The same biological state — "CD pathobiont-module dominance at the species level + the metabolic consequences thereof" — is the dominant cross-modality structure in HMP2 and is detectable as a single canonical correlation pair at r=0.96.

Methodological notes

• 2-modality vs 3-modality: pathway modality dropped per data-scope constraint. Adding pathways via cMD_IBD reprocessing on cohort-aligned subjects (NB07a) would give 3-modality MOFA but requires sample-level pairing (cMD_IBD pathway slice exists but is NOT paired with HMP2 metabolomics).
• CCA vs MOFA+: CCA on PC scores captures the same canonical-correlation signal that MOFA+ would on the 2-modality case. MOFA+ would additionally model modality-specific factors (factors loading only on one modality), which CCA does not. The 4 PCs we used are joint factors; modality-specific structure remains in the per-modality PCA components.
• Ecotype as covariate (per plan v1.7 N13): ecotype is implicit in the species-loading structure (E1-Bact2-transitional + E3-Bacteroides-expanded species both load positive on CC1; E0-commensal species load negative). Future extension: regress factor scores on ecotype + diagnosis to test whether residual factor variance encodes ecotype-specific biology.

Limitations

• Pathway modality not in the analysis (HMP2 fact_pathway_abundance unavailable; v1.9 no-raw-reads constraint).
• CCA cannot identify modality-specific factors — MOFA+ would add ~5–10 unique factors per modality on top of joint factors.
• Pilot-scale — not a substitute for a full MOFA+ analysis with 3 modalities and proper factor-relevance gating. The 4 canonical correlations are all very strong (>0.88), which suggests the per-modality PCAs are over-aligned; a proper MOFA+ would impose sparsity to gate factor count.
• n=106 subjects, modest for a 4-factor decomposition. The CC2/CC3/CC4 signal beyond CC1 is suggestive but not strongly disease-discriminative.

Output artifacts:
- data/nb07d_cca_loadings.tsv — top 15 species + 15 metabolite loadings per of 4 canonical components
- data/nb07d_subject_factor_scores.tsv — 106 subjects × {CC1-4 species score, metab score, joint score} + diagnosis
- data/nb07d_mofa_pilot_verdict.json — formal pilot verdict
- figures/NB07d_mofa_pilot.png — joint factor space + CC1 species + CC1 metabolite loadings

(Script: run_nb07d.py. Per plan v1.7 NB07d (scope-adjusted to 2-modality per v1.9) + N13 multi-modal QC prerequisite (NB04e + NB07a + NB09a all passed).)

19. NB12 — Pathobiont × phage targetability matrix (Pillar 4 opener)

The Pillar 4 framework starts with a per-pathobiont phage-availability profile (Tier-B in the project's 4-tier rubric) layered on top of the NB05 actionable Tier-A scoring + Pillar-3 mechanism profile. NB12 builds the foundational matrix from ref_phage_biology (12-organism literature-curated synthesis from indicator_taxa_literature_review).

Phage-availability score (Tier-B), 0-3 ordinal scale: 0 = no known phages, 1 = temperate / limited, 2 = lytic literature, 3 = clinical trial / commercial cocktail.

Per-actionable Tier-A profile:

Stratification — 4 phage-availability classes among 6 actionable Tier-A: Class 3 clinical (1: E. coli), Class 2 lytic literature (2: E. lenta, E. bolteae), Class 1 temperate-limited (1: M. gnavus), Class 0 gap (2: H. hathewayi, F. plautii).

Critical: the 2 highest-NB05-scored species (H. hathewayi 4.0, M. gnavus 3.8) have the WEAKEST phage availability. F. plautii additionally has the HIGHEST BA-coupling cost — phage targeting may be deprioritized in favor of BA-binding co-therapy.

Pillar 4 → Pillar 5 hand-off framework — 3 design strategies for the 6 actionable Tier-A:
1. Direct phage targeting (Tier-1): E. coli (AIEC subset, EcoActive clinical-trial cocktail; require strain-resolution diagnostic per NB07b/NB08a).
2. Phage targeting with monitoring (Tier-2): E. lenta (PMBT5), E. bolteae (PMBT24); M. gnavus if lytic-locked engineering succeeds.
3. Phage GAP — alternative strategies needed:
- H. hathewayi: highest priority for external DB query; if no phages, fall back to GAG-degrading enzyme inhibitors.
- F. plautii: lowest Pillar-5 priority due to highest BA-cost — consider deprioritizing or co-administering UDCA / BA-binding agent.
- M. gnavus: biochemical glucorhamnan-synthesis targets (Henke 2019) if lytic-locked engineering fails.

Output artifacts:
- data/nb12_phage_targetability_matrix.tsv — per-pathobiont scoring matrix with NB05 + phage + Pillar-5 class
- data/nb12_phage_targetability_verdict.json — formal verdict + Pillar-4/5 hand-off + limitations
- figures/NB12_phage_targetability.png — 2-panel: Tier-A × Tier-B scatter + Pillar-5 priority bar

Quantitative-augmentation data found during NB12 execution (deferred to NB13):
- HMP2 fact_viromics (3,039 sample-rows): direct phage observations in patient stool. 16 unique E. coli phages + 13 unique Klebsiella phages observed. NO observations of phages targeting H. hathewayi / M. gnavus / F. plautii / E. lenta / E. bolteae — confirms the gut-anaerobe phage coverage gap from the literature-curation side.
- PhageFoundry phagefoundry_strain_modelling (BERDL): 96 phages × 188 E. coli strains × experimentally-tested susceptibility (17,672 pairs; 3,929 susceptible / 13,743 resistant; 22 % susceptibility rate). The definitive E. coli AIEC phage-cocktail design dataset for this project.
- ref_viromics_cd_vs_nonibd (22 viruses pre-computed CD-vs-nonIBD DA): top hit Gokushovirus WZ-2015a CD-DOWN log2FC=-2.7, FDR=1e-11 (microvirus depleted in CD).

(Script: run_nb12.py. Per plan v1.7 NB12 + v1.9 no-raw-reads. BERDL Spark auth restored mid-execution; PhageFoundry quantitative analysis promoted to NB13.)

20. NB13 — PhageFoundry quantitative E. coli phage-cocktail design

NB12 established the qualitative literature-curated foundation; NB13 adds the quantitative cocktail-design layer at strain-resolution level using BERDL phagefoundry_strain_modelling (96 phages × 188 E. coli strains × 17,672 experimentally-tested susceptibility pairs from Gaborieau 2025-10-02 phage-prediction experiment, AUC=0.88). 3,929 / 17,672 = 22 % susceptibility rate across the matrix.

Greedy minimum-set-cover cocktail design — smallest cocktail covering ≥X % of 188 strains:

A 5-phage cocktail covers 94.7 % of 188 E. coli strains in PhageFoundry. Comparable in size to the EcoActive 7-phage clinical-trial cocktail (per NB12 §1) — both are 5-7 phages. The top-ranked phage DIJ07_P2 alone (genus Phapecoctavirus) lyses 63.8 % of strains — remarkable single-phage breadth.

Phage host-phylogroup analysis — AIEC relevance:

AIEC strains are predominantly phylogroup B2 (~80 %) and D (~20 %) per Dogan 2014/Dubinsky 2022. PhageFoundry distribution:

65 of 94 phages (69 %) are isolated against B2/D phylogroup hosts — strongly AIEC-relevant. Top broadest-host-range B2/D phages include several isolated against canonical AIEC reference strains: LF82_P8 (60.1 %, Mosigvirus) — LF82 is THE AIEC reference (Darfeuille-Michaud 2004); LF73_P1 (62.8 %) — LF73 is also a CD-associated AIEC strain; 536_P7 / 536_P9 — E. coli 536 is a UPEC/B2 archetype.

HMP2 viromics × PhageFoundry overlap = 0 — the 7 unique E. coli phages observed in HMP2 fact_viromics (D108, EC6, ECML-117, Murica, slur16, vB_EcoM-VpaE1, vB_EcoM_AYO145A) do NOT name-overlap with PhageFoundry phages. The two datasets are complementary, not overlapping: PhageFoundry = research/clinical isolates with experimental susceptibility; HMP2 viromics = natural phages observed in patient stool. PhageFoundry is the primary source for cocktail-design.

Pillar 5 hand-off — concrete E. coli AIEC phage-cocktail recommendation:
1. Tier-1 cocktail (5 phages, 95 % strain coverage): DIJ07_P2 + LF73_P1 + AL505_Ev3 + 55989_P2 + LF110_P2
2. Extended Tier-1+ (8 phages, 99 % strain coverage): above + NIC06_P2 + LF73_P4 + BCH953_P4 — comparable in size to EcoActive
3. Strain-level diagnostic requirement: AIEC-vs-commensal E. coli discrimination via per-patient pks-island / Yersiniabactin / Enterobactin gene-presence detection (per NB07b + NB08a) is required before cocktail selection — the cocktail covers broad E. coli strain diversity but does not by itself distinguish AIEC from commensal E. coli.

Limitations:
- PhageFoundry has no explicit AIEC-vs-commensal strain annotation. The 188 strains include canonical AIEC strains (LF82, LF73, 536, NIC06, BCH953 — named as phage isolation hosts) but Pillar-5 cross-reference to literature isolate annotations is needed for accurate AIEC coverage estimation.
- Susceptibility matrix is binary (no titer/burst size). Real-world cocktail dosing depends on phage burst size + receptor binding kinetics + gastric passage survival — out of project scope.
- 26 strains (14 %) are phage-resistant at ≤5 % susceptibility — clinical significance of these escape strains needs cross-reference to AIEC pks-island annotation.
- HMP2 viromics × PhageFoundry overlap = 0 means in-vivo gut delivery of cocktail phages is not validated by the project's data; would require in-vivo testing.

Output artifacts:
- data/nb13_phage_host_range.tsv — 96 phages × {host_range_pct, Family, Genus, Phage_host_phylo, …}
- data/nb13_strain_phage_susceptibility.tsv — 188 strains × phage_susceptibility_pct
- data/nb13_phagefoundry_cocktail_verdict.json — formal verdict + cocktail compositions at 50/75/90/95/99 % coverage
- figures/NB13_phagefoundry_cocktail.png — 3-panel: phage host-range + strain susceptibility + cocktail-coverage curve

(Script: run_nb13.py. BERDL Spark Connect via fresh KBASE_AUTH_TOKEN from /home/aparkin/.env. Per plan v1.7 NB13 + v1.9 no-raw-reads.)

21. NB14 — HMP2 endogenous phageome × ecotype × diagnosis

NB14 closes Pillar 4 by adding the in-vivo phage-community lens to the curated-literature (NB12) + experimental-susceptibility (NB13) layers. HMP2 fact_viromics (3,039 sample-rows × 273 viruses × VirMAP taxonomic profile) cross-referenced with NB04h ecotype projections gives 630 of 648 viromics samples (97 %) with ecotype calls — 27 (E0) + 484 (E1) + 55 (E2) + 50 (E3) samples spanning CD/UC/healthy.

Per-ecotype × per-virus CD-vs-nonIBD Mann-Whitney: 85 DA tests across 4 ecotypes; 3 pass strict FDR<0.10:

Gokushovirus WZ-2015a is robustly CD-DOWN across multiple ecotypes, with the strongest signal in E1 (cliff=-0.36, FDR=5e-7). This independently rediscovers the precomputed ref_viromics_cd_vs_nonibd top hit (log2fc=-2.7, FDR=1e-11). Gokushovirus is a Microviridae member (single-stranded DNA, ~5 kb genome; Gokushovirinae subfamily includes the canonical "crassphage-like" lineage that infects Bacteroides / Prevotella gut commensals) — CD-DOWN consistent with Norman 2015 / Clooney 2019 IBD-Microviridae depletion.

Tier-A pathobiont species × phage-family Spearman ρ (n=630 paired viromics+metaphlan3 samples):

E. coli correlates positively with Podoviridae (+0.18) and Myoviridae (+0.13) — both Caudovirales families that include T7-like + T4-like E. coli phages. This is a plausible endogenous phage-host correlation: when E. coli is abundant, E. coli phages tend to also be abundant (commensal phage carriage / lysogenic-state co-occurrence). All correlations |ρ|≤0.18 — no strong (|ρ|>0.30) endogenous phage candidates targeting Tier-A pathobionts.

H. hathewayi and *M. gnavus* correlate NEGATIVELY with the "Unknown" phage family (which captures 80 % of HMP2 viromics observations that VirMAP couldn't classify) — pathobiont blooms displace some unclassified phages, consistent with reduced phage diversity in CD dysbiosis.

Per-ecotype phage-family abundance (Panel A) shows modest ecotype-specific variation: Anelloviridae E1-specific; Parvoviridae E2-elevated; Unknown family E2-dominant. The dominant signal is the "Unknown" classification (80 % of observations) — VirMAP family-level classification gap is the methodological limit, not biology.

Pillar 4 closure synthesis — three complementary phage-evidence layers:

Combined Pillar-4 verdict: phage-therapy feasibility for E. coli is high (clinical cocktail + experimental susceptibility + in-vivo phage correlation all aligned); for M. gnavus / H. hathewayi / F. plautii the gap is confirmed across all three layers (no clinical phages, no PhageFoundry coverage, no HMP2 in-vivo Podoviridae/Myoviridae signal). External DB queries (INPHARED + IMG/VR) for the 3 gut-anaerobe gaps remain the highest-priority Pillar-4 follow-up before Pillar 5 cocktail drafts.

Limitations:
- VirMAP family-level classification gap: 80 % of phages classified as "Unknown" family — limits per-family DA + correlation power.
- Per-ecotype DA asymmetric power: E1 has 231 CD samples, but E0/E2/E3 have 12-30 — limits within-ecotype detection in smaller ecotypes (Gokushovirus E0 cliff=-0.55 is biologically interesting but doesn't pass strict FDR<0.10 due to power).
- No CRISPR-spacer-derived phage-host predictions in HMP2 viromics — would require IMG/VR cross-reference for Tier-A-host phage discovery.
- "Unknown" family aggregation hides species-level signal; per-virus DA is more informative than per-family for IBD-relevant phages.

Output artifacts:
- data/nb14_viromics_da_per_ecotype.tsv — 85 (ecotype × virus) DA tests
- data/nb14_pathobiont_phage_family_corr.tsv — 126 (Tier-A species × phage family) Spearman ρ
- data/nb14_endogenous_phageome_verdict.json — formal verdict
- figures/NB14_endogenous_phageome.png — 3-panel: ecotype × family abundance heatmap + top 12 DA viruses + species × family ρ heatmap

(Script: run_nb14.py. Per plan v1.7 NB14 endogenous phageome stratification + v1.9 no-raw-reads.)

22. NB15 — UC Davis per-patient profile + cocktail draft (Pillar 5 opener)

NB15 assembles per-patient profiles for 23 UC Davis CD patients combining all Pillar 1-4 evidence and produces concrete per-patient cocktail drafts with components + caveats per patient. The Pillar-5 hand-off framework from REPORT § Pillar 4 closure (3 design strategies) is applied as a per-patient rule book.

Inputs per patient: NB02 ecotype (E0/E1/E3/mixed) + demographics (Montreal, calprotectin, medication) + Kuehl_WGS Kaiju Tier-A pathobiont presence (≥0.001 relative abundance) + NB05 Tier-A score + Pillar-3 mechanism profile + Pillar-4 phage availability.

Per-target prescribing rate (across 23 patients):

E. coli present in only 8/23 patients (35 %) — interesting given AIEC is canonically CD-associated. Implication: NB13's 5-phage AIEC cocktail directly applicable to ~35 % of cohort; H. hathewayi / M. gnavus near-universal carriage means their phage GAP (NB12) is the dominant unmet need.

Per-ecotype cocktail summary:

14 of 23 patients have concrete phage cocktail drafts. All 9 E1 patients receive concrete cocktails (PMBT24 + PMBT5 + AIEC 5-phage if E. coli present); 4 of 6 E3 patients; 1 mixed (patient 6967); E0 patients lack concrete components because their priority targets (H. hathewayi, M. gnavus) are in Pillar-4 GAP / temperate-only.

Per-patient stratification — 4 cocktail-design categories:

F. plautii BA-cost is the dominant E1 design constraint: F. plautii present in 78 % of patients (18/23) AND has the HIGHEST BA-coupling cost (NB09c §13 active 7α-dehydroxylator). All 9 E1 patients carry F. plautii. Phage targeting shifts BA pool toward inflammatory primary tauro-conjugated forms; alternative is deprioritize F. plautii from cocktail + co-administer UDCA / BA-binding agent. F. plautii is also Pillar-4 GAP (not in ref_phage_biology), making it a "leave alone" target for Pillar 5.

E1 hybrid cocktail (9 patients): 3 strategies blended:
- Direct phage: E. coli (if present, NB13 5-phage), E. bolteae (PMBT24), E. lenta (PMBT5)
- Alternative: H. hathewayi (GAG-degrading enzyme inhibitors), F. plautii (BA-binding co-therapy or deprioritize)
- Limited: M. gnavus (temperate-only; lytic-locked engineering OR biochemical glucorhamnan-synthesis target)

Pure phage cocktail is not feasible for any E1 patient — only a hybrid with non-phage alternatives. The 3 GAP species require non-phage strategies; INPHARED + IMG/VR external DB queries for new gut-anaerobe phages are the highest-priority Pillar-4 follow-up before clinical translation.

Pillar 5 hand-off — concrete deliverables:

Limitations:
- 23-patient cohort is small — results are exemplars/templates, not statistically robust per-patient validation.
- Kuehl_WGS uses Kaiju (not MetaPhlAn3) — Tier-A presence calls have lower confidence vs CMD analyses (NB02 classifier-mismatch asymmetry).
- No per-patient bile-acid measurements; BA-cost is ecotype-level, not per-patient.
- No patient-specific AIEC strain-resolution diagnostic for the 8 E. coli-positive patients; cocktail recommendation assumes AIEC subset prevalence per Dogan 2014.
- F. plautii deprioritization is precautionary based on NB09c mechanism — clinical validation requires per-patient BA panels.

Output artifacts:
- data/nb15_patient_profile.tsv — 23 patients × full per-patient profile
- data/nb15_per_patient_cocktail_draft.tsv — long-format per-patient × per-target cocktail breakdown
- data/nb15_pillar5_cocktail_verdict.json — formal verdict + per-ecotype summary
- figures/NB15_patient_cocktail_draft.png — 3-panel: presence heatmap + ecotype prescribing rate + calp × n_targets scatter

(Script: run_nb15.py. Pillar 5 opener; integrates NB02 + NB05 + NB06 + NB09c + NB10a + NB12-NB14 evidence into per-patient cocktail recommendations. Per plan v1.9 no-raw-reads.)

23. NB16 — Patient 6967 longitudinal stability + state-dependent dosing strategy

NB16 is the central per-patient longitudinal stability test for the project, focused on patient 6967 — the only multi-timepoint UC Davis CD patient with documented E1↔E3 ecotype drift across 2 visits. The notebook also uses patient 1112 (2 reseq replicates of the same biological sample) as a Kaiju reliability validation.

Patient 6967 longitudinal — clear E1→E3 transition with M. gnavus dominant:

M. gnavus 14× expansion is the dominant signature of the E1→E3 transition. All other Tier-A pathobionts also expand (1.3–3.1×) reflecting general dysbiosis worsening, but the M. gnavus fold-change dominates. E. coli remains absent in both visits — patient 6967 is not an AIEC carrier. This matches NB01b ecotype biology: E3 = severe-Bacteroides-expanded with M. gnavus as a dominant expansion axis; E1 = Bact2-transitional with milder pathobiont burden.

Per-visit cocktail composition (would the cocktail change?):

• Visit 1 (E1): 5 priority targets (E1_CD module 0), all 5 present → cocktail = {H. hathewayi, M. gnavus, E. bolteae, E. lenta, F. plautii}.
• Visit 2 (E3): 4 priority targets (E3_CD module 1), 3 present (E. coli absent) → cocktail = {E. lenta, H. hathewayi, M. gnavus}.
• Shared (both visits): H. hathewayi, M. gnavus, E. lenta — universal Tier-1 trio.
• Visit-1-only: E. bolteae, F. plautii — E1-specific (would be deprioritized on E3 transition).
• Cocktail Jaccard (visit 1 × visit 2) = 0.60 — moderate overlap; ecotype drift implies non-trivial cocktail re-design, not a complete rewrite.

Patient 1112 technical replicate concordance — Kaiju reliability validation:

Spearman ρ = 1.000 (p < 0.001) on 6 Tier-A — perfect rank concordance across reseq replicates. Kaiju calls are highly reliable for the actionable Tier-A pathobionts in UC Davis samples.

State-dependent dosing rule (5 concrete recommendations):

1. Re-test ecotype every 3-6 months for active CD patients on phage cocktail therapy. Patient 6967 demonstrates ecotype is dynamic, not static.
2. F. plautii inclusion is E1-specific — if patient transitions E1→E3, deprioritize F. plautii from cocktail. This also reduces BA-coupling-cost concern (NB09c). Per NB07c, F. plautii is an E1_CD module 0 species, not an E3_CD module 1 species.
3. E. coli inclusion is E3-specific — if patient transitions E3→E1, deprioritize E. coli (subject to AIEC strain detection per NB07b/NB08a). E. coli is in E3_CD module 1 priority list, not E1_CD.
4. Universal Tier-1 trio (M. gnavus, H. hathewayi, E. lenta) span both E1 and E3 — these don't need re-evaluation on ecotype shift; they are the cocktail backbone for any active-disease CD patient.
5. M. gnavus qPCR as cheap ecotype-state indicator — the 14× expansion in patient 6967's E3 transition suggests M. gnavus abundance via qPCR could serve as a clinical proxy for ecotype shift, avoiding the need for full metagenomics at every visit. A 5-fold change in M. gnavus abundance might be the threshold for triggering full ecotype re-evaluation.

Pillar 5 hand-off — clinical workflow recommendation:

Initial visit:
  └─ Stool metagenomics → ecotype assignment
      ├─ E0: limited cocktail, flare reserve
      ├─ E1: full hybrid 3-strategy cocktail (NB13 5-phage E. coli if AIEC+;
      │       PMBT24; PMBT5; alternatives for H. hathewayi + F. plautii)
      └─ E3: focused cocktail (E. coli if present; PMBT5; alternatives)

Follow-up visits (3-6 month):
  ├─ Calprotectin: assess disease activity
  ├─ M. gnavus qPCR: cheap ecotype-state indicator
  │   └─ if 5-fold change → trigger full ecotype re-test
  └─ If full re-test shows ecotype shift:
      ├─ E1 → E3: drop F. plautii; consider adding E. coli cocktail
      ├─ E3 → E1: add F. plautii alternative; reassess E. coli targeting
      └─ Stable ecotype: continue current cocktail

Limitations:
- Patient 6967 is the only multi-timepoint patient with biological replicate samples in the UC Davis cohort — ecotype-drift conclusions are based on n=1 longitudinal trajectory.
- Visit 1 confidence (0.64) and visit 2 confidence (0.41) are both moderate; the ecotype call shift could partly reflect classifier uncertainty. However, the 14× M. gnavus expansion and 3× E. lenta expansion are quantitative biological observations independent of the ecotype label.
- No timing information between the 2 patient 6967 visits — duration of E1→E3 drift is unknown, limiting clinical-workflow timing recommendations.
- State-dependent dosing rule is theoretical; not yet validated in clinical practice.
- Cocktail Jaccard 0.60 is moderate; whether 3-of-5 component overlap implies clinically-meaningful cocktail re-design depends on patient response heterogeneity (out of project scope).

Output artifacts:
- data/nb16_p6967_tier_a_longitudinal.tsv — patient 6967 per-visit Tier-A abundance + fold change
- data/nb16_longitudinal_verdict.json — formal verdict + state-dependent dosing rules
- figures/NB16_longitudinal_dosing.png — 3-panel: per-visit Tier-A + tech replicate scatter + cocktail composition shift

(Script: run_nb16.py. Pillar 5 second notebook — central longitudinal-stability test + state-dependent dosing rule for Pillar 5 hand-off. Per plan v1.9 no-raw-reads.)