bam: optional htslib writer behind the htslib-bam feature

[BenjaminDEMAILLE]

Optional htslib BAM writer behind the htslib-bam feature.

What changed

src/io/bam_htslib.rs adds a second unsorted-BAM-to-file backend using rust-htslib, selected by the htslib-bam feature. htslib's hts_set_threads gives a genuinely parallel deflate; worker count is --runThreadN - 1, so a single-threaded run stays serial.

Only the unsorted-BAM-to-file path is taken over, the one where BAM writing can dominate a run. Stdout and the sorted path keep the default writer.

Why

BAM writing was measured as the dominant cost at high thread counts, and noodles_bgzf::MultithreadedWriter was tried first and rejected (~37% slower; numbers on #143 and #151). htslib is a different mechanism, not a layer over the same compressor. Whether it is actually faster is the open question below.

The measurement, and why yeast cannot settle it

The earlier revision reported five interleaved A/B pairs at 16 threads, medians 2.51 s for noodles against 2.08 s for htslib, "~17% faster", with no command, input or platform. test/bench_bam_backends.sh is now in the branch and does that comparison reproducibly.

First I ran it badly and posted the result anyway. I did not check whether the machine was idle, and it runs unrelated 15-core jobs. The tail of my series drifted from 2.7 s to 5.1 s on both backends, which I attributed to the machine warming up; that is what an unrelated job arriving mid-series looks like. Those numbers are withdrawn. The script now samples CPU idle before building and refuses below 80%, printing the idle figure beside every pair, so a contaminated series is visible in its own output.

Re-run on an idle machine (90-93% idle throughout), yeast index, 200k reads of ERR12389696, --outSAMtype BAM Unsorted --runThreadN 16, eight pairs with the order flipped on even ones:

	median	range
noodles	2.46 s	2.33 to 2.56
htslib	2.50 s	2.41 to 2.69

But that comparison cannot resolve the claim, and neither could the original. On this data BAM writing is a rounding error in the run:

workload	--outSAMtype None	BAM Unsorted	BAM share
200k reads	2.29 s	2.36 s	~3%
2M reads	21.9, 22.4 s	22.8, 22.6 s	1 to 4%

A 17% improvement in BAM writing would move total wall clock by about half a percent at 200k, well under the run-to-run spread. The 0.04 s between the two medians above is inside that noise, so it is not evidence either way.

So the honest position is that the headline is unverified, not refuted, and it is out of the title until it is measured on a workload where BAM writing is a meaningful share of the run. On yeast it never is, at any scale reachable here. The premise in the section above, that BAM writing dominates at high thread counts, comes from human-scale data with far larger records and output; that is the workload this benchmark needs, and the revision that quotes a number should say which one it used.

Reproduce:

test/bench_bam_backends.sh "$RUSTAR_IDX" "$DATA/reads/subset.fastq" 16 8

How the backends are kept honest

Both are handed the same rendered header, and every record goes through the same SAM rendering before htslib parses it, so they cannot drift in field encoding. backends_agree_on_the_decoded_record_stream states the contract: BGZF block boundaries may differ, that is framing, and htslib chooses them differently, but the decoded record stream may not.

Verification

• cargo test --release (default features): 561 lib + integration, green
• cargo test --release --features htslib-bam: green
• cargo clippy --all-targets -- -D warnings, with and without the feature: clean
• cargo fmt --check: clean

Trade-off on record

This pulls htslib (C) into a project that is otherwise pure Rust apart from libdeflater, on a crate published to crates.io. Optional and platform-gated limits the blast radius but does not remove the dependency-story change. If the dependency is unwanted, the measurement still stands on its own: it locates where the BAM write cost is.

[Psy-Fer]

I don't think it's worth pulling in htslib and it's C for such a small performance change.
If it matters, it would be better to find the bottleneck in noodles and put a PR there, then everyone gets the benefit.

Thoughts?

James

[BenjaminDEMAILLE]

Agreed, and my own measurement since opening this says the same thing more strongly than your objection does.

I tried to verify the 17% headline and could not, because on yeast BAM writing is 1-4% of the run at any reachable scale: --outSAMtype None 21.9 s against BAM Unsorted 22.8 s at 2M reads. A 17% improvement on that slice moves total wall clock by roughly half a percent, below the run-to-run spread. So the number in the original description is unverified rather than refuted, and this workload cannot settle it either way. (I also had to retract a "~2% slower" figure I posted here: it was taken with an unrelated job at 1522% CPU.)

That removes the case for a C dependency, since there is no demonstrated win to weigh against it.

Your alternative is the right one and it turned out not to need noodles-side work at all: noodles-bgzf already ships MultithreadedWriter, and noodles-bgzf with the libdeflate feature plus libdeflater were already dependencies here. Only the wiring was missing. That is #164, which switches make_bgzf_writer to it and also splits record encoding across the pool, byte-identical to the current output and with no new dependency.

Closing this in favour of #164.