fix(cram): read CRAM files that use Huffman coding for the bases series (#429)

[mwiewior]

Fixes #429.

Draft — depends on biodatageeks/datafusion-bio-formats#225. This branch pins datafusion-bio-formats by revision so it can be built and tested now; that pin should become a released tag before merge.

What was broken

Reading a CRAM whose bases series uses Huffman coding aborted the query with not yet implemented, raised from noodles's Byte::decode_take. It affected scan_cram/read_cram and pb.depth, which is what @J-Wall was calling in #429.

Only the single-byte Byte::decode path handled Huffman; the bulk decode_take path used for the bases series did not. Upstream implemented it in zaeleus/noodles#393, so the fix is to stop being 727 commits behind: our noodles fork has been rebased onto upstream master, and datafusion-bio-formats repointed at it.

What's in the rebase

Two of our patches were dropped as obsolete — upstream fixed the CRAM container landmark calculation itself, and migrated xz2 → lzma-rust2. Everything still needed was replayed, three patches by hand against changed upstream APIs. Details in biodatageeks/datafusion-bio-formats#225.

The rebase also surfaced a second bug, found because this repo's CRAM tests started failing: crai::io::reader::read_index() reused its line buffer without clearing it, concatenating each record onto the previous one, so any CRAI with more than one record failed with invalid digit found in string. Fixed at the source with regression tests; it deserves an upstream PR too.

All noodles crates now resolve to a single revision. Two divergent fork revisions used to be pinned at once, so cargo was building duplicate copies of noodles-core, noodles-csi and noodles-bgzf.

Tests

tests/test_io_cram_huffman.py, with a 2.5 KB fixture. Five of its six tests fail on master with the reported panic — I checked, because a regression test that passes before the fix is worthless. The sixth pins the CRAI behaviour below and passes either way by design.

Getting a fixture that reproduces the bug took some care, and the test module documents why: the unmapped reads are required (the bases series is only read for unmapped records), and the fixture must ship without a .crai (an indexed scan skips unmapped reads and never reaches the Huffman block).

Full suite: 1074 passed, 3 skipped against the lockfile-pinned environment, versus 1068 on master — the six new tests, no regressions.

Documentation — CRAI limitations

While building the fixture I found something worth writing down, so docs/features/reading.md gains a CRAM index (CRAI) limitations section:

# sample.cram holds 300 mapped reads and 200 unmapped ones
pb.scan_cram("sample.cram").collect().height   # 300 — with sample.cram.crai
pb.scan_cram("sample.cram").collect().height   # 500 — without it

A CRAI records only placed alignment slices and, unlike BAI/CSI, cannot address the unmapped tail. So a whole-file scan of an indexed CRAM silently drops unmapped reads. BAM is not affected — an indexed scan_bam returns them, which I verified rather than assumed. Region queries are unaffected, since they ask for placed reads by definition.

This is pre-existing behaviour, not something this PR changes, but it is a quiet way to lose records and was undocumented.

Not covered here

pb.depth reaches CRAM through the pileup crate rather than the table provider. The test_depth_on_huffman_cram test above covers it end to end from Python, but a Rust-level test belongs in datafusion-bio-functions.

[mwiewior]

Updated: the CRAI limitation this PR originally documented is fixed, not documented, and the pin now points at datafusion-bio-formats f4e0854.

I need to correct something I got wrong in the first version of this PR. I wrote that a CRAI "cannot address the unmapped tail" and documented the dropped reads as a property of the format. That is not true. Dumping the index of the test fixture:

0   1  1550  327   196  503     ← chr1
1   1  1550  1047  196  505     ← chr2
-1  0  1     1770  185  494     ← the unmapped slice

The unmapped slice is right there, and noodles has exposed Reader::query_unmapped for it since 0.90.0. The dropped reads were a gap in our own CRAM reader, so the honest fix was to close it rather than write the limitation down.

Behaviour now: a whole-file scan_cram/read_cram returns the same records with or without a .crai. Unmapped reads surface with a null chrom. Region queries are unchanged — they ask for placed reads by definition, so they never return the tail. BAM was never affected.

# sample.cram: 300 mapped + 200 unmapped
pb.scan_cram("sample.cram").collect().height                              # 500, indexed or not
pb.scan_cram("sample.cram").filter(pl.col("chrom").is_null()).collect()   # the 200
pb.scan_cram("sample.cram").filter(pl.col("chrom") == "chr1").collect()   # 150, no tail

The docs section is renamed from "CRAM index (CRAI) limitations" to "Unmapped reads and indexed scans", since it now describes behaviour rather than a caveat, and carries a short note for anyone on an older version who worked around it by deleting their index. The capability-matrix footnote is gone.

TestCRAMIndexUnmappedReads is replaced by TestCRAMIndexedUnmappedReads, which asserts the indexed and sequential scans agree, that unmapped reads have a null chrom, that no record is duplicated, and that a region query still excludes the tail.

Full suite: 1076 passed, 3 skipped against the lockfile-pinned environment.

Still draft, for the same reason as before: the datafusion-bio-formats pin is a revision and should become a released tag once #225 merges.

[mwiewior]

@codex pls review

[mwiewior]

@claude pls review

[claude]

Claude finished @mwiewior's task in 2m 23s —— View job

---

Code Review

• Gather context and understand the request
• Read changed files (CHANGELOG.md, Cargo.toml, docs, tests)
• Review test coverage (tests/test_io_cram_huffman.py)
• Review documentation changes (docs/features/reading.md)
• Review Cargo dependency changes
• Post comprehensive review

---

This is a well-structured PR that fixes a real bug, closes a silent data-loss gap in the indexed CRAM path, and documents the corrected behaviour. The test design is notably good. A few observations below.

Hard Blocker

Revision pin in Cargo.toml (rev = "f4e08540...") must become a released tag before merge. You've already called this out in the PR description, so this is just confirming it's the only thing blocking merge readiness — nothing can be done until datafusion-bio-formats#225 lands and gets a tag.

Tests (tests/test_io_cram_huffman.py)

Strengths

• The module docstring explaining why specific fixture properties are load-bearing (unmapped reads required; no .crai intentional) is exactly what future maintainers need. Rare to see this done right.
• test_read_matches_pysam cross-validates against an independent reader — the correct way to write a decode-correctness test.
• The explicit constant assert set(df["sequence"].unique().to_list()) == {"N" * READ_LENGTH} pins the fixture expectation clearly and would catch a decoding bug that produced wrong bytes rather than a panic.
• TestCRAMIndexedUnmappedReads covers all four behaviours worth pinning: agreement (indexed == unindexed), null chrom for unmapped, no duplicates, region query excludes the tail.

One observation

test_indexed_scan_does_not_duplicate_records asserts df["name"].n_unique() == TOTAL_READS, which assumes all 500 read names are unique. If the fixture holds paired-end reads, mates share a name and this test would fail spuriously. The module docstring doesn't say one way or the other. If these are paired-end reads (300 mapped pairs + 200 unmapped pairs is a natural shape), this assertion is fragile. Worth either documenting that the fixture is single-end, or switching to df.height == TOTAL_READS and df.is_duplicated().sum() == 0 on a subset of columns that are stable per alignment.

test_depth_on_huffman_cram asserts only result.height > 0. Since the goal is confirming that pb.depth no longer panics, that's acceptable as the primary assertion. Adding {"chr1", "chr2"} <= set(result["contig"].unique().to_list()) (note the direction — both contigs should be present, not just a subset of them) would make it slightly stronger, but this is minor.

Documentation (docs/features/reading.md)

The new "Unmapped reads and indexed scans" section is well-placed, clear, and the !!! note "Versions before 0.34.0" admonition is exactly the right pattern for migration guidance. The example numbers (300/200/150) match the test fixture, which is a nice consistency.

One cosmetic nit: the diff introduces an extra blank line between the capability matrix table and ## Performance features (line 35 in the diff). Not a problem but it looks unintentional.

CHANGELOG.md

Well-organised across Added / Fixed / Changed / Documentation sections. The Fixed bullet for the CRAI line-buffer bug is clear. One suggestion: the sentence BAM was never affected at the end of the Added bullet reads slightly awkward without a period. Fix this →

Cargo / Dependency (Cargo.toml, Cargo.lock)

• Deduplicating the three noodles crates (noodles-core, noodles-bgzf, noodles-csi) from two divergent fork revisions down to one is a clean improvement — fewer symbols, faster builds, no stale artefacts.
• New transitive deps (borsh, borsh-derive, bit-vec 0.10.1) look like they come from updated noodles internals (likely the Byte::decode_take Huffman implementation). Nothing alarming.
• datafusion-bio-format-{bam,bbi,...} minor version bumps in the lockfile (1.8.7 → 1.8.8, 1.6.6 → 1.6.7) are consistent with the upstream crate revisions.

What's Not Here (per PR description)

A Rust-level test for pb.depth in datafusion-bio-functions is called out as out-of-scope. The Python end-to-end test_depth_on_huffman_cram does cover the regression path from the user's perspective, so this is a reasonable deferral — just worth tracking as a follow-up.

Summary

The code itself looks correct and the test design is high quality. The only merge-blocking item is replacing the git revision pin with a proper released tag once upstream PR #225 lands. Once that's done, this is ready to merge.

---

[chatgpt-codex-connector]

Codex Review: Didn't find any major issues. 👍

Reviewed commit: 2d4822bf4a

ℹ️ About Codex in GitHub

Your team has set up Codex to review pull requests in this repo. Reviews are triggered when you

• Open a pull request for review
• Mark a draft as ready
• Comment "@codex review".

If Codex has suggestions, it will comment; otherwise it will react with 👍.

Codex can also answer questions or update the PR. Try commenting "@codex address that feedback".