feat: add 1-indel search via bidirectional DFS in Rust#24
Conversation
|
Really nice work, Roman. It looks good for the most part. I traced both test cases by hand and the logic is correct for 1 indel. Not mergeable yet though. Here is what needs to happen before it lands, roughly in priority order. 1. Rebase onto current master.
2. Remove 3. Give indels their own column. 4. Strengthen the tests + add a brute-force oracle in tests. 5. k-selection should log loudly. 6. CLI + docs. |
Implements insertion/deletion matching as a new search path in the Rust engine. Uses the pigeonhole principle for seed selection, bidirectional DFS with full budget splitting, and two-level deduplication across seeds. Adds max_indels parameter to Matcher with mutual exclusion against max_mismatches. max_indels > 1 raises ValueError pending validation.
Tests verified by hand against Q8V336 (NCAP_DUGBA) in the test proteome.
QNALVEATRFC is constructed so the only match in the test proteome requires deleting Q at query position 0 (terminal). Verified that removing the guard produces a spurious NALVEATRFC hit; the guard correctly blocks it. Also adds *.pepidx to .gitignore to prevent index files from being staged.
Deletions at query position 0 or query_len-1 are now permitted when the protein position is not at a sequence boundary, allowing valid biological indel matches that were previously blocked.
Roman-Young
force-pushed
the
feature/indel-matching
branch
from
b90ce8a to
09a6757
Compare
Not imported anywhere in the package; relocating per PR IEDB#24 review (item 2) since its role is validating the Rust indel implementation, not shipping with it.
Previously the indel count rode in the Mismatches slot. Split it into its own column (right after Mismatches in FINAL_COLUMNS) so indel rows report Mismatches=0/Indels=N and non-indel rows report Indels=0, without disturbing existing Mismatches semantics for exact/mismatch/ discontinuous search. Per PR IEDB#24 review item 3.
is_terminal_deletion treated p_idx == 0 and p_idx == protein_len - 1 as boundary artifacts, but both are real, existing residues, not out-of-bounds. This required 2 residues of protein context around a deletion's gap when 1 is sufficient (a real adjacent residue already confirms the query residue is genuinely absent, not just missing because the recorded sequence happens to end there). Only strictly out-of-bounds references (p_idx < 0 or p_idx >= protein_len) are actually ambiguous.
- Replace indel.py's seed/DFS port with a genuinely independent brute-force oracle (direct window enumeration, no seeding/recursion) so it can catch bugs the DFS and its Python port would share - Add ValueError path tests for max_indels>1 and indel+mismatch combos - Add peptide_len < k edge case test - Add multi-hit test covering two different proteins in one query - Add randomized property test comparing Rust output against the brute-force oracle across 25 seeded random peptide/proteome pairs
Silently kicking off a full-proteome preprocessing pass with a terse message is a surprise for large proteomes. Per PR IEDB#24 review item 5.
Per PR IEDB#24 review item 6.
Trailing constraint note (after the code block) rather than folded into the intro sentence, mirroring how Counts-Only Mode documents its best_match incompatibility.
…ts_only Both combinations previously produced silent wrong behavior: best_match was ignored (indel_search() results were never passed through _best_match_filter), and counts_only silently returned exact-match counts (max_mismatches=0) rather than indel counts. Neither is implemented yet, so both now fail loudly like the existing max_indels/max_mismatches mutual-exclusivity check.
_with_deletion/_with_insertion are named for the edit applied to the query string, which is the opposite of the Indels label the search engine reports for the resulting match (removing a query residue means the protein has extra content there, an insertion match, and vice versa). Noted during audit to prevent misreading.
A deletion at the query's own first or last residue has no query-side context to confirm the residue is genuinely absent, as opposed to the alignment simply starting/ending one residue short -- so it's now blocked regardless of protein-side buffer, reverting the earlier relaxation. Applied identically to match.rs's is_terminal_deletion and the brute-force oracle, so both stay in lockstep. Added test_terminal_insertion_blocked and test_insertion_at_second_to_last_position_found to cover the analogous insertion case: a padded leading/trailing insertion is unverifiable and blocked, while an interior insertion one position before the query's own boundary (e.g. LPDGVWEESS) has real query-side context on both sides and remains valid. Renamed indel.py to indel_brute_force.py for clarity ahead of future 2-indel brute-force work.
dmx2
left a comment
dmx2
left a comment
There was a problem hiding this comment.
Very good Roman, thank you. We are so close to merging. One extremely minor point and then one change in pepmatch/matcher.py and we're good to go.
-
Fix the PR description
PR summary still says “deletions at query position 0 or query_len-1 are now permitted” but we decided on this already. So just fix the description in case anyone is following this repo and they don't get confused on what we changed. -
Separate
MismatchesandIndelscolumns in each mode
Right now, a mismatch search returns an Indels column (always 0) and an indel search returns a Mismatches column (always 0). I think that will confuse users, as they will see an Indels column on a normal mismatch run and wonder if it also searched for indels.match()method concatenates the linear result with the discontinuous result for mixed queries but I believe this is disabled for indel searches (we do not do discontinuous indel searches). Check for that and the easiest would be to leave mismatching and discontinuous functions alone, but just replaced "Mismatches" with "Indels" and populate it correctly. This should be easy since the column will stay in one place and alternate for either mode.
I'll do one last review and then merge if these changes are made. Thanks!!!
| FINAL_COLUMNS = [ | ||
| 'Query ID','Query Sequence','Matched Sequence','Protein ID','Protein Name','Species', | ||
| 'Taxon ID','Gene','Mismatches','Mutated Positions','Index start','Index end', | ||
| 'Taxon ID','Gene','Mismatches','Indels','Mutated Positions','Index start','Index end', |
There was a problem hiding this comment.
I know this might be annoying for the code readability, but we should split the output when doing indel vs mismatching mode. See my final comments.
| ] | ||
|
|
||
| def _to_dataframe(self, cols): | ||
| def _to_dataframe(self, cols, is_indels=False): |
There was a problem hiding this comment.
Yeah, we already have a flag here, so popping the indel column in should be easy.
…ntinuous Each mode now carries a single edit-count column instead of always showing both: Indels for indel search, Mismatches for every other mode. The prior schema surfaced an all-zero Indels column on mismatch runs (and an all-zero Mismatches column on indel runs), which reads as a search that did not happen; dropping the twin also restores the original schema for the mismatch, exact, discontinuous, and best_match outputs. _to_dataframe now names the single column by mode via _final_columns(is_indels), replacing the fixed FINAL_COLUMNS list. Since a single match() call can concat a linear result frame with a discontinuous one, and the two must share columns, max_indels combined with discontinuous epitopes is now rejected in __init__ -- a discontinuous indel search is undefined anyway (a position-anchored epitope has no contiguous window to seed or extend). Adds test_indel_mode_emits_indels_column_only, test_mismatch_mode_emits_mismatches_column_only, and test_max_indels_and_discontinuous_raises.
|
Merged! |
2 participants
Summary
max_indels=1search in the Rust engine via exhaustive bidirectional DFS anchored on pigeonhole-guaranteed k-mer seedsrun_indelinmatch.rswith theindel_search_peptide→extend_bidirectional→dfscall chainIndelsfor indel search,Mismatchesfor every other mode — never both; indel rows leaveMutated Positionsempty-i/--max_indelsCLI flag and a README "Indel Searching" sectionKey design decisions
Indelscolumn on a mismatch run (and vice versa);max_indelsis rejected together with discontinuous epitopes, which have no contiguous window to alignTest plan
test_indel_search: 1-insertion match (NALVEARFC → NALVEATRFC) and 1-deletion match (EQLVPSYQQA → EQLVSYQQA) against Q8V336test_query_terminal_deletion_blocked: confirms QNALVEATRFC does not match NALVEATRFC (deleting Q at query position 0 is a barred query-terminal deletion)ValueErrorpaths,peptide_len < k, multi-hit across proteins, and per-mode column teststests/indel_brute_force.py)