fix(dist): clamp the graph's reference window at a contig end (#189) - #201
Merged
Conversation
The alignment graph's trailing node took its sequence from a substr that clamps at the contig end, but its end coordinate from an unclamped ref_end + 1. At a contig end the two disagreed: on a 10 bp contig with a SNP on the final base the trailing node spanned [2, 4) while holding zero bases, and graph->ref was two characters that the coordinate 4 overran. Clamp the window once, in win_end, and use it for the trailing query node, the trailing truth node, this->truth and this->ref. Every node's coordinate span now equals the length of the sequence it holds, with zero-width insertion variant and bypass nodes the deliberate exception. At a contig end the trailing node becomes a zero-width sink, the mirror of the zero-width entry node #177 added at position 0. No behavior change: nothing reads the trailing node's end coordinate. Both alignment endpoints are built from qseqs/tseqs lengths rather than coordinates, edit distances go through the already-clamped graph->ref and graph->truth, the trailing node is the sink so its end coordinate is never matched against any node's begin coordinate, and the insertion-leap rule additionally requires tidxs >= 0 || tskips >= 0, which the trailing node never satisfies. Its sequence was already empty at a contig end, so the endpoint cell was already the sentinel. A chr20 before/after run is byte-identical across every scored output file, which covers the loop restructure; the fixture's last variant sits 158 kb from the contig end, so the clamp itself is covered by the new tests rather than by it. Retire the TODO asking why the +1 is needed "compared to generate_ptrs_strs()". That function was removed in 232036e and reached exactly as far, via `for (int ref_pos = beg_pos; ref_pos <= end_pos; )` and `ref_end = end_pos+1`, so there was never an asymmetry between them to explain. The remaining +1 is the window's right flank and is left alone: unlike the coordinate, the trailing node's sequence is read, so shortening it would be a scoring change. Add the coverage the issue asks for, mirroring #177's: five GraphCtor tests over a variant on a contig's final base and on the second-to-last (all five fail without the clamp), three PrecRecall tests pinning TP, FN and FP there, and a contig_end_snp integration scenario placing a SNP on the final base of sc1. Refs #189
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Note
Authorship: the content below was drafted by Claude Opus 5 (an AI coding agent) and
filed via
ghunder @TimD1-bot, a bot account operated by @TimD1. It reflects theagent's analysis, not a statement authored by @TimD1.
What
The alignment graph's trailing node took its sequence from a
substrthat clamps at the contig end, but its end coordinate from an unclampedref_end + 1. At a contig end the two disagreed. On a 10 bp contig with a SNP on the final base, the trailing node spanned[2, 4)while holding zero bases, andgraph->refwas two characters that the coordinate4overran.Clamp the window once, in
win_end, and use it for the trailing query node, the trailing truth node,this->truth, andthis->ref. Every node's coordinate span now equals the length of the sequence it holds, with zero-width insertion variant and bypass nodes the deliberate exception. At a contig end the trailing node becomes a zero-width sink -- the mirror of the zero-width entry node #177 added at position 0.Why this is behavior-neutral
Nothing reads the trailing node's end coordinate:
src/dist.cpp:178-180andsrc/dist.cpp:305-307useqseqs[qnodes-1].length()-1. Edit distances go throughgraph->ref/graph->truth(src/dist.cpp:309-310,src/dist.cpp:361-362), which were alreadysubstr-clamped.qends[n2] == qbegs[n1], and no node begins after the last one.tidxs >= 0 || tskips >= 0(src/dist.cpp:1066-1069), which the trailing node never satisfies, so its width is irrelevant there.Its sequence was already empty at a contig end, so the endpoint cell was already the sentinel. The clamp changes the coordinate and nothing else.
Verification
precision-recall{,-summary}.tsv,query.tsv,truth.tsv,summary.vcf,switchflips.tsv,phase-blocks.tsv,phasing-summary.tsv,genotype-errors.tsv). Onlyruntime.tsvand stderr timings differ.The retired
TODOThe
TODOasked why the+1is needed "compared togenerate_ptrs_strs()". That function was removed in232036e, and read from that commit's parent it reached exactly as far, viafor (int ref_pos = beg_pos; ref_pos <= end_pos; )andref_end = end_pos+1. There was never an asymmetry between them to explain.The remaining
+1is the window's right flank and is left alone -- unlike the coordinate, the trailing node's sequence is read, so shortening it would be a scoring change. That asymmetry (one base of left flank, two on the right) is out of scope here.Tests
Mirroring #177's coverage:
GraphCtortests over a variant on a contig's final base and on the second-to-last. All five fail without the clamp -- verified by revertingsrc/dist.cppand rebuilding.PrecRecalltests pinning TP, FN, and FP on the final base. These pass with and without the clamp, which is the issue's own point that labeling was already correct -- they close the coverage gap rather than guard the fix.contig_end_snpintegration scenario placing a SNP on the final base ofsc1(with aPOS 350mid-contig control), which does exercise the clamp end-to-end.Closes #189