Skip to content

repair silently chooses one of multiple overlapping sequence placements #679

Description

@SuhasSrinivasan

repair silently chooses one of multiple overlapping sequence placements

Summary

modkit repair checks whether the acceptor sequence occurs uniquely within the donor sequence before projecting MM/ML calls. The current search misses overlapping occurrences, so an acceptor such as ACAC is treated as uniquely placed in donor ACACAC even though it starts at offsets 0 and 2.

The command then silently chooses the first placement and can project modification calls to the wrong acceptor positions.

Severity

Severity: Medium — scientific correctness

Rationale: The command accepts the input, reports the record as repaired, and emits plausible MM/ML tags derived from an arbitrary placement. This silently changes modification coordinates for the affected record. The trigger is limited to repetitive donor sequence in which the complete acceptor occurs at multiple, overlapping starts, so exposure is data-dependent and record-local; the safe correction is to reject that record rather than guess a placement.

User and scientific impact

  • Affected result or workflow: MM/ML/MN tags produced by repair for trimmed or otherwise changed reads.
  • Direction of error: modification probabilities can be assigned to the wrong positions within the acceptor read.
  • Likely exposure: data-dependent, particularly for short acceptor reads and repetitive or low-complexity sequence.
  • Detectability or workaround: users must independently search the donor sequence for overlapping occurrences. A normal substring search that reports only non-overlapping matches is not a sufficient check.

Affected versions and environment

  • Released version: modkit 0.6.4.
  • Development revision: 5cecc3fb3a9336068d9e3c68d5c08d678153dd2c.
  • Operating system and architecture used for verification: macOS on Apple Silicon.
  • Input format: query-name-sorted donor and acceptor BAM.
  • Related tool used to generate and inspect the compact fixture: samtools 1.23.1.

Steps to reproduce

Create one donor whose modified second A is compatible with two overlapping acceptor placements:

cat > donor.sam <<'EOF'
@HD	VN:1.6	SO:queryname	SS:queryname:lexicographical
repeat	4	*	0	0	*	*	0	0	ACACAC	IIIIII	MM:Z:A+a.,1;	ML:B:C,200
EOF

cat > acceptor.sam <<'EOF'
@HD	VN:1.6	SO:queryname	SS:queryname:lexicographical
repeat	4	*	0	0	*	*	0	0	ACAC	IIII
EOF

samtools view -b -o donor.bam donor.sam
samtools view -b -o acceptor.bam acceptor.sam

modkit repair \
  --donor donor.bam \
  --acceptor acceptor.bam \
  --output repaired.bam \
  --threads 4

samtools view repaired.bam

Observed behavior

modkit 0.6.4 exits successfully and reports one repaired record:

repeat  4  *  0  0  *  *  0  0  ACAC  IIII  MM:Z:A+a.,1;  ML:B:C,200  MN:i:4

The MM:Z:A+a.,1; result corresponds to choosing the occurrence at donor offset 0. Choosing the equally valid occurrence at offset 2 would place the same donor modification at the first acceptor A instead.

Control or independent oracle

ACAC occurs in ACACAC at zero-based offsets 0 and 2. Because the two placements imply different projected coordinates, neither placement is uniquely justified by sequence identity alone.

A corrected implementation reports zero repaired and one failed record for this fixture, emits no repeat alignment, and continues processing unrelated records.

Expected behavior

  • Detect overlapping as well as non-overlapping acceptor occurrences in the donor's forward sequence.
  • Treat more than one occurrence as a per-record repair rejection rather than selecting one silently.
  • Continue processing other records according to the existing nonfatal record-rejection policy.
  • Preserve exact MM/ML/MN behavior for unique forward and reverse-strand placements.

Root-cause evidence

  • Relevant code: repair_tags.rs:327-332.
  • str::match_indices returns non-overlapping matches. After accepting the match at offset 0, it does not report the overlapping match at offset 2.
  • The focused parent regression accepts ACACAC/ACAC; the corrected search rejects it while preserving unique, absent, non-overlapping-duplicate, one-sided empty-sequence, and reverse-strand controls.

Proposed fix scope

Find the first occurrence, then search again beginning one byte after that start. Any second occurrence is ambiguous, including one that overlaps the first. This requires no collection of every match and keeps the existing unique-placement projection unchanged.

Non-goals

  • No heuristic choice among multiple placements and no alignment-based disambiguation.
  • No change to query-name merging, parallel output ordering, fatal I/O handling, or atomic publication.
  • No change to the policy for duplicate-flagged donor records.
  • No change to the existing per-record rejection contract.

Acceptance criteria

  • ACAC in ACACAC and equivalent reverse-strand cases are rejected as ambiguous.
  • Separate non-overlapping repeated occurrences remain rejected.
  • A unique occurrence produces byte-equivalent MM/ML/MN tags.
  • An absent occurrence and the existing one-sided empty-sequence controls return bounded per-record errors rather than producing a repaired record.
  • Focused repair tests and the applicable workspace suite pass.

Reproduction artifacts

The complete minimal SAM fixtures are generated inline above. With samtools 1.23.1, the donor BAM was 294 bytes with SHA-256 3343bc265731d6f17536d47c2470295d156f337ef87a6c8983d39aaa15836802; the acceptor BAM was 275 bytes with SHA-256 13625fe3ae3f071a40403956f35ad5e153268668cd497446889dcf66e83702c7.

Related work

  • Query-name merge ordering and fatal/atomic output lifecycle are tracked separately.
  • Proposed PR: one small standalone repair, to be linked after publication.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions