feat: add the rustqc dna subcommand (depth of coverage, samtools, preseq)

[BenjaminDEMAILLE]

Part of #157.

Second of the DNA stack. Stacked on #152, which extracts the shared modules this builds on; review that one first, and the diff here will shrink once it merges.

Adds rustqc dna <BAM>...: a single-pass DNA QC pipeline that needs no annotation.

What it does

One rayon worker per contig, each holding its own depth array and feeding a depth accumulator, a BamStatAccum and a PreseqAccum from the same record stream, so the alignment is read once. A separate pass over unmapped records feeds the counters flagstat and idxstats report.

Outputs land under mosdepth/, samtools/ and preseq/, or flat with --flat-output, alongside rustqc_summary.json and CITATIONS.md.

Parity

Every output is compared against the real upstream tools on a committed fixture, a real public human chr22 slice from nf-core/test-datasets, duplicate-marked locally:

Output	Result
mosdepth.summary.txt	identical
mosdepth.global.dist.txt (1094 lines)	identical
mosdepth.region.dist.txt	identical
per-base.bed.gz (721 intervals)	identical
regions.bed.gz, thresholds.bed.gz	identical
flagstat, idxstats	identical
samtools stats	all 1889 data lines identical

samtools stats headers differ by design: RustQC writes its own # banner naming itself rather than reproducing samtools' version string. That is pre-existing behaviour shared with the rna pipeline. This PR does fix a bug in it, though: the banner hardcoded rustqc rna, so DNA output was labelled as RNA output.

The fixture generation script pins mosdepth 0.3.14 and samtools 1.24 and refuses to run against other versions, recording them in a VERSIONS.txt the test suite asserts against, so fixtures and tool versions cannot drift apart.

Two things worth knowing, both measured rather than assumed

Mate-overlap correction is the dominant behaviour, not an edge case. mosdepth counts a base once when both mates of a pair cover it, unless --fast-mode is given. On this dataset that is the difference between 469875 and 247878 total covered bases, a factor of two. It has its own commit so it can be reviewed on its own.

The distribution files have a non-obvious emission rule, reverse-engineered from the fixtures and checked against every row: every depth from 0 up to min(300, max) gets a row whether or not any base sits at it, and above 300 only depths that occur and lie strictly below the maximum do. So the maximum gets a row when it falls inside the dense range and none when it does not. The global distribution tops out at 866 with a maximum of 867, while the region distribution does emit its maximum of 204. The region distribution is also over windows and their rounded mean depth, not over bases.

Defaults

-Q/--mapq defaults to 0 for dna, matching mosdepth, rather than the 30 rna uses. Every other shared flag keeps its rna name, short form and RUSTQC_* environment variable.

--max-depth-workers bounds how many contig depth arrays are live at once. Each worker costs four bytes per base, about 1 GB for GRCh38 chr1, so the default is a 4 GB budget divided by the largest contig. There is no portable way to read free memory from std, so that budget is a constant; the flag is the escape hatch.

Not in this PR

Picard metrics (CollectWgsMetrics, CollectInsertSizeMetrics, CollectGcBiasMetrics, CollectHsMetrics) and Qualimap bamqc land in the following PRs. --targets is accepted and warns that targeted metrics are not implemented yet. The .csi companion indexes for the bgzf outputs are not written yet either.

One inconsistency to flag: CITATIONS.md cites samtools v1.22.1, the version constant the RNA pipeline was validated against, while these fixtures were generated with 1.24. I left the constant alone rather than silently claiming the RNA pipeline was revalidated.

Tests

283 green, up from 233. Includes 15 parity tests, and unit tests covering the paths the fixture cannot reach, notably the MAPQ filter: mosdepth -Q 30 returns the same totals as the default on this data, so that path is covered on synthetic records instead.

🤖 Generated with Claude Code