Skip to content
Merged
2 changes: 1 addition & 1 deletion .coderabbit.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -329,7 +329,7 @@ reviews:
require it to be called out. Do NOT file wording-only nits on help text
or doc comments unless the text states a behavior the code does not
implement.
- path: "src/lib/{grouper,template,mi_group,read_info}.rs"
- path: "src/lib/{grouper,template,template_filter,mi_group,read_info}.rs"
instructions: >-
Core read-grouping and template-assembly logic that sits outside the
crates but feeds the same fgbio-validated output as `fgumi-umi` and
Expand Down
18 changes: 18 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,16 @@ All notable changes to this project will be documented in this file.

### Changed

- **Breaking:** `dedup --metrics` gains eleven columns: one `filtered_*` column per template
filter reason, `filtered_templates` for the total, and `passthrough_templates` for the
`--include-unmapped` templates that bypass the filter. Consumers that parse the file
positionally must be updated; consumers keyed on column names are unaffected. Note that
`fgumi compare metrics` treats a column-set difference as `DIFFER`, so comparing a
pre-change dedup metrics file against a post-change one will fail by design.

These columns count *templates*, whereas the `discarded_*` columns of `group`'s
fgbio-shaped `.grouping_metrics.txt` count *primary records*. `group`'s output is unchanged.

- **Breaking:** consensus downsampling now selects which reads to retain by a hash of the read
name rather than by shuffling a seeded random number generator. This affects `simplex
--max-reads`, `duplex --max-reads-per-strand`, and `codec --max-reads`.
Expand All @@ -28,6 +38,14 @@ All notable changes to this project will be documented in this file.

### Bug Fixes

