Skip to content

known_differences.md

Latest commit

 

History

History
133 lines (92 loc) · 4.94 KB

known_differences.md

File metadata and controls

133 lines (92 loc) · 4.94 KB
Raw

Known Differences Between vcfkit and bcftools

This document tracks intentional behavioral differences between vcfkit normalize and bcftools norm. Each entry notes the version where the limitation was introduced, a concrete example, and the planned fix.

1. Multi-allelic indel left-alignment (v0.1.x — planned fix in v0.2)

Behavior

vcfkit normalize passes multi-allelic indel records through unchanged when --no-split is used (or when splitting is disabled). It does not attempt to left-align multi-allelic indels.

bcftools norm (without -m) will left-align multi-allelic indels jointly — i.e. it applies the Tan 2015 algorithm treating all ALTs together as a group, choosing the leftmost position where all alleles remain valid.

Concrete Example

Reference (chr22, hg19 b37) around position 16404838 is a poly-A run:

...GGGAAAAAAAAAAAAAAAA...
         ^16404838

Input VCF record that is not fully left-aligned (shifted right by 2):

22  16404840  .  AA  AAA,A  100  PASS  DP=100

bcftools norm output — shifts the record left by 2 into the poly-A run:

22  16404838  .  GA  GAA,G  100  PASS  DP=100

vcfkit normalize --no-split output — record passes through unchanged:

22  16404840  .  AA  AAA,A  100  PASS  DP=100

Note: the 1000 Genomes Project chr22 data contains this exact record at position 16404838 (GA -> GAA,G), which is already fully left-aligned. The divergence only manifests when the input record is not yet at the leftmost valid position.

Adversarial test

tests/corpus/synthetic/multiallelic_polyA.vcf contains a multi-allelic indel at position 9 of a poly-A tract (ref NNNNGAAAAAAAAAAAGCNNNN, both ALTs shiftable leftward). Tests:

  • normalize_test::multiallelic_polya_fixture_passes_through_unchanged — always runs; confirms vcfkit passes the record through unchanged.
  • normalize_test::diff_multiallelic_polya_matches_bcftools_no_split — #[ignore], requires bcftools; confirms vcfkit and bcftools norm (without -m) produce the same output on this adversarial input.

Root Cause

The biallelic Tan 2015 algorithm in vcfkit-core/src/normalize.rs operates on a single (REF, ALT) pair. Applying it to only the first ALT of a multi-allelic record would silently drop the remaining ALTs — a data-corrupting bug fixed in v0.1.4 by skipping left-alignment for multi-allelic records entirely.

Proper joint multi-allelic left-alignment requires selecting a single anchor position where all ALT alleles can be simultaneously represented, which is a non-trivial extension of the algorithm.

Planned Fix

v0.2 will implement joint Tan 2015 left-alignment across all ALTs: find the leftmost position P such that trim_and_extend(REF, ALT_k, P) is valid for every k simultaneously, then rewrite the record at P.

2. --ask low-confidence translations (v0.3.x)

Behavior

When vcfkit filter --ask produces a translation with confidence below 50%, the --yes (non-interactive) mode exits non-zero with an error. This prevents scripts from silently running ambiguous queries.

Add --accept-low-confidence to override:

vcfkit filter -a "missense variants in BRCA1" --yes --accept-low-confidence input.vcf

Low confidence typically indicates the query references fields not present in the VCF header (e.g. CSQ when the file has no VEP annotation), or is phrased in a way the model cannot translate with high certainty.

No equivalent in bcftools

bcftools does not have a natural-language filter mode; this difference is vcfkit-only.

3. --ask: schema-based grounding (v0.3.0-alpha.2+)

Translations use the VCF header schema (field names, types, descriptions) but not example values from the data. Two failure modes result:

String value mismatch. A translation like INFO/varType == 'SNP' may return zero records if the field is defined in the header but not populated in any record, or if the file uses different casing or naming (e.g. 'snp', 'SNV', 'snv'). The expression is syntactically valid and grounded in the schema; only the value guess is wrong. Confirmed on GiAB HG001: varType is in the header but absent from all records.

Data sparsity. A query like "variants with QUAL above 100" returns zero records if the VCF uses a fixed QUAL value ≤ 100. GiAB HG001 sets QUAL = 50 for all records; QUAL > 100 correctly returns nothing. Similarly, a truth-set VCF may have all records with FILTER = PASS, making FILTER != 'PASS' return zero.

Workaround: inspect example values before trusting low-record-count results. Use -e with a hand-written expression once you know the actual data values. A future release may surface example values to the LLM.

4. No other known differences

The differential tests in tests/normalize_test.rs compare vcfkit output against bcftools norm on the synthetic corpus (basic SNPs, multi-allelic SNPs, unnormalized indels) and all pass. No additional behavioral differences have been identified.