From 4945fefcaa958a5c61a0ea6c2e0978c0984b4df2 Mon Sep 17 00:00:00 2001 From: Mariia Zueva Date: Mon, 31 Aug 2026 23:09:22 +0200 Subject: [PATCH] Pool replicate samples instead of refusing the run --- .changeset/pool-replicate-samples.md | 27 +++ model/src/columns.ts | 12 ++ model/src/index.ts | 39 +++- model/src/types.ts | 15 ++ software/src/errors.py | 7 +- software/src/main.py | 10 + software/src/pipeline.py | 16 +- software/src/pooling.py | 117 +++++++++++ software/src/validate.py | 56 ++---- software/tests/conftest.py | 18 +- software/tests/integration/test_cli.py | 107 ++++++++-- software/tests/unit/test_pooling.py | 258 +++++++++++++++++++++++++ software/tests/unit/test_validate.py | 67 +------ ui/src/components/RunStatistics.vue | 58 ++++++ ui/src/components/SettingsDrawer.vue | 3 +- 15 files changed, 667 insertions(+), 143 deletions(-) create mode 100644 .changeset/pool-replicate-samples.md create mode 100644 software/src/pooling.py create mode 100644 software/tests/unit/test_pooling.py diff --git a/.changeset/pool-replicate-samples.md b/.changeset/pool-replicate-samples.md new file mode 100644 index 0000000..5a231da --- /dev/null +++ b/.changeset/pool-replicate-samples.md @@ -0,0 +1,27 @@ +--- +'@platforma-open/milaboratories.sort-seq-analysis.software': patch +'@platforma-open/milaboratories.sort-seq-analysis.model': patch +'@platforma-open/milaboratories.sort-seq-analysis.ui': patch +'@platforma-open/milaboratories.sort-seq-analysis.block': patch +--- + +Pool replicate samples instead of refusing the run + +A condition-and-gate group holding more than one sample used to fail the run. Its reads are +now summed, matching what `titeseq-analysis` and `clonotype-enrichment` already do. + +The pooling is reported rather than silent, which is the objection the original refusal was +raised against: the manifest gains a `pooledGroups` list, the run log names every merged +group and its samples, and the Run statistics dialog shows a warning listing them. The dialog +resolves sample ids to sample labels through the sample axis's label column, so the alert +names samples the way the user does. Only retained conditions are reported — an excluded +condition is not part of the run. Where replicates supplied different sort fractions, the +non-null values are averaged and that is flagged separately. + +Two things broke quietly on replicate samples and are fixed at the source, by pooling before +anything reads the table: the read-distribution file emitted duplicate `(variantKey, gate)` +keys, and the sort-fraction check counted a twice-collected gate's fraction twice and refused +runs that were entirely valid. + +This is temporary. It does not decide whether replicates should be pooled at all rather than +scored separately, nor what to do when they disagree on the sort fraction. diff --git a/model/src/columns.ts b/model/src/columns.ts index 16f7e38..b559fe7 100644 --- a/model/src/columns.ts +++ b/model/src/columns.ts @@ -31,6 +31,8 @@ export const PColumnName = { MutationCount: "pl7.app/repertoire/mutationCount", /** The variant axis's label column — shown as "Variant Id". */ VariantLabel: "pl7.app/label", + /** The **sample** axis's label column. Same name as `VariantLabel`; the axis picks one. */ + SampleLabel: "pl7.app/label", /** The per-variant mutation list, shown as "Mutations". */ Mutations: "pl7.app/repertoire/mutations", } as const; @@ -104,6 +106,16 @@ export const metadataSelector: AnchoredPColumnSelector = { name: PColumnName.Metadata, }; +/** + * The sample axis's label column. `pl7.app/isLabel` is required alongside the name, which the + * variant label column shares. + */ +export const sampleLabelSelector: AnchoredPColumnSelector = { + axes: [{ anchor: "main", idx: 0 }], + name: PColumnName.SampleLabel, + annotations: { "pl7.app/isLabel": "true" }, +}; + /** * The variants' mutation count. Resolved in the workflow rather than here — this constant * exists so the model and the workflow state the same predicate, and so a reader can see diff --git a/model/src/index.ts b/model/src/index.ts index 9f74e0f..a1e82bc 100644 --- a/model/src/index.ts +++ b/model/src/index.ts @@ -12,7 +12,14 @@ import { type SUniversalPColumnId, } from "@platforma-sdk/model"; import { kind } from "@platforma-open/milaboratories.sort-seq-analysis.kind"; -import { Annotation, FacsBin, isAbundanceAnchor, metadataSelector, PColumnName } from "./columns"; +import { + Annotation, + FacsBin, + isAbundanceAnchor, + metadataSelector, + PColumnName, + sampleLabelSelector, +} from "./columns"; import { blockDataModel } from "./dataModel"; import type { BlockArgs, BlockData, RunManifest } from "./types"; @@ -24,10 +31,9 @@ export * from "./types"; * Every configuration rule, checked here and **nowhere else**. These six are decidable from the * arguments and snapshotted column values alone, so they are refused before the run starts. * - * The two data-value rules — sort fractions, one sample per condition-and-gate group — belong - * to the computation and are deliberately not approximated here: a duplicated rule is one that - * will disagree, and it fails by drifting looser, so the settings pass and the run fails anyway - * with a different message. + * The one data-value rule — sort fractions — belongs to the computation and is deliberately not + * approximated here: a duplicated rule is one that will disagree, and it fails by drifting + * looser, so the settings pass and the run fails anyway with a different message. */ export function settingsIssues(data: BlockData): string[] { const issues: string[] = []; @@ -233,6 +239,29 @@ export const platforma = BlockModelV3.create({ dataModel: blockDataModel, kind } return ctx.createPFrame(columns as PColumn[]); }) + /** + * The sample label column, so the UI can resolve the pooling report's `PlId`s to sample names. + * Kept out of `metadataColumnsPframe`, whose columns are the roles the pickers offer. + */ + .output("sampleLabelPframe", (ctx) => { + const anchor = ctx.data.abundanceRef; + if (!anchor) return undefined; + + const columns = ctx.resultPool.getAnchoredPColumns({ main: anchor }, [sampleLabelSelector]); + if (!columns || columns.length === 0) return undefined; + + return ctx.createPFrame(columns as PColumn[]); + }) + + /** The `PObjectId` for reading that column out of the frame above. */ + .output("sampleLabelColumnId", (ctx) => { + const anchor = ctx.data.abundanceRef; + if (!anchor) return undefined; + + const columns = ctx.resultPool.getAnchoredPColumns({ main: anchor }, [sampleLabelSelector]); + return columns?.[0]?.id; + }) + // --------------------------------------------------------------------------- // The three views. // --------------------------------------------------------------------------- diff --git a/model/src/types.ts b/model/src/types.ts index 4ae7ece..f914c4e 100644 --- a/model/src/types.ts +++ b/model/src/types.ts @@ -105,6 +105,19 @@ export type GateCollected = { depth: number; }; +/** One condition-and-gate group whose reads were pooled across several samples. */ +export type PooledGroup = { + condition: string; + gate: string; + /** Sample **ids**, sorted. The UI resolves them to labels via `sampleLabelPframe`. */ + samples: string[]; + /** + * The replicates supplied different fractions and the non-null values were averaged. Absent + * in the uncorrected mode. + */ + sortFractionsDiffer?: boolean; +}; + /** The manifest's per-condition entry. */ export type ConditionSummary = { /** Verbatim, exactly as it appears in the metadata column. */ @@ -141,5 +154,7 @@ export type RunManifest = { * no mutation-count table at all, which is a third state rather than a reason. */ parentAbsenceReason: string | null; + /** Retained conditions only; empty on a run with no replicates. */ + pooledGroups: PooledGroup[]; conditions: ConditionSummary[]; }; diff --git a/software/src/errors.py b/software/src/errors.py index 47a1c35..ec81ef0 100644 --- a/software/src/errors.py +++ b/software/src/errors.py @@ -6,10 +6,9 @@ condition excluded, a negative floor — are refused by the block model before the run starts and never arrive here. -What arrives here is the data-value class: a sort-fraction set breaking -`sort-fraction-values`, and more than one sample in a condition-and-gate group. Those -raise `Refusal`, which `main` turns into a non-zero exit that names the offending -values and writes no file. +What arrives here is the data-value class, which is a single rule: a sort-fraction set +breaking `sort-fraction-values`. It raises `Refusal`, which `main` turns into a non-zero +exit that names the offending values and writes no file. A `Refusal` is therefore *not* an internal assertion. Anything raised as some other exception type is a bug in this package or a caller violating the interface, and is diff --git a/software/src/main.py b/software/src/main.py index 7542b1b..f4c412b 100644 --- a/software/src/main.py +++ b/software/src/main.py @@ -90,6 +90,16 @@ def _report(manifest: dict) -> None: detail = reason if reason is not None else "no mutation-count table supplied" print(f"Parent row: not identified ({detail})") + for entry in manifest["pooledGroups"]: + samples = " + ".join(entry["samples"]) + line = ( + f"Pooled condition {entry['condition']!r} gate {entry['gate']!r}: " + f"{len(entry['samples'])} samples merged ({samples})" + ) + if entry.get("sortFractionsDiffer"): + line += " — WARNING: replicates supplied different sort fractions; averaged" + print(line) + for entry in manifest["conditions"]: gates = ", ".join(f"{gate['gate']}={gate['depth']}" for gate in entry["gatesCollected"]) summary = ( diff --git a/software/src/pipeline.py b/software/src/pipeline.py index a459ba1..7da2e71 100644 --- a/software/src/pipeline.py +++ b/software/src/pipeline.py @@ -32,18 +32,22 @@ ) from io_layer import distribution_file_name, score_file_name, write_manifest, write_table from params import Params -from validate import check_one_sample_per_group, check_sort_fractions +from pooling import pool_replicates +from validate import check_sort_fractions def run(reads: pl.DataFrame, variants: pl.DataFrame | None, params: Params, out_dir: Path) -> dict: """Score every retained condition and write every file. Returns the manifest.""" in_scope = selected_gates(reads, params) + + # Before pooling: the report is confined to retained conditions. retained = retained_conditions(in_scope, params) - # Both refusals run over the whole run before anything is written, so a failure leaves - # nothing partial behind. Over the in-scope rows only: a second sample in a gate the run - # does not cover, or a missing sort fraction on one, is not this run's problem. - check_one_sample_per_group(in_scope, retained) + # Must run before anything else reads the table — see `pooling`. Over the in-scope rows + # only, so a second sample in a gate the run does not cover is not pooled. + in_scope, pooled_groups = pool_replicates(in_scope, params.sort_fraction_column, retained) + + # Over the whole run before anything is written, so a failure leaves nothing partial. if params.sort_fraction_column is not None: check_sort_fractions(in_scope, params.sort_fraction_column, retained) @@ -61,6 +65,8 @@ def run(reads: pl.DataFrame, variants: pl.DataFrame | None, params: Params, out_ # Null where there was no mutation-count table at all: the two reasons below are # for a table that exists and does not identify a single parent. "parentAbsenceReason": parent.absence_reason, + # Retained conditions only; empty on a run with no replicates. + "pooledGroups": pooled_groups, "conditions": conditions, } write_manifest(manifest, out_dir) diff --git a/software/src/pooling.py b/software/src/pooling.py new file mode 100644 index 0000000..5ca65ff --- /dev/null +++ b/software/src/pooling.py @@ -0,0 +1,117 @@ +"""Pool replicate samples that share a condition-and-gate group. + +Pooling here, once, before anything else reads the table keeps the rest of the package on its +one-row-per-(gate, variant) grain. Two places break silently otherwise: +`scoring.read_distribution` emits duplicate (variantKey, gate) keys, and +`check_sort_fractions` counts a twice-collected gate's fraction twice and over-sums. + +Temporary. Whether replicates should be pooled at all rather than scored separately is open, +as is what to do when they disagree on the sort fraction (here: averaged and flagged). +""" + +from __future__ import annotations + +import polars as pl +from constants import ( + COL_CONDITION, + COL_GATE, + COL_READS, + COL_SAMPLE, + COL_VARIANT, +) + +# A single-sample group keeps its bare sample id, so pooling is a no-op there. +_LABEL_SEPARATOR = "+" + +_SAMPLES = "_samples" +_FRACTION_VALUES = "_fractionValues" + + +def pool_replicates( + reads: pl.DataFrame, sort_fraction_column: str | None, retained_conditions: list[str] +) -> tuple[pl.DataFrame, list[dict]]: + """Collapse to one row per (condition, gate, variant): reads summed, sort fraction the mean + of the group's non-null values, other metadata the first in sample-id order, sample ids + joined. Returns the pooled table and the groups that held more than one sample. + + Metadata resolves per (condition, gate), not per variant: a variant missing from one + replicate would give its gate a shorter label, and `check_sort_fractions` groups by sample + and would then see that gate twice. + + A `sort_fraction_column` the table does not carry is left for `check_sort_fractions` to + refuse, rather than pre-empted here with a polars traceback. + + Only retained conditions are reported. Excluded rows are still pooled; nothing reads them. + """ + fraction = sort_fraction_column if sort_fraction_column in reads.columns else None + + # Read off the frame rather than listed: the workflow decides which metadata columns the + # reads table carries, and one named here but not there would vanish from the run. + fixed = {COL_CONDITION, COL_GATE, COL_VARIANT, COL_SAMPLE, COL_READS} + passthrough = [name for name in reads.columns if name not in fixed and name != fraction] + + # Deduplicated to one row per sample first: metadata repeats across the sample's variant + # rows, so aggregating the raw rows would weight a sample by how many variants it detected. + metadata_columns = [COL_CONDITION, COL_GATE, COL_SAMPLE, *passthrough] + if fraction is not None: + metadata_columns.append(fraction) + + by_sample = pl.col(COL_SAMPLE).sort() + meta_aggs = [pl.col(COL_SAMPLE).unique().sort().alias(_SAMPLES)] + if fraction is not None: + meta_aggs.append(pl.col(fraction).alias(_FRACTION_VALUES)) + # Sample-id order, so the pick is deterministic rather than row-order dependent. + meta_aggs += [pl.col(name).sort_by(by_sample).first().alias(name) for name in passthrough] + + group_meta = ( + reads.select(metadata_columns) + .unique() + .group_by(COL_CONDITION, COL_GATE) + .agg(*meta_aggs) + .sort(COL_CONDITION, COL_GATE) + ) + + report = _report(group_meta, fraction, retained_conditions) + + resolved = [pl.col(_SAMPLES).list.join(_LABEL_SEPARATOR).alias(COL_SAMPLE), *passthrough] + if fraction is not None: + resolved.append(pl.col(_FRACTION_VALUES).list.drop_nulls().list.mean().alias(fraction)) + resolved_meta = group_meta.select(COL_CONDITION, COL_GATE, *resolved) + + pooled = ( + reads.group_by(COL_CONDITION, COL_GATE, COL_VARIANT) + .agg(pl.col(COL_READS).sum().alias(COL_READS)) + .join(resolved_meta, on=[COL_CONDITION, COL_GATE], how="left") + .select(reads.columns) + .sort(COL_CONDITION, COL_GATE, COL_VARIANT) + ) + return pooled, report + + +def _report( + group_meta: pl.DataFrame, sort_fraction_column: str | None, retained_conditions: list[str] +) -> list[dict]: + """The replicated groups of retained conditions, in condition-and-gate order. + + `sortFractionsDiffer` means the group's replicates supplied different fractions, which is + either two separate sorts or wrong metadata — the fraction is a property of the gate. + """ + replicated = group_meta.filter( + (pl.col(_SAMPLES).list.len() > 1) + & pl.col(COL_CONDITION).is_in(retained_conditions) + ) + if replicated.height == 0: + return [] + + entries = [] + for row in replicated.iter_rows(named=True): + entry = { + "condition": row[COL_CONDITION], + "gate": row[COL_GATE], + "samples": list(row[_SAMPLES]), + } + if sort_fraction_column is not None: + supplied = [value for value in row[_FRACTION_VALUES] if value is not None] + entry["sortFractionsDiffer"] = len(set(supplied)) > 1 + entries.append(entry) + return entries diff --git a/software/src/validate.py b/software/src/validate.py index c69b1a6..eb0b41b 100644 --- a/software/src/validate.py +++ b/software/src/validate.py @@ -1,8 +1,8 @@ -"""The two data-value refusals this side owns. +"""The one data-value refusal this side owns. -`validation-boundary` assigns exactly these to the computation, because both need the -column *values* rather than a picked column reference: whether a set of fractions sums to -1, and whether two samples share a gate, are properties of the project's data. +`validation-boundary` assigns it to the computation because it needs the column *values* +rather than a picked column reference: whether a set of fractions sums to 1 is a property +of the project's data. Every configuration rule — a required argument absent, an anchor resolving to nothing or to more than one column, the three metadata roles not distinct, an empty gate order, @@ -12,7 +12,11 @@ specific and bad, a model-side check drifting looser so the settings pass and the run fails anyway with a different message. -Both refusals name the offending values and are raised before any file is written, so +`pooling` has already collapsed the table to one row per (condition, gate, variant), which +is what keeps `_check_sum_per_condition` correct: it sums one value per sample, and a gate +collected twice would otherwise contribute its fraction twice. + +The refusal names the offending values and is raised before any file is written, so nothing partial is produced. """ @@ -31,45 +35,6 @@ from errors import Refusal -def check_one_sample_per_group(reads: pl.DataFrame, retained_conditions: list[str]) -> None: - """v1 admits at most one sample per condition-and-gate group. - - Two samples in a group **fail the run**. Their reads are not pooled and neither sample - is preferred: summing them would move every depth, frequency and weighted mean in the - run, which is exactly the shape of failure `input-defaults` bars — output of ordinary - length and plausible content, with nothing saying an aggregation was chosen. - - Replicate support is a v2 addition, and what v2 must supply is an aggregation rule for - the reads plus an agreement rule for the per-sample values the group then has more - than one of. Until it does, refusing is the only honest option. - - No sample for a pair is **not** an error — that gate was not collected at that - condition, which clauses 1 and 2 already accommodate. - - Only retained conditions are checked; an excluded value is not part of the run. The - caller likewise hands over only the selected gates' rows, so two samples sharing a gate - the run does not cover pass unremarked — that pair is not a group of this run. - """ - offenders = ( - reads.filter(pl.col(COL_CONDITION).is_in(retained_conditions)) - .group_by(COL_CONDITION, COL_GATE) - .agg(pl.col(COL_SAMPLE).unique().sort().alias("samples")) - .filter(pl.col("samples").list.len() > 1) - .sort(COL_CONDITION, COL_GATE) - ) - if offenders.height == 0: - return - - detail = "; ".join( - f"condition {row[COL_CONDITION]!r} gate {row[COL_GATE]!r}: samples {', '.join(row['samples'])}" - for row in offenders.iter_rows(named=True) - ) - raise Refusal( - f"more than one sample in a condition-and-gate group, which v1 does not support " - f"(reads are not pooled): {detail}" - ) - - def check_sort_fractions(reads: pl.DataFrame, column: str, retained_conditions: list[str]) -> None: """The three requirements of `sort-fraction-values`. @@ -104,7 +69,8 @@ def check_sort_fractions(reads: pl.DataFrame, column: str, retained_conditions: in_run = reads.filter(pl.col(COL_CONDITION).is_in(retained_conditions)) # One supplied value per sample: the sort fraction is per-sample metadata, repeated - # across the sample's variant rows. + # across the sample's variant rows. After `pooling` this is also one value per + # (condition, gate), which is the grain the sum below must be taken on. per_sample = in_run.group_by(COL_SAMPLE, COL_CONDITION, COL_GATE).agg( pl.col(column).first().alias("fraction"), pl.col(COL_READS).sum().alias("sample_reads"), diff --git a/software/tests/conftest.py b/software/tests/conftest.py index d2feecc..a3b01c4 100644 --- a/software/tests/conftest.py +++ b/software/tests/conftest.py @@ -26,8 +26,9 @@ def reads_frame(rows: list[tuple[str, str, int]], condition: str = "pH7", fractions: dict[str, float] | None = None): """Build a reads table from (gate, variantKey, reads) triples. - One sample per (condition, gate) — the v1 grain `sample-gate-grouping` requires — with - the sample id derived from the pair so a caller cannot accidentally collide two. + One sample per (condition, gate) — the grain the arithmetic is defined on — with the + sample id derived from the pair so a caller cannot accidentally collide two. Use + `replicate_frame` for the second sample of a gate. """ frame = pl.DataFrame( { @@ -46,6 +47,19 @@ def reads_frame(rows: list[tuple[str, str, int]], condition: str = "pH7", fracti return frame +def replicate_frame( + rows: list[tuple[str, str, int]], + sample: str, + condition: str = "pH7", + fractions: dict[str, float] | None = None, +): + """A second sample for gates `reads_frame` already covered. The sample id is given rather + than derived, so it collides on (condition, gate); concat the two frames.""" + return reads_frame(rows, condition=condition, fractions=fractions).with_columns( + pl.lit(sample).alias("sampleId") + ) + + # The base table. Every gate's depth is exactly 100. # # gate | P | A | B | C | depth diff --git a/software/tests/integration/test_cli.py b/software/tests/integration/test_cli.py index ac0de5b..90471bd 100644 --- a/software/tests/integration/test_cli.py +++ b/software/tests/integration/test_cli.py @@ -15,6 +15,7 @@ BASE_MUTATION_COUNTS, BASE_ROWS, reads_frame, + replicate_frame, variants_frame, write_params, write_tsv, @@ -412,30 +413,100 @@ def test_over_summing_fractions_exits_non_zero_and_writes_nothing(tmp_path, caps assert "REFUSED" in captured.err -def test_replicate_samples_exit_non_zero_and_write_nothing(tmp_path, capsys): +def test_a_run_with_no_replicates_reports_no_pooling(tmp_path): + """The empty list is the statement: one sample per gate.""" + _, _, manifest = invoke(tmp_path, reads_frame(BASE_ROWS), variants_frame(BASE_MUTATION_COUNTS)) + + assert manifest["pooledGroups"] == [] + + +def test_replicate_samples_are_pooled_and_the_pooling_is_reported(tmp_path, capsys): + """g1 is replicated with its own read counts, so its depth doubles (100 -> 200) along with + every read count in it — leaving every frequency and every weighted mean unchanged. + """ + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + reads = pl.concat([reads_frame(BASE_ROWS), replicate_frame(g1_rows, "replicate_of_g1")]) + + code, out_dir, manifest = invoke(tmp_path, reads, variants_frame(BASE_MUTATION_COUNTS)) + + assert code == 0 + assert manifest["pooledGroups"] == [ + {"condition": "pH7", "gate": "g1", "samples": ["replicate_of_g1", "s_pH7_g1"]} + ] + + entry = manifest["conditions"][0] + # The pooled depth: the number the arithmetic used, which is what a read floor is set against. + depths = {gate["gate"]: gate["depth"] for gate in entry["gatesCollected"]} + assert depths == {"g1": 200, "g2": 100, "g3": 100} + assert read_scores(out_dir, entry["gateRankMeanFile"], "gateRankMean") == pytest.approx(BASE_MEANS, rel=REL) + + captured = capsys.readouterr() + assert "Pooled condition 'pH7' gate 'g1'" in captured.out + assert "replicate_of_g1" in captured.out + + +def test_replicates_disagreeing_on_the_sort_fraction_warn_rather_than_refuse(tmp_path, capsys): + """0.5 and 0.3 average to 0.4, so the condition sums to 0.9 and the run proceeds.""" + fractions = {"g1": 0.5, "g2": 0.3, "g3": 0.2} + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + reads = pl.concat( + [ + reads_frame(BASE_ROWS, fractions=fractions), + replicate_frame(g1_rows, "replicate_of_g1", fractions={**fractions, "g1": 0.3}), + ] + ) + + code, _, manifest = invoke( + tmp_path, reads, variants_frame(BASE_MUTATION_COUNTS), sort_fraction_column="sortFraction" + ) + + assert code == 0 + assert manifest["pooledGroups"][0]["sortFractionsDiffer"] is True + assert manifest["conditions"][0]["sortFractionSum"] == pytest.approx(0.9, rel=REL) + assert "different sort fractions" in capsys.readouterr().out + + +def test_excluded_conditions_are_not_reported_as_pooled(tmp_path): + """A factor column with three values, one selected: only that one's groups are the run's.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + reads = pl.concat( + [ + reads_frame(BASE_ROWS, condition=condition) + for condition in ("specificity", "affinity", "polyspecificity") + ] + + [ + replicate_frame(g1_rows, f"rep_{condition}", condition=condition) + for condition in ("specificity", "affinity", "polyspecificity") + ] + ) + + code, _, manifest = invoke( + tmp_path, + reads, + variants_frame(BASE_MUTATION_COUNTS), + excluded=["affinity", "polyspecificity"], + ) + + assert code == 0 + assert [entry["condition"] for entry in manifest["conditions"]] == ["specificity"] + assert [entry["condition"] for entry in manifest["pooledGroups"]] == ["specificity"] + + +def test_a_replicate_in_an_unselected_gate_is_not_pooled(tmp_path): + """A gate the run does not rank is dropped whole, replicate included, rather than merged.""" reads = pl.concat( [ reads_frame(BASE_ROWS), - pl.DataFrame( - { - "sampleId": ["replicate"], - "variantKey": ["P"], - "reads": [5], - "condition": ["pH7"], - "gate": ["g1"], - }, - schema_overrides={"reads": pl.Int64}, - ), + reads_frame([("unsorted", "P", 40)]), + replicate_frame([("unsorted", "P", 40)], "replicate_of_unsorted"), ] ) - code, out_dir, manifest = invoke(tmp_path, reads, variants_frame(BASE_MUTATION_COUNTS)) - assert code == 1 - assert manifest is None - assert not out_dir.exists() - captured = capsys.readouterr() - assert "replicate" in captured.out - assert "replicate" in captured.err + code, _, manifest = invoke(tmp_path, reads, variants_frame(BASE_MUTATION_COUNTS)) + + assert code == 0 + assert manifest["pooledGroups"] == [] + assert [gate["gate"] for gate in manifest["conditions"][0]["gatesCollected"]] == ["g1", "g2", "g3"] def test_sort_fraction_column_missing_from_the_reads_table_fails(tmp_path, capsys): diff --git a/software/tests/unit/test_pooling.py b/software/tests/unit/test_pooling.py new file mode 100644 index 0000000..8bf4725 --- /dev/null +++ b/software/tests/unit/test_pooling.py @@ -0,0 +1,258 @@ +"""Pooling replicate samples that share a condition-and-gate group. + +Every case asserts a hand-computed number or an expected absence, per +`computation-test-suite`. +""" + +from __future__ import annotations + +import polars as pl +import pytest +from conftest import BASE_MEANS, BASE_ROWS, GATE_RANKS, means_as_dict, reads_frame, replicate_frame + +import scoring +import validate +from pooling import pool_replicates + +FRACTIONS = {"g1": 0.5, "g2": 0.3, "g3": 0.2} +RETAINED = ["pH7"] + + +# --------------------------------------------------------------------------- +# The no-replicate case: pooling must not move anything. +# --------------------------------------------------------------------------- + + +def test_a_run_with_no_replicates_is_unchanged(): + """Only the row order differs, which nothing reads — every emitted file sorts by its own + keys.""" + frame = reads_frame(BASE_ROWS, fractions=FRACTIONS) + pooled, report = pool_replicates(frame, "sortFraction", RETAINED) + + assert report == [] + key = ["condition", "gate", "variantKey"] + assert pooled.sort(key).equals(frame.select(pooled.columns).sort(key)) + + +# --------------------------------------------------------------------------- +# Pooling the reads. +# --------------------------------------------------------------------------- + + +def test_pooling_an_exact_duplicate_of_a_gate_leaves_every_score_unchanged(): + """Replicating g1 with its own read counts doubles its depth (100 -> 200) and every + variant's reads in it, so every frequency in g1 — and every weighted mean — is unchanged. + Pooling merges depth without tilting the profile. + """ + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat([reads_frame(BASE_ROWS), replicate_frame(g1_rows, "rep_g1")]) + + pooled, _ = pool_replicates(frame, None, RETAINED) + means = scoring.gate_rank_means(scoring.per_gate_frequencies(pooled, None), GATE_RANKS) + + assert means_as_dict(means) == pytest.approx(BASE_MEANS, rel=1e-12) + + +def test_pooled_reads_are_summed_and_the_weighted_mean_follows(): + """Two gates, two samples in g1, hand-computed end to end. + + g1 base P 20 A 30 replicate P 30 A 20 pooled P 50 A 50 depth 100 + g2 P 80 A 20 depth 100 + + freq P: g1 .5 g2 .8 A: g1 .5 g2 .2 + P: den .5 + .8 = 1.3 num 1(.5) + 2(.8) = 2.1 mean 21/13 + A: den .5 + .2 = .7 num 1(.5) + 2(.2) = .9 mean 9/7 + """ + ranks = {"g1": 1, "g2": 2} + base = reads_frame([("g1", "P", 20), ("g1", "A", 30), ("g2", "P", 80), ("g2", "A", 20)]) + replicate = replicate_frame([("g1", "P", 30), ("g1", "A", 20)], "rep_g1") + + pooled, report = pool_replicates(pl.concat([base, replicate]), None, RETAINED) + + assert pooled.filter((pl.col("gate") == "g1") & (pl.col("variantKey") == "P")).item(0, "reads") == 50 + assert pooled.filter((pl.col("gate") == "g1") & (pl.col("variantKey") == "A")).item(0, "reads") == 50 + assert [entry["gate"] for entry in report] == ["g1"] + + means = means_as_dict(scoring.gate_rank_means(scoring.per_gate_frequencies(pooled, None), ranks)) + assert means["P"] == 21 / 13 + assert means["A"] == 9 / 7 + + +def test_a_gate_collected_once_keeps_one_row_per_variant(): + """One row per (condition, gate, variant), replicated gate or not — the grain + `read_distribution` and `check_sort_fractions` rely on.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat([reads_frame(BASE_ROWS), replicate_frame(g1_rows, "rep_g1")]) + + pooled, _ = pool_replicates(frame, None, RETAINED) + + duplicates = pooled.group_by("condition", "gate", "variantKey").len().filter(pl.col("len") > 1) + assert duplicates.height == 0 + + +def test_conditions_are_pooled_independently(): + """A sample id repeated across conditions is two groups, not one.""" + frame = pl.concat( + [ + reads_frame([("g1", "P", 10)], condition="pH7"), + reads_frame([("g1", "P", 40)], condition="pH5"), + ] + ) + pooled, report = pool_replicates(frame, None, ["pH5", "pH7"]) + + assert report == [] + reads = dict(zip(pooled["condition"].to_list(), pooled["reads"].to_list(), strict=True)) + assert reads == {"pH7": 10, "pH5": 40} + + +# --------------------------------------------------------------------------- +# The pooled label. +# --------------------------------------------------------------------------- + + +def test_the_pooled_label_joins_the_sample_ids_and_is_the_same_on_every_variant(): + """One label per (condition, gate), not per variant: a variant missing from one replicate + would otherwise give its gate a shorter label, and `check_sort_fractions` groups by + sample.""" + base = reads_frame([("g1", "P", 10), ("g1", "A", 20)]) + # The replicate detected P and not A, which is the case that splits the label. + replicate = replicate_frame([("g1", "P", 5)], "rep_g1") + + pooled, _ = pool_replicates(pl.concat([base, replicate]), None, RETAINED) + + labels = set(pooled["sampleId"].to_list()) + assert labels == {"rep_g1+s_pH7_g1"} + + +# --------------------------------------------------------------------------- +# The sort fraction. This is where dropping the refusal alone goes wrong. +# --------------------------------------------------------------------------- + + +def test_a_valid_two_replicate_run_is_no_longer_refused_for_over_summing(): + """`check_sort_fractions` sums one value per sample, so g1 collected twice would contribute + 0.5 twice and over-sum to 1.5. Pooled, g1 supplies one value and the sum is 1.0.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat( + [ + reads_frame(BASE_ROWS, fractions=FRACTIONS), + replicate_frame(g1_rows, "rep_g1", fractions=FRACTIONS), + ] + ) + + pooled, report = pool_replicates(frame, "sortFraction", RETAINED) + + assert report[0]["sortFractionsDiffer"] is False + fractions = pooled.group_by("gate").agg(pl.col("sortFraction").first()) + assert fractions["sortFraction"].sum() == 1.0 + # No refusal. + validate.check_sort_fractions(pooled, "sortFraction", ["pH7"]) + + +def test_replicates_disagreeing_on_the_fraction_are_averaged_and_flagged(): + """0.5 and 0.3 average to 0.4, and the disagreement is reported rather than resolved.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat( + [ + reads_frame(BASE_ROWS, fractions=FRACTIONS), + replicate_frame(g1_rows, "rep_g1", fractions={**FRACTIONS, "g1": 0.3}), + ] + ) + + pooled, report = pool_replicates(frame, "sortFraction", RETAINED) + + assert report[0]["gate"] == "g1" + assert report[0]["sortFractionsDiffer"] is True + g1 = pooled.filter(pl.col("gate") == "g1") + assert g1["sortFraction"].unique().to_list() == [0.4] + + +def test_a_null_fraction_on_one_replicate_does_not_hide_a_value_on_the_other(): + """The mean is over the non-null values, so a group that supplied a fraction keeps it. + All-null still averages to null, so the presence refusal can still fire.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + replicate = replicate_frame(g1_rows, "rep_g1", fractions=FRACTIONS).with_columns( + pl.lit(None, dtype=pl.Float64).alias("sortFraction") + ) + frame = pl.concat([reads_frame(BASE_ROWS, fractions=FRACTIONS), replicate]) + + pooled, _ = pool_replicates(frame, "sortFraction", RETAINED) + + assert pooled.filter(pl.col("gate") == "g1")["sortFraction"].unique().to_list() == [0.5] + + +def test_a_group_null_on_every_replicate_stays_null(): + frame = pl.concat( + [ + reads_frame([("g1", "P", 10)], fractions=FRACTIONS), + replicate_frame([("g1", "P", 10)], "rep_g1", fractions=FRACTIONS), + ] + ).with_columns(pl.lit(None, dtype=pl.Float64).alias("sortFraction")) + + pooled, _ = pool_replicates(frame, "sortFraction", RETAINED) + + assert pooled["sortFraction"].to_list() == [None] + + +# --------------------------------------------------------------------------- +# The distribution file. The second thing dropping the refusal alone breaks. +# --------------------------------------------------------------------------- + + +def test_the_distribution_carries_one_row_per_variant_and_gate(): + """`read_distribution` joins `per_gate` onto a variant x gate grid, so unpooled a + replicated gate matches twice. The workflow imports this file as a p-column keyed on + exactly those two axes.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat([reads_frame(BASE_ROWS), replicate_frame(g1_rows, "rep_g1")]) + + pooled, _ = pool_replicates(frame, None, RETAINED) + per_gate = scoring.per_gate_frequencies(pooled, None) + scored = scoring.gate_rank_means(per_gate, GATE_RANKS) + distribution = scoring.read_distribution(per_gate, scored, GATE_RANKS, 20) + + duplicates = distribution.group_by("variantKey", "gate").len().filter(pl.col("len") > 1) + assert duplicates.height == 0 + # g1's pooled reads for P: 10 + 10. + p_g1 = distribution.filter((pl.col("variantKey") == "P") & (pl.col("gate") == "g1")) + assert p_g1.item(0, "gateReads") == 20 + + +# --------------------------------------------------------------------------- +# The report is confined to the run. +# --------------------------------------------------------------------------- + + +def test_a_replicated_group_in_an_excluded_condition_is_not_reported(): + """An excluded condition is not part of the run. Its rows are still pooled; nothing reads + them.""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat( + [ + reads_frame(BASE_ROWS, condition="specificity"), + replicate_frame(g1_rows, "rep_specificity", condition="specificity"), + reads_frame(BASE_ROWS, condition="affinity"), + replicate_frame(g1_rows, "rep_affinity", condition="affinity"), + ] + ) + + _, report = pool_replicates(frame, None, ["specificity"]) + + assert [(entry["condition"], entry["gate"]) for entry in report] == [("specificity", "g1")] + + +def test_every_retained_condition_is_reported(): + """The filter is the exclusion list, not "one condition".""" + g1_rows = [row for row in BASE_ROWS if row[0] == "g1"] + frame = pl.concat( + [ + reads_frame(BASE_ROWS, condition="specificity"), + replicate_frame(g1_rows, "rep_specificity", condition="specificity"), + reads_frame(BASE_ROWS, condition="affinity"), + replicate_frame(g1_rows, "rep_affinity", condition="affinity"), + ] + ) + + _, report = pool_replicates(frame, None, ["affinity", "specificity"]) + + assert [entry["condition"] for entry in report] == ["affinity", "specificity"] diff --git a/software/tests/unit/test_validate.py b/software/tests/unit/test_validate.py index 8977f0a..3bbebd7 100644 --- a/software/tests/unit/test_validate.py +++ b/software/tests/unit/test_validate.py @@ -1,7 +1,7 @@ -"""The two data-value refusals. +"""The sort-fraction refusal — the one data-value refusal. -Each asserts the refusal fires **and** that its message names the offending values — the -message is the whole of what the user gets, since nothing partial is produced. +Each case asserts the refusal fires **and** that its message names the offending values — +the message is the whole of what the user gets, since nothing partial is produced. """ from __future__ import annotations @@ -11,70 +11,11 @@ from conftest import BASE_ROWS, reads_frame from errors import Refusal -from validate import check_one_sample_per_group, check_sort_fractions +from validate import check_sort_fractions RETAINED = ["pH7"] -def test_one_sample_per_group_passes_on_the_v1_grain(): - check_one_sample_per_group(reads_frame(BASE_ROWS), RETAINED) - - -def test_two_samples_in_one_group_refuses_and_names_them(): - """Reads are not pooled and neither sample is preferred: summing them would move every - depth, frequency and weighted mean in the run, invisibly.""" - frame = pl.concat( - [ - reads_frame(BASE_ROWS), - pl.DataFrame( - { - "sampleId": ["replicate_of_g1"], - "variantKey": ["P"], - "reads": [5], - "condition": ["pH7"], - "gate": ["g1"], - }, - schema_overrides={"reads": pl.Int64}, - ), - ] - ) - - with pytest.raises(Refusal) as excinfo: - check_one_sample_per_group(frame, RETAINED) - - message = str(excinfo.value) - assert "replicate_of_g1" in message - assert "'g1'" in message - assert "not pooled" in message - - -def test_a_group_with_no_sample_is_not_an_error(): - """That gate was not collected at that condition — clauses 1 and 2 accommodate it.""" - rows = [row for row in BASE_ROWS if row[0] != "g3"] - check_one_sample_per_group(reads_frame(rows), RETAINED) - - -def test_duplicate_in_an_excluded_condition_is_ignored(): - """An excluded value is not part of the run.""" - frame = pl.concat( - [ - reads_frame(BASE_ROWS), - reads_frame([("g1", "P", 5)], condition="pH5"), - pl.DataFrame( - { - "sampleId": ["other"], - "variantKey": ["P"], - "reads": [5], - "condition": ["pH5"], - "gate": ["g1"], - }, - schema_overrides={"reads": pl.Int64}, - ), - ] - ) - check_one_sample_per_group(frame, RETAINED) - - # --------------------------------------------------------------------------- # Sort fractions. # --------------------------------------------------------------------------- diff --git a/ui/src/components/RunStatistics.vue b/ui/src/components/RunStatistics.vue index 5a39d6a..041b282 100644 --- a/ui/src/components/RunStatistics.vue +++ b/ui/src/components/RunStatistics.vue @@ -12,6 +12,7 @@ * Open state is a local `ref`, not `BlockData`: in `data` it is shared, so one person opening * the dialog would open it for everyone with the project open. */ +import { getSingleColumnData, type PObjectId } from "@platforma-sdk/model"; import { PlAgOverlayLoading, PlAlert, @@ -19,6 +20,7 @@ import { PlDialogModal, PlLogView, PlMaskIcon24, + useWatchFetch, } from "@platforma-sdk/ui-vue"; import { computed, ref } from "vue"; import { useApp } from "../app"; @@ -49,6 +51,49 @@ const parentAbsence = computed(() => { return "No mutation-count column was available, so no bin score was produced."; }); +/** + * `sampleId` -> sample name. The label column has one axis, so `axesData`'s single entry lines + * up with `data` by position. Empty until resolved, or where upstream published no labels. + */ +const sampleLabels = useWatchFetch( + () => ({ + pframe: app.model.outputs.sampleLabelPframe, + columnId: app.model.outputs.sampleLabelColumnId, + }), + async ({ pframe, columnId }) => { + const out: Record = {}; + if (!pframe || !columnId) return out; + + const column = await getSingleColumnData(pframe, columnId as PObjectId); + const ids = Object.values(column?.axesData ?? {})[0] ?? []; + const labels = column?.data ?? []; + for (const [index, id] of ids.entries()) { + const label = labels[index]; + if (id !== null && id !== undefined && label !== null && label !== undefined) { + out[String(id)] = String(label); + } + } + return out; + }, +); + +/** The sample's name, falling back to the raw id. */ +function labelFor(sampleId: string): string { + return sampleLabels.value?.[sampleId] ?? sampleId; +} + +/** The pooled groups, or `undefined` where nothing was pooled — see `parentAbsence` on why. */ +const pooling = computed(() => { + const groups = manifest.value?.pooledGroups; + if (!groups || groups.length === 0) return undefined; + return { + lines: groups.map( + (group) => `${group.condition} / ${group.gate}: ${group.samples.map(labelFor).join(" + ")}`, + ), + fractionsDiffer: groups.some((group) => group.sortFractionsDiffer === true), + }; +}); + /** Already ordered by declared gate rank, so this reads along the binding axis. */ function gateList(gates: { gate: string; depth: number }[]): string { return gates.map((entry) => `${entry.gate} (${entry.depth})`).join(", "); @@ -104,6 +149,19 @@ function binScoreCell(entry: { {{ parentAbsence }} + + + These gates were collected more than once, and their reads were summed. Depths, frequencies + and every score below are over the pooled reads. +
    +
  • {{ line }}
  • +
+ +
+ diff --git a/ui/src/components/SettingsDrawer.vue b/ui/src/components/SettingsDrawer.vue index 0f1d15e..4bc7fb6 100644 --- a/ui/src/components/SettingsDrawer.vue +++ b/ui/src/components/SettingsDrawer.vue @@ -139,7 +139,8 @@ function setGateColumn(ref: SUniversalPColumnId | undefined) { @update:model-value="setGateColumn" >