- `dedup --metrics` now reports the templates dropped by its template filter, broken out by
reason. `dedup` is a read filter as well as a duplicate marker, but the drop counts were
collected, merged across workers, and then discarded at the serialization boundary: a run
that filtered out every input template wrote a row of zeros and logged only "Deduplication
complete: 0 templates", with nothing to indicate records had been dropped or why. A
per-reason summary is now also logged, so the drops are visible without `--metrics`
([#739](https://github.com/fulcrumgenomics/fgumi/issues/739)).

- `duplex --max-reads-per-strand 0` is now rejected at startup. It previously exited 0 having
written an empty BAM, because a cap of zero empties every strand. `simplex` and `codec`
already validated their equivalents.
Expand Down
164 changes: 164 additions & 0 deletions crates/fgumi-metrics/src/dedup.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,164 @@
//! Metrics for the `dedup` command.

use serde::{Deserialize, Serialize};

use crate::Metric;

/// Metrics written by `fgumi dedup --metrics`.
///
/// Three sections, in column order:
///
/// 1. **Filtering** — `filtered_templates` through `passthrough_templates`.
/// `dedup` is a read *filter* as well as a duplicate marker: templates failing
/// the pre-grouping filter are dropped and never reach the output. The
/// `filtered_*` columns give the count per reason in **template** units — note
/// that `group`'s fgbio-shaped `.grouping_metrics.txt` counts primary
/// *records*, so the two files' discard counts are not comparable.
/// 2. **Deduplication** — `total_templates` through `duplicate_reads`, covering
/// only what survived filtering. These count what *passed the filter*, not what
/// reaches the output file: under `--remove-duplicates` the duplicates are
/// dropped at write time, after they have been counted here.
/// 3. **Diagnostics** — `secondary_reads`, `supplementary_reads`,
/// `missing_tc_tag`.
///
/// Templates read from the input are `filtered_templates + total_templates`; that
/// sum is not emitted as its own column because it is derivable, matching
/// `UmiGroupingMetrics`' treatment of its derivable fields.
// NB: plain code spans above, never intra-doc links — this doc is rendered verbatim
// into docs/src/metrics/deduplication-metrics.md, where a rustdoc link would come
// out as literal markdown. Keep field docs to a single short line: they become the
// Description column of that page's table.
#[derive(Debug, Clone, Default, Serialize, Deserialize)]
pub struct DeduplicationMetrics {
/// Templates dropped by the filter (any reason)
pub filtered_templates: u64,
/// Templates dropped (record shorter than the minimum BAM record length)
pub filtered_malformed_record: u64,
/// Templates dropped (no primary R1 or R2)
pub filtered_no_primary_reads: u64,
/// Templates dropped (no mapped primary read)
pub filtered_unmapped: u64,
/// Templates dropped (QC-fail flag)
pub filtered_not_passing_filter: u64,
/// Templates dropped (mapping quality below `--min-map-q`)
pub filtered_low_mapping_quality: u64,
/// Templates dropped (mate `MQ` below `--min-map-q`)
pub filtered_low_mate_mapping_quality: u64,
/// Templates dropped (no UMI tag)
pub filtered_missing_umi: u64,
/// Templates dropped (N base in the UMI)
pub filtered_ns_in_umi: u64,
/// Templates dropped (UMI shorter than `--min-umi-length`)
pub filtered_umi_too_short: u64,
/// Templates emitted untouched by `--include-unmapped`
pub passthrough_templates: u64,
/// Templates that passed the filter
pub total_templates: u64,
/// Templates kept as their family's representative
pub unique_templates: u64,
/// Templates marked as duplicates
pub duplicate_templates: u64,
/// `duplicate_templates` / `total_templates`
#[serde(with = "crate::float")]
pub duplicate_rate: f64,
/// Reads in templates that passed the filter
pub total_reads: u64,
/// Reads not marked as duplicates
pub unique_reads: u64,
/// Reads marked as duplicates
pub duplicate_reads: u64,
/// Secondary alignments that passed the filter
pub secondary_reads: u64,
/// Supplementary alignments that passed the filter
pub supplementary_reads: u64,
/// Secondary/supplementary reads lacking the `tc` tag
pub missing_tc_tag: u64,
}

impl DeduplicationMetrics {
/// Creates a metrics struct with all counts initialized to zero.
#[must_use]
pub fn new() -> Self {
Self::default()
}
}

impl Metric for DeduplicationMetrics {
fn metric_name() -> &'static str {
"deduplication"
}
}

impl crate::ProcessingMetrics for DeduplicationMetrics {
fn total_input(&self) -> u64 {
self.filtered_templates + self.total_templates
}

fn total_output(&self) -> u64 {
self.total_templates
}

fn total_filtered(&self) -> u64 {
self.filtered_templates
}
}

#[cfg(test)]
mod tests {
use super::*;
use crate::template_filter::TemplateFilterReason;
use fgoxide::io::DelimFile;
use tempfile::NamedTempFile;

fn header_of(metrics: DeduplicationMetrics) -> String {
let file = NamedTempFile::new().expect("temp file");
DelimFile::default().write_tsv(file.path(), [metrics]).expect("write");
std::fs::read_to_string(file.path())
.expect("read")
.lines()
.next()
.expect("header")
.to_string()
}

/// Every filter reason must have a column. This is the test that would have
/// caught fgumi#739: the counts existed and were merged correctly, then were
/// dropped at the serialization boundary because no column held them.
#[test]
fn every_filter_reason_has_a_column() {
let header = header_of(DeduplicationMetrics::default());
let columns: Vec<&str> = header.split('\t').collect();
for reason in TemplateFilterReason::ALL {
let column = format!("filtered_{}", reason.column_suffix());
assert!(columns.contains(&column.as_str()), "missing {column} for {reason:?}");
}
}

/// Column order is the file's contract; pin it whole so a field reorder is a
/// deliberate, reviewed change rather than a silent schema break.
#[test]
fn serializes_expected_columns_in_order() {
assert_eq!(
header_of(DeduplicationMetrics::default()),
"filtered_templates\tfiltered_malformed_record\tfiltered_no_primary_reads\t\
filtered_unmapped\tfiltered_not_passing_filter\tfiltered_low_mapping_quality\t\
filtered_low_mate_mapping_quality\tfiltered_missing_umi\tfiltered_ns_in_umi\t\
filtered_umi_too_short\tpassthrough_templates\ttotal_templates\tunique_templates\t\
duplicate_templates\tduplicate_rate\ttotal_reads\tunique_reads\tduplicate_reads\t\
secondary_reads\tsupplementary_reads\tmissing_tc_tag"
);
}

#[test]
fn processing_metrics_reconciles_input_from_the_columns() {
use crate::ProcessingMetrics;
let metrics = DeduplicationMetrics {
filtered_templates: 10,
total_templates: 90,
..DeduplicationMetrics::default()
};
assert_eq!(metrics.total_input(), 100);
assert_eq!(metrics.total_output(), 90);
assert_eq!(metrics.total_filtered(), 10);
}
}
94 changes: 94 additions & 0 deletions crates/fgumi-metrics/src/group.rs
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
use serde::{Deserialize, Serialize};

use crate::Metric;
use crate::template_filter::{TemplateFilterCounts, TemplateFilterReason};

/// Build a size distribution from (size, count) pairs.
///
Expand Down Expand Up @@ -107,6 +108,45 @@ impl UmiGroupingMetrics {
pub fn new() -> Self {
Self::default()
}

/// Builds the fgbio filter columns from [`TemplateFilterCounts`].
///
/// fgbio's `UmiGroupingMetric` counts **primary records**, not templates, so
/// this reads the primary-read side of `counts`. The exhaustive match makes
/// adding a new [`TemplateFilterReason`] a compile error until it is assigned
/// a column, so a rejection can never be silently dropped from the output.
///
/// Six reasons share `discarded_poor_alignment` because fgbio's schema has no
/// finer bucket for them. That collapse is lossy and is retained only for
/// fgbio parity; `dedup` reports each reason separately.
///
/// Leaves the family-size and summary fields at their defaults; the caller
/// fills those.
#[must_use]
pub fn from_filter_counts(counts: &TemplateFilterCounts) -> Self {
let mut metrics = Self {
total_records: counts.total_primary_reads(),
accepted_records: counts.accepted_primary_reads(),
..Self::default()
};

for reason in TemplateFilterReason::ALL {
let reads = counts.rejected_primary_reads(reason);
let column = match reason {
TemplateFilterReason::NotPassingFilter => &mut metrics.discarded_non_pf,
TemplateFilterReason::NsInUmi => &mut metrics.discarded_ns_in_umi,
TemplateFilterReason::UmiTooShort => &mut metrics.discarded_umi_too_short,
TemplateFilterReason::MalformedRecord
| TemplateFilterReason::NoPrimaryReads
| TemplateFilterReason::Unmapped
| TemplateFilterReason::LowMappingQuality
| TemplateFilterReason::LowMateMappingQuality
| TemplateFilterReason::MissingUmi => &mut metrics.discarded_poor_alignment,
};
*column += reads;
}
metrics
}
}

impl Metric for UmiGroupingMetrics {
Expand Down Expand Up @@ -258,6 +298,60 @@ mod tests {
assert_eq!(metrics.unique_molecule_ids, 0);
}

/// fgbio columns are in **primary-read** units, and the six reasons fgbio has
/// no bucket for must all land in `discarded_poor_alignment`. This mapping is
/// what keeps group's output byte-identical across the `TemplateFilterCounts`
/// refactor.
#[rstest::rstest]
#[case::malformed(TemplateFilterReason::MalformedRecord, 0, 2, 0, 0)]
#[case::no_primary(TemplateFilterReason::NoPrimaryReads, 0, 2, 0, 0)]
#[case::unmapped(TemplateFilterReason::Unmapped, 0, 2, 0, 0)]
#[case::non_pf(TemplateFilterReason::NotPassingFilter, 2, 0, 0, 0)]
#[case::low_mapq(TemplateFilterReason::LowMappingQuality, 0, 2, 0, 0)]
#[case::low_mate_mapq(TemplateFilterReason::LowMateMappingQuality, 0, 2, 0, 0)]
#[case::missing_umi(TemplateFilterReason::MissingUmi, 0, 2, 0, 0)]
#[case::ns_in_umi(TemplateFilterReason::NsInUmi, 0, 0, 2, 0)]
#[case::umi_too_short(TemplateFilterReason::UmiTooShort, 0, 0, 0, 2)]
fn set_filter_counts_maps_reasons_to_fgbio_columns(
#[case] reason: TemplateFilterReason,
#[case] expected_non_pf: u64,
#[case] expected_poor_alignment: u64,
#[case] expected_ns_in_umi: u64,
#[case] expected_umi_too_short: u64,
) {
let mut counts = TemplateFilterCounts::new();
counts.record_accepted(2);
counts.record_rejected(reason, 2);

let metrics = UmiGroupingMetrics::from_filter_counts(&counts);

assert_eq!(metrics.accepted_records, 2, "accepted is in primary-read units");
assert_eq!(metrics.total_records, 4);
assert_eq!(metrics.discarded_non_pf, expected_non_pf);
assert_eq!(metrics.discarded_poor_alignment, expected_poor_alignment);
assert_eq!(metrics.discarded_ns_in_umi, expected_ns_in_umi);
assert_eq!(metrics.discarded_umi_too_short, expected_umi_too_short);
}

/// No rejection may be dropped: the four fgbio columns must sum to the total.
#[test]
fn set_filter_counts_conserves_every_rejection() {
let mut counts = TemplateFilterCounts::new();
counts.record_accepted(2);
for reason in TemplateFilterReason::ALL {
counts.record_rejected(reason, 2);
}

let metrics = UmiGroupingMetrics::from_filter_counts(&counts);

let discarded = metrics.discarded_non_pf
+ metrics.discarded_poor_alignment
+ metrics.discarded_ns_in_umi
+ metrics.discarded_umi_too_short;
assert_eq!(discarded, counts.total_rejected_primary_reads());
assert_eq!(metrics.total_records, metrics.accepted_records + discarded);
}

#[test]
fn test_umi_grouping_metrics_serializes_fgbio_five_columns() {
use fgoxide::io::DelimFile;
Expand Down
4 changes: 4 additions & 0 deletions crates/fgumi-metrics/src/lib.rs
Original file line number Diff line number Diff line change
Expand Up @@ -12,12 +12,14 @@
pub mod clip;
pub mod consensus;
pub mod correct;
pub mod dedup;
pub mod duplex;
pub mod float;
pub mod group;
pub mod rejection;
pub mod shared;
pub mod simplex;
pub mod template_filter;
pub mod writer;

use serde::{Deserialize, Serialize};
Expand Down Expand Up @@ -98,6 +100,7 @@ pub trait ProcessingMetrics {
pub use clip::{ClipCounts, ClippingMetrics, ClippingMetricsCollection, ReadType};
pub use consensus::{ConsensusCallerKind, ConsensusKvMetric, ConsensusMetrics};
pub use correct::UmiCorrectionMetrics;
pub use dedup::DeduplicationMetrics;
pub use duplex::{
DuplexFamilySizeMetric, DuplexMetricsCollector, DuplexUmiMetric, DuplexYieldMetric,
FamilySizeMetric,
Expand All @@ -106,6 +109,7 @@ pub use group::{FamilySizeMetrics, PositionGroupSizeMetrics, UmiGroupingMetrics}
pub use rejection::{RejectionReason, format_count};
pub use shared::UmiMetric;
pub use simplex::{SimplexFamilySizeMetric, SimplexMetricsCollector, SimplexYieldMetric};
pub use template_filter::{TemplateFilterCounts, TemplateFilterReason};
pub use writer::{read_metrics, read_metrics_auto, write_metrics};

#[cfg(test)]
Expand Down
Loading
Loading