Repository navigation
docs: add TAPs guide and update EM-Seq guide - #174
Conversation
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## main #174 +/- ##
=======================================
Coverage 89.24% 89.24%
=======================================
Files 120 120
Lines 59642 59642
=======================================
+ Hits 53228 53229 +1
+ Misses 6414 6413 -1 ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
|
@coderabbitai review |
✅ Actions performedReview triggered.
|
📝 WalkthroughWalkthroughUpdated methylation pipeline documentation and related build configuration. The guide was generalized to cover both EM-Seq and TAPs/Illumina chemistries, with chemistry-specific instructions for conversion targets, alignment tools, and re-alignment parameters. Example commands updated to accept chemistry mode as a parameter rather than hardcoding EM-Seq. The mdBook configuration now references the unified methylation guide instead of EM-Seq-specific content. The documentation footer branding image was also changed from dark to light variant. 🚥 Pre-merge checks | ✅ 3✅ Passed checks (3 passed)
✏️ Tip: You can configure your own custom pre-merge checks in the settings. ✨ Finishing Touches📝 Generate docstrings
🧪 Generate unit tests (beta)
Thanks for using CodeRabbit! It's free for OSS, and your support helps us grow. If you like it, consider giving us a shout-out. Comment |
|
@coderabbitai review |
✅ Actions performedReview triggered.
|
|
@coderabbitai review |
✅ Actions performedReview triggered.
|
There was a problem hiding this comment.
Actionable comments posted: 8
🧹 Nitpick comments (8)
crates/fgumi-consensus/src/codec_caller.rs (1)
456-459: Consider avoiding per-read CIGAR reparse here.
get_cigar_ops+simplify_cigar_from_rawon each call adds avoidable work in a hot path. You can thread precomputed/simplified original CIGAR throughClippedRecordInfoand reuse it.🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@crates/fgumi-consensus/src/codec_caller.rs` around lines 456 - 459, The code reparses the CIGAR per read by calling bam_fields::get_cigar_ops and bam_fields::simplify_cigar_from_raw to build original_cigar in the hot path; instead, compute and simplify the CIGAR once earlier and thread that simplified CIGAR through the ClippedRecordInfo struct so callers reuse it. Add a field (e.g., simplified_original_cigar) to ClippedRecordInfo, populate it where records are first constructed, remove the per-call get_cigar_ops/simplify_cigar_from_raw usage in codec_caller.rs (the original_cigar block), and update any functions that construct ClippedRecordInfo to pass the precomputed simplified CIGAR.src/commands/common.rs (1)
64-88: LGTM: Methylation reference loading.Validates
--refrequirement, builds contig name mapping from header.One note:
ref_timer.log_completion(0)logs 0 records processed—consider logging contig count instead for more useful metrics.🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@src/commands/common.rs` around lines 64 - 88, The ref_timer currently logs 0 records; update load_methylation_reference to log the actual contig count by passing ref_names.len() (cast to the expected numeric type, e.g., as u64 or usize as required) into ref_timer.log_completion instead of 0 so the timer reports the number of contigs processed; reference symbols: load_methylation_reference, ref_timer.log_completion, ref_names.src/commands/filter.rs (1)
859-872: Unwrap is safe here but depends on validation fix above.The
methylation_tags.as_ref().unwrap()at line 866 is guarded byneeds_methylation_tagslogic. However,methylation_modecould still beDisabledif user omits--methylation-mode- the validation fix above would complete the safety guarantee.🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@src/commands/filter.rs` around lines 859 - 872, The unwrap on methylation_tags (methylation_tags.as_ref().unwrap()) in the conversion-fraction block is only safe if upstream validation guarantees methylation_mode is not Disabled when needs_methylation_tags is true; fix by either enforcing that validation or adding a local guard: before calling check_conversion_fraction_raw_with_ref_bases_and_tags (the function used here), ensure methylation_mode != Disabled and methylation_tags.is_some(), and skip the call or set pass=false if those conditions are not met; alternatively update the validation logic referenced above to require a non-Disabled methylation_mode whenever needs_methylation_tags is set so the unwrap is provably safe.crates/fgumi-consensus/src/vanilla_caller.rs (1)
1393-1415: Verify strand defaulting behavior on emptyoriginal_raws.Line 1397-1398:
is_none_ordefaultsis_top = truewhenoriginal_rawsis empty. While this shouldn't occur in practice (methylation requires source reads), the implicit default could mask bugs if the invariant ever breaks.Consider explicit handling:
-let is_top = original_raws - .first() - .is_none_or(|raw| crate::methylation::is_top_strand(fgumi_raw_bam::flags(raw))); +let is_top = original_raws + .first() + .map(|raw| crate::methylation::is_top_strand(fgumi_raw_bam::flags(raw))) + .unwrap_or(true); // Default to top strand; original_raws should never be empty hereNot blocking — current code is functionally correct given callers always have reads.
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@crates/fgumi-consensus/src/vanilla_caller.rs` around lines 1393 - 1415, When methylation is present we currently compute is_top using original_raws.first().is_none_or(...), which silently defaults to true for an empty original_raws; instead explicitly handle the empty case to avoid masking invariants: inside the methylation block check original_raws.first() and if it is None either early-return/log an error or assert (e.g., return Err / debug_assert!) before calling crate::methylation::is_top_strand, then use the real boolean for is_top when calling crate::methylation::build_mm_ml_tags and appending MM/ML and cu/ct tags (referencing original_raws, is_top, crate::methylation::is_top_strand, crate::methylation::build_mm_ml_tags, and self.options.methylation_mode).crates/fgumi-consensus/src/filter.rs (2)
1697-1700: Consider explicit match arms for clarity.The wildcard
_catches bothEmSeqandDisabled. While functionally acceptable (caller guards for enabled mode), explicit arms document intent and prevent silent surprises if a new variant is added.Suggested diff
let numerator = match methylation_mode { crate::MethylationMode::Taps => u64::from(cu), - _ => u64::from(ct), + crate::MethylationMode::EmSeq | crate::MethylationMode::Disabled => { + u64::from(ct) + } };🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@crates/fgumi-consensus/src/filter.rs` around lines 1697 - 1700, The match on methylation_mode currently uses a wildcard for non-Taps variants; change it to explicit arms so each variant is handled (e.g., match methylation_mode { crate::MethylationMode::Taps => u64::from(cu), crate::MethylationMode::EmSeq => u64::from(ct), crate::MethylationMode::Disabled => u64::from(ct), } ), making the intent clear and preventing future silent behavior if new variants are added (update the match in the same scope where numerator is computed).
3154-3310: Missing TAPs mode coverage.All
check_conversion_fraction_rawtests useEmSeq. Add a test withMethylationMode::Tapsto verify the inverted numerator logic.Want me to draft a TAPs-specific test case?
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@crates/fgumi-consensus/src/filter.rs` around lines 3154 - 3310, Add a TAPs-mode unit test mirroring the existing EmSeq tests that verifies the inverted numerator logic for TAPS: create a new test (e.g., test_conversion_fraction_taps_mode_inverts_counts) that builds a reference and a read with cu/ct tags similar to the high-conversion and low-conversion EmSeq cases but flipped for Taps, call check_conversion_fraction_raw with crate::MethylationMode::Taps, and assert the opposite pass/fail expectations (i.e., values that passed under EmSeq should fail under TAPS and vice versa); reference the test helper check_conversion_fraction_raw and the RecordBuilder/add_i16_array_tag usage to locate where to mimic setup.docs/em-seq-guide.md (1)
233-246: Add explicit--allow-c-to-tguidance in the correction step.The doc implies the chemistry expectation, but it should explicitly state when
--allow-c-to-tis needed and that it applies to all UMI segments.✏️ Proposed addition
### Step 2: UMI Correction Correct UMIs against the known whitelist: ```bash fgumi correct \ --input unmapped.bam \ --output corrected.bam \ --umi-files known_umis.txt \ --max-mismatches 1 \ --min-distance 1 \ --metrics correction_metrics.txt \ --threads 8
+If your UMI design includes unmethylated cytosines, add
--allow-c-to-t.
+This flag is applied uniformly across all UMI segments.</details> Based on learnings: In the fgumi EM-seq correct workflow (PR `#166`, fulcrumgenomics/fgumi), both R1 and R2 UMI segments are in forward orientation. Therefore only C→T conversion tolerance is needed; G→A tolerance is not required and was reverted. The `--allow-c-to-t` flag applies uniformly to all UMI segments regardless of read pair index. <details> <summary>🤖 Prompt for AI Agents</summary>Verify each finding against the current code and only fix it if needed.
In
@docs/em-seq-guide.mdaround lines 233 - 246, Add explicit guidance to the
"Step 2: UMI Correction" section to instruct users to include the --allow-c-to-t
flag when their UMI design contains unmethylated cytosines; state that in the
fgumi correct EM-seq workflow (PR#166) both R1 and R2 UMI segments are
forward-oriented so only C→T tolerance is required (no G→A), and make clear that
--allow-c-to-t is applied uniformly across all UMI segments regardless of
read-pair index.</details> </blockquote></details> <details> <summary>crates/fgumi-raw-bam/src/testutil.rs (1)</summary><blockquote> `116-158`: **Guard aux parsing against truncated data.** This helper can panic on malformed/truncated aux buffers due to unchecked indexing/slicing. Prefer early `None` on short buffers for safer test diagnostics. <details> <summary>🛡️ Suggested hardening pattern</summary> ```diff match typ { b'B' => { + if i + 5 > aux.len() { + return None; + } let sub = aux[i]; i += 1; let count = u32::from_le_bytes([aux[i], aux[i + 1], aux[i + 2], aux[i + 3]]) as usize; i += 4; if t == tag && sub == b'C' { + if i + count > aux.len() { + return None; + } return Some(aux[i..i + count].to_vec()); } let elem_size = match sub { b's' | b'S' => 2, b'i' | b'I' | b'f' => 4, _ => 1, }; + if i + count.saturating_mul(elem_size) > aux.len() { + return None; + } i += count * elem_size; }🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed. In `@crates/fgumi-raw-bam/src/testutil.rs` around lines 116 - 158, The function find_u8_array_tag_in_aux can panic on malformed/truncated aux by indexing/slicing without bounds checks; update it to guard every read/slice by checking that aux has enough remaining bytes before each access (e.g., ensure at least 3 bytes before reading t and typ, at least 1+4 bytes before reading sub and count when typ==b'B', verify i + count <= aux.len() before returning aux[i..i+count].to_vec(), and verify i + count*elem_size <= aux.len() before advancing i), return None early on any short buffer, and also guard the Z string loop by ensuring i < aux.len() before each increment; apply these checks inside find_u8_array_tag_in_aux so no unchecked indexing or slicing occurs.
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.
Inline comments:
In `@crates/fgumi-consensus/src/duplex_caller.rs`:
- Around line 1240-1280: The code wrongly assumes ab_consensus is always
top-strand and ba_consensus bottom by passing hardcoded true/false into
build_mm_tag_no_ml and build_mm_ml_tags; instead persist and read the
per-simplex strand orientation on the consensus structs (e.g., add or use a
field like strand_is_top or simplex_strand on
ab_consensus/ba_consensus/consensus) and pass that boolean into
build_mm_tag_no_ml and build_mm_ml_tags when building am/bm and MM/ML tags; also
ensure any logic that copies BA-only consensuses into ab_consensus preserves
that stored strand orientation so tags remain correct for reverse-oriented or
BA-only outputs.
- Around line 924-938: The conversion-artifact detection and handling in
is_conversion_artifact and its branch is incomplete: change the predicate to
check methylation evidence on both reads (a.methylation and b.methylation) for
is_ref_c at index i rather than only a.methylation, e.g. require either side has
evidence.get(i).is_ref_c; when the artifact branch is taken (in the
is_conversion_artifact true path that calls Self::unconverted_base and sums
qualities using Self::cap_quality) also ensure you do not increment error
counters or penalize cE/errors — treat the two observations as agreement for
downstream error accounting (where errors and cE are updated), so update the
code paths that update errors/cE to skip or decrement their penalties when
is_conversion_artifact is true.
In `@crates/fgumi-consensus/src/methylation.rs`:
- Around line 25-34: Change the MethylationEvidence counters unconverted_count
and converted_count from i16 to a wider integer (e.g., i32 or i64) and update
all code that writes into or sums these fields (the simplex-accumulation and
duplex-summation sites that currently assign directly into
unconverted_count/converted_count) to use the wider type; keep clamping to the
i16 range only at the SAM-array emission/serialization boundary so emitted cu/ct
remain safe. Also update any other similar methylation counter fields mentioned
in the review (the other occurrences flagged in the diff) to the same wider type
and adjust any arithmetic, casts, and serialization points accordingly so
accumulation uses the wider type and only the final output is clamped.
In `@crates/fgumi-consensus/src/tags.rs`:
- Around line 119-137: The new methylation tag constants UNCONVERTED_COUNT,
CONVERTED_COUNT, AB_UNCONVERTED_COUNT, AB_CONVERTED_COUNT, BA_UNCONVERTED_COUNT,
and BA_CONVERTED_COUNT are declared but not included in the tag lists used for
reversal/iteration; add the string constants ("cu","ct","au","at","bu","bt") to
both T AFGS_TO_REVERSE and ALL_TAGS (where those arrays/list constants are
defined) so methylation per-base arrays are handled on reverse-strand code paths
and included in full tag iteration. Ensure ordering/duplication matches existing
conventions in TAGS_TO_REVERSE and ALL_TAGS.
In `@crates/fgumi-raw-bam/src/cigar.rs`:
- Around line 15-25: The match over raw_op in the CIGAR opcode decoder currently
treats opcode 6 and any unknown values via the catch-all arm; change it to
explicitly handle 6 => Kind::Pad and replace the fallback arm with a panic (or
unreachable!()) that reports the invalid opcode so malformed BAMs fail fast;
update the module doc comment above the decoder (the doc comment near the CIGAR
parsing function) to state that invalid opcodes will cause a panic consistent
with noodles_compat.rs behavior.
In `@docs/em-seq-guide.md`:
- Around line 26-32: Add a language specifier to the fenced pipeline block to
satisfy markdownlint MD040; change the opening fence from ``` to something like
```text (or ```none) so the block becomes ```text followed by the pipeline lines
and a closing ```; update the fenced block surrounding the Phase 1/Phase 2
diagram accordingly.
In `@src/commands/filter.rs`:
- Around line 1068-1076: When self.min_conversion_fraction is Some, also
validate that self.methylation_mode is set: in the same block that checks range
and self.reference, add a check that self.methylation_mode.is_some() (or
equivalent) and bail with a clear message like "--min-conversion-fraction
requires --methylation-mode" if not present; keep the existing range check for
min_conversion_fraction and the existing reference check so all required flags
(min_conversion_fraction, methylation_mode, reference) are validated together.
In `@src/commands/zipper.rs`:
- Around line 629-643: The loop in the
Kind::Match/SequenceMatch/SequenceMismatch arm indexes ref_bases[ref_offset + i]
with no bounds check; add a guard like if ref_offset + i >= ref_bases.len() {
continue; } before accessing ref_bases to avoid panics at contig edges, keeping
the rest of the logic (comparison to ref_target and writing to seq[read_pos +
i]) unchanged; ensure you still reference the existing symbols ref_bases,
ref_offset, read_pos, seq, converted_base, unconverted_base, modified, and keep
the post-loop increments of read_pos and ref_offset as-is.
---
Nitpick comments:
In `@crates/fgumi-consensus/src/codec_caller.rs`:
- Around line 456-459: The code reparses the CIGAR per read by calling
bam_fields::get_cigar_ops and bam_fields::simplify_cigar_from_raw to build
original_cigar in the hot path; instead, compute and simplify the CIGAR once
earlier and thread that simplified CIGAR through the ClippedRecordInfo struct so
callers reuse it. Add a field (e.g., simplified_original_cigar) to
ClippedRecordInfo, populate it where records are first constructed, remove the
per-call get_cigar_ops/simplify_cigar_from_raw usage in codec_caller.rs (the
original_cigar block), and update any functions that construct ClippedRecordInfo
to pass the precomputed simplified CIGAR.
In `@crates/fgumi-consensus/src/filter.rs`:
- Around line 1697-1700: The match on methylation_mode currently uses a wildcard
for non-Taps variants; change it to explicit arms so each variant is handled
(e.g., match methylation_mode { crate::MethylationMode::Taps => u64::from(cu),
crate::MethylationMode::EmSeq => u64::from(ct), crate::MethylationMode::Disabled
=> u64::from(ct), } ), making the intent clear and preventing future silent
behavior if new variants are added (update the match in the same scope where
numerator is computed).
- Around line 3154-3310: Add a TAPs-mode unit test mirroring the existing EmSeq
tests that verifies the inverted numerator logic for TAPS: create a new test
(e.g., test_conversion_fraction_taps_mode_inverts_counts) that builds a
reference and a read with cu/ct tags similar to the high-conversion and
low-conversion EmSeq cases but flipped for Taps, call
check_conversion_fraction_raw with crate::MethylationMode::Taps, and assert the
opposite pass/fail expectations (i.e., values that passed under EmSeq should
fail under TAPS and vice versa); reference the test helper
check_conversion_fraction_raw and the RecordBuilder/add_i16_array_tag usage to
locate where to mimic setup.
In `@crates/fgumi-consensus/src/vanilla_caller.rs`:
- Around line 1393-1415: When methylation is present we currently compute is_top
using original_raws.first().is_none_or(...), which silently defaults to true for
an empty original_raws; instead explicitly handle the empty case to avoid
masking invariants: inside the methylation block check original_raws.first() and
if it is None either early-return/log an error or assert (e.g., return Err /
debug_assert!) before calling crate::methylation::is_top_strand, then use the
real boolean for is_top when calling crate::methylation::build_mm_ml_tags and
appending MM/ML and cu/ct tags (referencing original_raws, is_top,
crate::methylation::is_top_strand, crate::methylation::build_mm_ml_tags, and
self.options.methylation_mode).
In `@crates/fgumi-raw-bam/src/testutil.rs`:
- Around line 116-158: The function find_u8_array_tag_in_aux can panic on
malformed/truncated aux by indexing/slicing without bounds checks; update it to
guard every read/slice by checking that aux has enough remaining bytes before
each access (e.g., ensure at least 3 bytes before reading t and typ, at least
1+4 bytes before reading sub and count when typ==b'B', verify i + count <=
aux.len() before returning aux[i..i+count].to_vec(), and verify i +
count*elem_size <= aux.len() before advancing i), return None early on any short
buffer, and also guard the Z string loop by ensuring i < aux.len() before each
increment; apply these checks inside find_u8_array_tag_in_aux so no unchecked
indexing or slicing occurs.
In `@docs/em-seq-guide.md`:
- Around line 233-246: Add explicit guidance to the "Step 2: UMI Correction"
section to instruct users to include the --allow-c-to-t flag when their UMI
design contains unmethylated cytosines; state that in the fgumi correct EM-seq
workflow (PR `#166`) both R1 and R2 UMI segments are forward-oriented so only C→T
tolerance is required (no G→A), and make clear that --allow-c-to-t is applied
uniformly across all UMI segments regardless of read-pair index.
In `@src/commands/common.rs`:
- Around line 64-88: The ref_timer currently logs 0 records; update
load_methylation_reference to log the actual contig count by passing
ref_names.len() (cast to the expected numeric type, e.g., as u64 or usize as
required) into ref_timer.log_completion instead of 0 so the timer reports the
number of contigs processed; reference symbols: load_methylation_reference,
ref_timer.log_completion, ref_names.
In `@src/commands/filter.rs`:
- Around line 859-872: The unwrap on methylation_tags
(methylation_tags.as_ref().unwrap()) in the conversion-fraction block is only
safe if upstream validation guarantees methylation_mode is not Disabled when
needs_methylation_tags is true; fix by either enforcing that validation or
adding a local guard: before calling
check_conversion_fraction_raw_with_ref_bases_and_tags (the function used here),
ensure methylation_mode != Disabled and methylation_tags.is_some(), and skip the
call or set pass=false if those conditions are not met; alternatively update the
validation logic referenced above to require a non-Disabled methylation_mode
whenever needs_methylation_tags is set so the unwrap is provably safe.
🪄 Autofix (Beta)
Fix all unresolved CodeRabbit comments on this PR:
- Push a commit to this branch (recommended)
- Create a new PR with the fixes
ℹ️ Review info
⚙️ Run configuration
Configuration used: Organization UI
Review profile: CHILL
Plan: Pro
Run ID: cb9f853b-0428-4424-aef9-2fcca40c2d52
📒 Files selected for processing (19)
crates/fgumi-consensus/src/codec_caller.rscrates/fgumi-consensus/src/duplex_caller.rscrates/fgumi-consensus/src/filter.rscrates/fgumi-consensus/src/lib.rscrates/fgumi-consensus/src/methylation.rscrates/fgumi-consensus/src/tags.rscrates/fgumi-consensus/src/vanilla_caller.rscrates/fgumi-raw-bam/src/builder.rscrates/fgumi-raw-bam/src/cigar.rscrates/fgumi-raw-bam/src/tags.rscrates/fgumi-raw-bam/src/testutil.rsdocs/em-seq-guide.mddocs/taps-guide.mdsrc/commands/common.rssrc/commands/duplex.rssrc/commands/filter.rssrc/commands/simplex.rssrc/commands/zipper.rssrc/lib/reference.rs
ad1ca97 to
73f3636
Compare
Add docs/taps-guide.md covering TAPs pipeline workflows, output tags, filter options, and key differences from EM-Seq. Update em-seq-guide.md for --methylation-mode em-seq syntax.
There was a problem hiding this comment.
Actionable comments posted: 1
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.
Inline comments:
In `@docs/src/guide/methylation.md`:
- Around line 95-98: The fgumi fastq examples are missing the --no-read-suffix
flag, which causes /1 and /2 to be appended to FASTQ read names and later breaks
fgumi zipper's read-name matching; update both example command invocations that
start with "fgumi fastq --input unmapped.bam" to include --no-read-suffix (i.e.,
"fgumi fastq --no-read-suffix --input unmapped.bam ...") so the produced FASTQ
names match the original unmapped BAM when passed through "bwa mem" and "fgumi
zipper".
🪄 Autofix (Beta)
Fix all unresolved CodeRabbit comments on this PR:
- Push a commit to this branch (recommended)
- Create a new PR with the fixes
ℹ️ Review info
⚙️ Run configuration
Configuration used: Organization UI
Review profile: CHILL
Plan: Pro
Run ID: 22baf8c6-4c02-44a4-89eb-25eedab49abd
📒 Files selected for processing (3)
crates/xtask/src/generate_summary.rsdocs/src/guide/methylation.mddocs/theme/sidebar-brand.js
✅ Files skipped from review due to trivial changes (1)
- docs/theme/sidebar-brand.js
Summary
docs/taps-guide.mdcovering TAPs pipeline workflows, output tags, filter options, recommended parameters, and key differences from EM-Seqdocs/em-seq-guide.mdfor--methylation-mode em-seqsyntax--allow-c-to-tis only needed when UMIs contain unmethylated cytosinesStacked on #173.
Test plan