diff --git a/docs/index.qmd b/docs/index.qmd index 4b20dbe..df4f92c 100644 --- a/docs/index.qmd +++ b/docs/index.qmd @@ -76,20 +76,82 @@ html`
```{ojs} { - const pools = data.pools.filter(p => !p.simulated); - const classes = data.candidates.by_evidence_class.filter(d => !d.simulated); - const poolText = pools.length - ? pools.map(p => `${p.studies.join(" and ")} pool on ${phen(p.phenotype)}: ${fmt(p.n_features_k2plus)} genes measured in both, ${fmt(p.n_features_k2plus_significant)} significant after multiple-testing correction`).join("; ") - : "no cross-study pool has formed yet"; - const classText = classes.map(d => `${fmt(d.n)} ${d.evidence_class.replace(/_/g, " ")}`).join(", "); - const gap = k.n_real_comparisons_with_resilience_phenotype === 0 - ? "No real comparison measures a resilience outcome yet." - : `${k.n_real_comparisons_with_resilience_phenotype} of ${k.n_comparisons_real} real comparisons measure a resilience outcome.`; - return html`

${poolText}. Ranked candidates by evidence class: ${classText}.

-

${gap} All real studies so far are RNA-seq pathogen challenges, so the ranked list reflects transcriptional - response to infection, not resistance. A candidate here is an association across the available evidence, not a - validated biomarker. The next curation step is a second stressor class and a study with a measured survival or - pathogen-load phenotype (candidate studies).

`; + // Every sentence below is computed from data.json, so it stays current as studies are added. + const ASSAY = {rnaseq: "RNA-seq", methylation: "methylation", proteomics: "proteomics", metabolomics: "metabolomics"}; + const assay = (a) => ASSAY[a] ?? a; + const list = (xs) => xs.length < 2 ? xs.join("") : `${xs.slice(0, -1).join(", ")} and ${xs[xs.length - 1]}`; + // Lower-case a label mid-sentence, leaving acronyms like "pH" alone. + const lc = (s) => s.split(" ").map(w => w.length > 1 && w[1] === w[1].toLowerCase() ? w[0].toLowerCase() + w.slice(1) : w).join(" "); + const relevance = (id) => data.vocab.phenotype_resilience_relevance?.[id]; + + const real = data.studies.filter(s => !s.simulated); + const withEvidence = real.filter(s => s.n_evidence_records > 0); + const assaysRegistered = [...new Set(real.flatMap(s => s.assay_types))]; + const assaysWithEvidence = new Set(withEvidence.flatMap(s => s.assay_types)); + const awaiting = real.filter(s => s.n_evidence_records === 0 && s.qc_status !== "failed" + && s.assay_types.some(a => !assaysWithEvidence.has(a))); + const failed = real.filter(s => s.qc_status === "failed"); + + const pools = data.pools.filter(p => !p.simulated && p.feature_type === "gene"); + const cross = pools.filter(p => p.max_k >= 2).sort((a, b) => b.studies.length - a.studies.length); + const single = pools.filter(p => p.max_k < 2); + + const classes = data.candidates.by_evidence_class.filter(d => !d.simulated).sort((a, b) => b.n - a.n); + const nRanked = classes.reduce((s, d) => s + d.n, 0); + const top = data.candidates.top_real; + const topPhens = [...new Set(top.map(r => r.phenotype))]; + + const coverage = html`

${real.length} real studies cover ${k.n_real_stressor_classes} stressor + ${k.n_real_stressor_classes === 1 ? "class" : "classes"} + (${list(k.real_stressor_classes.map(s => lc(stress(s))))}) and + ${list(assaysRegistered.map(assay))} assays. ${withEvidence.length + ? `${withEvidence.length} of them have harmonized evidence (${k.n_comparisons_real_harmonized} of + ${k.n_comparisons_real} comparisons), ${assaysWithEvidence.size === 1 ? "all" : "spanning"} + ${list([...assaysWithEvidence].map(assay))}.` + : "None has harmonized evidence yet."}

`; + + const poolItems = cross.map(p => html`
  • ${phen(p.phenotype)} + (${p.studies.length} studies: ${p.studies.join(", ")}): ${fmt(p.n_features_k2plus)} genes measured in at + least two studies, ${fmt(p.n_features_k2plus_significant)} with a pooled effect significant after + multiple-testing correction.
  • `); + const poolBlock = cross.length + ? html`

    Cross-study pools, where the same gene can be compared between independent studies:

    ` + : html`

    No cross-study pool has formed yet; every phenotype is supported by a single study.

    `; + const singleText = single.length + ? html`

    ${list(single.map((p, i) => i ? lc(phen(p.phenotype)) : phen(p.phenotype)))} ${single.length === 1 ? "has" : "have"} evidence from one + study only (${list(single.flatMap(p => p.studies))}), so ${single.length === 1 ? "its" : "their"} genes stay + emerging until a second study replicates them.

    ` + : ""; + + const tierText = html`

    ${fmt(nRanked)} gene × phenotype candidates are ranked: + ${fmt(k.n_candidates_real_high_priority)} high priority (replicated across studies), + ${fmt(k.n_candidates_real_multi_omics)} multi-omics, ${fmt(k.n_candidates_real_emerging)} emerging. + By evidence class: ${classes.map(d => `${fmt(d.n)} ${d.evidence_class.replace(/_/g, " ")}`).join(", ")}.

    `; + + const caveats = []; + caveats.push(`${k.n_real_comparisons_with_resilience_phenotype} of ${k.n_comparisons_real} real comparisons + measure a resilience outcome; the rest measure disease, stress response, or exposure alone.`); + if (topPhens.length === 1 && relevance(topPhens[0]) !== "resilience") + caveats.push(`The top-ranked candidates all come from the ${lc(phen(topPhens[0]))} pool, so they describe the response to ${relevance(topPhens[0]) === "disease" ? "infection" : "the stressor"}, + not resistance to it; each pool's own top candidates are listed below.`); + if (k.n_candidates_real_multi_omics === 0) + caveats.push(`The multi-omics tier is empty because no non-RNA-seq study has been harmonized yet.`); + caveats.push(`A candidate here is an association across the available evidence, not a validated biomarker.`); + + const next = []; + if (awaiting.length) + next.push(`harmonize the registered ${list([...new Set(awaiting.flatMap(s => s.assay_types) + .filter(a => !assaysWithEvidence.has(a)))].map(assay))} studies (${awaiting.map(s => s.study_id).join(", ")}) + to open the multi-omics tier`); + if (single.length) + next.push(`find a second study for ${list(single.map(p => lc(phen(p.phenotype))))}`); + if (failed.length) + next.push(`resolve the QC failures in ${failed.map(s => s.study_id).join(", ")}`); + + return html`${coverage}${poolBlock}${singleText}${tierText} + ${Object.assign(document.createElement("p"), {innerHTML: caveats.join(" ")})} + ${next.length ? html`

    Next curation steps: ${list(next)} + (candidate studies).

    ` : ""}`; } ``` @@ -109,8 +171,8 @@ html`

    Ranked by tier, then score${topPhen ? html`, all for ${topPhen} ``` ```{ojs} -Inputs.table(data.candidates.top_real, { - columns: ["feature_id_standardized", "gene_description", ...(topPhen ? [] : ["phenotype_label"]), +candidateTable = (rows, {showPhenotype = false, maxRows = 15} = {}) => Inputs.table(rows, { + columns: ["feature_id_standardized", "gene_description", ...(showPhenotype ? ["phenotype_label"] : []), "tier", "k_studies", "pooled_effect", "adjusted_p_value"], header: {feature_id_standardized: "Gene", gene_description: "Description", phenotype_label: "Phenotype", tier: "Tier", k_studies: "Studies", pooled_effect: "log2FC", adjusted_p_value: "Adj. p"}, @@ -121,10 +183,33 @@ Inputs.table(data.candidates.top_real, { pooled_effect: v => (v == null ? "—" : v.toFixed(2)), adjusted_p_value: fmtP }, - rows: 15, layout: "auto" + rows: maxRows, layout: "auto" }) ``` +```{ojs} +candidateTable(data.candidates.top_real, {showPhenotype: !topPhen}) +``` + +### Top candidates per pool + +```{ojs} +{ + // The overall list above is dominated by the pool with the most studies, so each + // cross-study pool also shows its own best candidates. + const byPhen = data.candidates.top_real_by_phenotype ?? {}; + const pools = data.pools + .filter(p => !p.simulated && p.feature_type === "gene" && p.max_k >= 2 && byPhen[p.phenotype]?.length) + .sort((a, b) => b.studies.length - a.studies.length || a.phenotype.localeCompare(b.phenotype)); + if (!pools.length) return html`

    No cross-study pool has ranked candidates yet.

    `; + return html`${pools.map(p => html` +

    ${phen(p.phenotype)}

    +

    ${p.studies.length} studies (${p.studies.join(", ")}); top ${byPhen[p.phenotype].length} + by tier, then score.

    + ${candidateTable(byPhen[p.phenotype], {maxRows: byPhen[p.phenotype].length})}`)}`; +} +``` + ## Studies ```{ojs} diff --git a/src/reporting/dashboard.py b/src/reporting/dashboard.py index 12e63cb..2d48c20 100644 --- a/src/reporting/dashboard.py +++ b/src/reporting/dashboard.py @@ -45,6 +45,9 @@ def _built_crosswalks() -> list[Path]: SIGNIFICANCE_THRESHOLD = 0.05 TOP_N_CANDIDATES = 15 +# The overall top list is dominated by whichever pool has the most studies, so +# each phenotype also gets its own short list. +TOP_N_PER_PHENOTYPE = 5 # --------------------------------------------------------------------------- # @@ -266,7 +269,7 @@ def _gene_annotations(gene_ids: list[str]) -> dict[str, dict]: def _candidate_summary(phenotype_labels: dict, cards_dir: Path | None = None) -> dict: if not CANDIDATES_TABLE.exists(): return {"by_tier": [], "by_evidence_class": [], "top_real": [], "top_simulated": [], - "n_total": 0} + "top_real_by_phenotype": {}, "top_simulated_by_phenotype": {}, "n_total": 0} df = pd.read_csv(CANDIDATES_TABLE, sep="\t", low_memory=False) df["simulated"] = df["simulated"].map(_bool) by_tier = [ @@ -290,8 +293,8 @@ def _candidate_summary(phenotype_labels: dict, cards_dir: Path | None = None) -> "n_supporting_layers", "mapping_confidences", "simulated", "card_file"] cols = [c for c in cols if c in ranked.columns] - def _top(frame: pd.DataFrame) -> list[dict]: - rows = frame.head(TOP_N_CANDIDATES)[cols].to_dict(orient="records") + def _top(frame: pd.DataFrame, n: int = TOP_N_CANDIDATES) -> list[dict]: + rows = frame.head(n)[cols].to_dict(orient="records") for r in rows: r["feature_id_standardized"] = str(r["feature_id_standardized"]) r["phenotype_label"] = phenotype_labels.get(r.get("phenotype"), r.get("phenotype")) @@ -302,11 +305,18 @@ def _top(frame: pd.DataFrame) -> list[dict]: top_real = _top(ranked[~ranked["simulated"]]) top_sim = _top(ranked[ranked["simulated"]]) - for r in top_real + top_sim: + # Ranked order is kept within each phenotype, so head() is that phenotype's top. + by_phen_real = {phen: _top(grp, TOP_N_PER_PHENOTYPE) + for phen, grp in ranked[~ranked["simulated"]].groupby("phenotype", sort=True)} + by_phen_sim = {phen: _top(grp, TOP_N_PER_PHENOTYPE) + for phen, grp in ranked[ranked["simulated"]].groupby("phenotype", sort=True)} + all_real = top_real + [r for rows in by_phen_real.values() for r in rows] + all_rows = all_real + top_sim + [r for rows in by_phen_sim.values() for r in rows] + for r in all_rows: r["card_page"] = _publish_card(r.pop("card_file", None), cards_dir) - ann = _gene_annotations([r["feature_id_standardized"] for r in top_real + ann = _gene_annotations([r["feature_id_standardized"] for r in all_real if r.get("feature_type") == "gene"]) - for r in top_real: + for r in all_real: r.update(ann.get(r["feature_id_standardized"], {})) return { @@ -314,6 +324,8 @@ def _top(frame: pd.DataFrame) -> list[dict]: "by_evidence_class": by_class, "top_real": top_real, "top_simulated": top_sim, + "top_real_by_phenotype": by_phen_real, + "top_simulated_by_phenotype": by_phen_sim, "n_total": int(len(df)), "n_real": int((~df["simulated"]).sum()), "n_simulated": int(df["simulated"].sum()), diff --git a/tests/test_dashboard.py b/tests/test_dashboard.py index aae12a1..4be6222 100644 --- a/tests/test_dashboard.py +++ b/tests/test_dashboard.py @@ -17,6 +17,11 @@ } +def dashboard_tier_order(): + from prioritize.rank import TIER_ORDER + return TIER_ORDER + + def _point_at_reports(monkeypatch, reports_dir): monkeypatch.setattr(dashboard, "EVIDENCE_TABLE", reports_dir / "evidence" / "evidence_table.tsv") monkeypatch.setattr(dashboard, "CANDIDATES_TABLE", reports_dir / "evidence_cards" / "candidates.tsv") @@ -41,6 +46,7 @@ def test_builder_reports_missing_outputs_instead_of_failing(tmp_path, monkeypatc assert data["kpis"]["n_studies_real"] >= 1 assert data["kpis"]["n_evidence_records_real"] == 0 assert data["candidates"]["top_real"] == [] + assert data["candidates"]["top_real_by_phenotype"] == {} assert all(not c["harmonized"] for c in data["comparisons"]) # Nothing fabricated: no cards were copied when there were none to copy. assert not (out.parent / "cards").exists() @@ -70,6 +76,18 @@ def test_builder_summarizes_a_demo_pipeline_run(isolated_reports, isolated_regis assert kpis["n_candidates_real_high_priority"] == 0 assert all(r["simulated"] for r in data["candidates"]["top_simulated"]) assert data["candidates"]["top_real"] == [] + assert data["candidates"]["top_real_by_phenotype"] == {} + + # Each phenotype gets its own short list, so a pool that never reaches the + # overall top list still shows its best candidates, in ranked order. + by_phen = data["candidates"]["top_simulated_by_phenotype"] + assert {"thermal_tolerance", "larval_viability"} <= set(by_phen) + tier_rank = {t: i for i, t in enumerate(dashboard_tier_order())} + for phen, rows in by_phen.items(): + assert 0 < len(rows) <= dashboard.TOP_N_PER_PHENOTYPE + assert all(r["phenotype"] == phen and r["simulated"] for r in rows) + keys = [(tier_rank[r["tier"]], -r["score"]) for r in rows] + assert keys == sorted(keys) by_study = {r["study_id"]: r for r in data["evidence"]["by_study"]} assert set(ALL_DEMO_STUDY_IDS) <= set(by_study)