Add a Samples page: PSMs, peptides and proteins per sample, by condition

app.py

@@ -30,6 +30,7 @@
             st.Page(Path("content", "results_rescoring.py"), title="Rescoring", icon="📈"),
             st.Page(Path("content", "results_filtered.py"), title="Filtered PSMs", icon="🎯"),
             st.Page(Path("content", "results_abundance.py"), title="Protein Abundance", icon="📋"),
+            st.Page(Path("content", "results_samples.py"), title="Samples", icon="🧪"),
             st.Page(Path("content", "results_download.py"), title="Download", icon="⬇️"),
         ],
         # Steps run in order; each one works on the previous step's output.

content/results_samples.py

@@ -0,0 +1,142 @@
+"""Samples Page: PSMs, peptides and proteins per sample, grouped by condition."""
+from pathlib import Path
+
+import plotly.express as px
+import streamlit as st
+
+from src.common.common import page_setup, show_fig
+from src.common.insight_tables import show_insight_table
+from src.common.results_helpers import get_abundance_data, get_sample_group_map, get_workflow_dir
+from src.common.samples import COUNT_COLUMNS, channel_samples, condition_summary, lfq_samples
+from src.workflow.ParameterManager import ParameterManager
+
+params = page_setup()
+st.title("Samples")
+
+st.markdown(
+    """
+How many PSMs, peptides and proteins each sample contributed, grouped by the
+conditions set on the Configure page. A sample far below the others in its
+condition points to a failed run or a loading problem.
+"""
+)
+
+if "workspace" not in st.session_state:
+    st.warning("Please initialize your workspace first.")
+    st.stop()
+
+workflow_dir = get_workflow_dir(st.session_state["workspace"])
+results_dir = workflow_dir / "results"
+cache_dir = results_dir / "insight_cache"
+workflow_params = ParameterManager(workflow_dir, "TOPP Workflow").get_parameters_from_json()
+analysis_mode = workflow_params.get("analysis-mode", "LFQ")
+fdr = float(workflow_params.get("psm-fdr-percent", 1.0))
+
+
+def _mtimes(paths) -> tuple:
+    return tuple(Path(p).stat().st_mtime if Path(p).exists() else 0.0 for p in paths)
+
+
+@st.cache_data(show_spinner="Counting per sample...")
+def _lfq_samples(workflow_dir: str, files: tuple, groups: tuple, mtimes: tuple):
+    """lfq_samples, cached on the group assignments and the results' modification times."""
+    return lfq_samples(workflow_dir, list(files), dict(groups))
+
+
+if analysis_mode == "LFQ":
+    filter_dir = results_dir / "psm_filter"
+    if not filter_dir.is_dir() or not any(filter_dir.glob("*.idXML")):
+        st.info("No results available yet. Please run the workflow first.")
+        st.page_link("content/workflow_run.py", label="Go to Run", icon="🚀")
+        st.stop()
+    mzml_files = [str(f) for f in workflow_params.get("mzML-files", [])]
+    group_map = {
+        key[len("mzML-group-"):]: value
+        for key, value in workflow_params.items()
+        if key.startswith("mzML-group-") and value
+    }
+    sources = sorted(filter_dir.glob("*.idXML")) + sorted((results_dir / "quant_results").glob("*.csv"))
+    samples = _lfq_samples(
+        str(workflow_dir), tuple(mzml_files), tuple(sorted(group_map.items())), _mtimes(sources)
+    )
+    st.caption(
+        f"PSMs, peptides and proteins are the targets identified at {fdr:g} % PSM FDR, from the "
+        "top hit of each spectrum. Quantified proteins and peptides have an intensity above 0 in "
+        "the sample's ProteomicsLFQ result."
+    )
+else:
+    abundance = get_abundance_data(st.session_state["workspace"])
+    if abundance is None:
+        st.info("No quantification result available yet. Please run the workflow first.")
+        st.page_link("content/workflow_run.py", label="Go to Run", icon="🚀")
+        st.stop()
+    pivot_df, _, group_map = abundance
+    sample_group_map = get_sample_group_map(st.session_state["workspace"], pivot_df, group_map)
+    # The report's first columns are the protein and its description; the channels follow
+    sample_columns = [c for c in pivot_df.columns[4:] if c != "PeptideSequence"]
+    samples = channel_samples(pivot_df, sample_columns, sample_group_map)
+    sources = sorted((results_dir / "quant_results").glob("*.csv"))
+    st.caption(
+        "Each sample is an isobaric channel. Quantified proteins have an intensity above 0 in the "
+        "channel. Identifications are made per file, not per channel: see the Overview page."
+    )
+
+if samples.empty:
+    st.info("No samples in the last run.")
+    st.stop()
+
+def _table_height(rows: int) -> int:
+    """Room for every row up to a page of 12, plus the header and pager."""
+    return 140 + 34 * min(rows, 12)
+
+
+counts = [c for c in COUNT_COLUMNS if samples[c].notna().any()]
+group_key = "|".join(f"{s}={c}" for s, c in zip(samples["Sample"], samples["Condition"]))
+
+st.subheader("Per sample")
+metric = st.selectbox("Count", counts, key="samples-metric")
+fig = px.bar(
+    samples,
+    x="Sample",
+    y=metric,
+    color="Condition",
+    category_orders={"Sample": samples["Sample"].tolist()},
+)
+# One full-width bar per sample: px gives each condition its own offset group,
+# which would leave room beside every bar for the other conditions
+fig.update_traces(offsetgroup=None, alignmentgroup=None)
+fig.update_layout(xaxis_title=None, legend_title_text="Condition", bargap=0.2)
+# Slant the sample names once they no longer fit side by side
+label_chars = sum(len(s) for s in samples["Sample"])
+fig.update_xaxes(tickangle=-35 if label_chars > 90 else 0, automargin=True)
+show_fig(fig, f"samples_{metric.lower().replace(' ', '_')}")
+
+show_insight_table(
+    samples[["Sample", "Condition"] + counts],
+    name="samples_per_sample",
+    cache_dir=cache_dir,
+    source_files=sources,
+    extra_key=group_key,
+    title="samples",
+    height=_table_height(len(samples)),
+    layout="fitColumns",
+)
+
+st.subheader("Per condition")
+summary = condition_summary(samples)
+# Means go as text, since the table stores floats as float32 and would show 36.8 as 36.79999923706055
+for col in summary.columns[2:]:
+    summary[col] = summary[col].map(lambda v: "" if v != v else f"{v:g}")
+show_insight_table(
+    summary,
+    name="samples_per_condition",
+    cache_dir=cache_dir,
+    source_files=sources,
+    extra_key=group_key,
+    title="conditions",
+    height=_table_height(len(summary)),
+    layout="fitColumns",
+)
+if (samples["Condition"] == "Unassigned").any():
+    st.caption("Samples without a condition are shown as Unassigned. Set conditions on the Configure page.")
+    st.page_link("content/workflow_configure.py", label="Go to Configure", icon="⚙️")

src/common/samples.py

@@ -0,0 +1,147 @@
+"""Per-sample counts of a finished run: PSMs, peptides and proteins, by condition.
+
+Reads only what the workflow already wrote: the IDFilter idXMLs (identifications
+at the run's PSM FDR) and the quantification result. It makes no Streamlit
+calls, so the page caches it and tests call it.
+"""
+from pathlib import Path
+
+import pandas as pd
+
+from src.common.results_helpers import extract_filename_from_idxml, load_idxml
+
+UNASSIGNED = "Unassigned"
+
+SAMPLE_COLUMNS = [
+    "Sample",
+    "Condition",
+    "Target PSMs",
+    "Peptides",
+    "Proteins",
+    "Quantified proteins",
+    "Quantified peptides",
+]
+
+# Columns the bar chart can show, in the order the page offers them
+COUNT_COLUMNS = SAMPLE_COLUMNS[2:]
+
+
+def _is_decoy(hit) -> bool:
+    if not hit.metaValueExists("target_decoy"):
+        return False
+    value = hit.getMetaValue("target_decoy")
+    if isinstance(value, bytes):
+        value = value.decode()
+    return str(value).startswith("decoy")
+
+
+def identification_counts(idxml_path) -> dict:
+    """Target PSMs, distinct peptides and distinct proteins of one idXML.
+
+    Only the top hit of each spectrum counts, and decoy hits are left out.
+    Proteins are the accessions the target peptides map to.
+    """
+    _, peptides = load_idxml(idxml_path)
+    psms = 0
+    sequences, proteins = set(), set()
+    for pep in peptides:
+        hits = pep.getHits()
+        if not hits or _is_decoy(hits[0]):
+            continue
+        hit = hits[0]
+        psms += 1
+        sequences.add(hit.getSequence().toString())
+        proteins.update(ev.getProteinAccession() for ev in hit.getPeptideEvidences())
+    return {"Target PSMs": psms, "Peptides": len(sequences), "Proteins": len(proteins)}
+
+
+def lfq_quantified_counts(quant_csv) -> dict:
+    """Quantified proteins and peptides per input file name, from ProteomicsLFQ's CSV.
+
+    A protein (group) or peptide counts for a file when its intensity there is above 0.
+    """
+    try:
+        df = pd.read_csv(quant_csv, usecols=["ProteinName", "PeptideSequence", "Reference", "Intensity"])
+    except (OSError, ValueError):
+        return {}
+    df = df[pd.to_numeric(df["Intensity"], errors="coerce") > 0]
+    counts = df.groupby("Reference").agg(
+        proteins=("ProteinName", "nunique"), peptides=("PeptideSequence", "nunique")
+    )
+    return {
+        Path(str(ref)).name: {"Quantified proteins": int(row.proteins), "Quantified peptides": int(row.peptides)}
+        for ref, row in counts.iterrows()
+    }
+
+
+def _condition(group_map: dict, *names) -> str:
+    for name in names:
+        value = group_map.get(name)
+        if value is not None and str(value).strip():
+            return str(value).strip()
+    return UNASSIGNED
+
+
+def _order(df: pd.DataFrame) -> pd.DataFrame:
+    """Samples grouped by condition (unassigned last), in their original order within one."""
+    conditions = sorted(c for c in df["Condition"].unique() if c != UNASSIGNED)
+    if UNASSIGNED in set(df["Condition"]):
+        conditions.append(UNASSIGNED)
+    df = df.copy()
+    df["Condition"] = pd.Categorical(df["Condition"], categories=conditions, ordered=True)
+    df = df.sort_values("Condition", kind="stable").reset_index(drop=True)
+    df["Condition"] = df["Condition"].astype(str)
+    return df
+
+
+def lfq_samples(workflow_dir, mzml_files, group_map: dict) -> pd.DataFrame:
+    """One row per input file (SAMPLE_COLUMNS); files the run did not reach show empty cells.
+
+    ``group_map`` maps an input file name to its condition, as set on the Configure page.
+    """
+    results = Path(workflow_dir) / "results"
+    names = [Path(f).name for f in mzml_files]
+    filter_dir = results / "psm_filter"
+    idxmls = {}
+    if filter_dir.is_dir():
+        idxmls = {extract_filename_from_idxml(p, names): p for p in sorted(filter_dir.glob("*.idXML"))}
+    quant_csvs = sorted((results / "quant_results").glob("*.csv")) if (results / "quant_results").is_dir() else []
+    quantified = lfq_quantified_counts(quant_csvs[0]) if quant_csvs else {}
+
+    rows = []
+    for name in names:
+        ids = identification_counts(idxmls[name]) if name in idxmls else {}
+        rows.append({
+            "Sample": Path(name).stem,
+            "Condition": _condition(group_map, name, Path(name).stem),
+            **{col: ids.get(col) for col in ("Target PSMs", "Peptides", "Proteins")},
+            **{col: quantified.get(name, {}).get(col) for col in ("Quantified proteins", "Quantified peptides")},
+        })
+    return _order(pd.DataFrame(rows, columns=SAMPLE_COLUMNS))
+
+
+def channel_samples(pivot_df: pd.DataFrame, sample_columns, sample_group_map: dict) -> pd.DataFrame:
+    """One row per isobaric channel: the proteins quantified in it (intensity above 0).
+
+    Identifications are per file, not per channel, so those columns stay empty.
+    """
+    rows = []
+    for col in sample_columns:
+        values = pd.to_numeric(pivot_df[col], errors="coerce")
+        rows.append({
+            "Sample": col,
+            "Condition": _condition(sample_group_map, col),
+            "Quantified proteins": int((values > 0).sum()),
+        })
+    return _order(pd.DataFrame(rows, columns=SAMPLE_COLUMNS))
+
+
+def condition_summary(samples: pd.DataFrame) -> pd.DataFrame:
+    """Per condition: the number of samples and the mean of each count over them."""
+    counts = [c for c in COUNT_COLUMNS if samples[c].notna().any()]
+    numeric = samples[["Condition"] + counts].copy()
+    numeric[counts] = numeric[counts].apply(pd.to_numeric, errors="coerce")
+    summary = numeric.groupby("Condition", sort=False).agg(
+        Samples=("Condition", "size"), **{f"Mean {c}": (c, "mean") for c in counts}
+    )
+    return summary.round(1).reset_index()