Skip to content

test: Add bcftools pipeline equivalence tests - #334

Open
mwiewior wants to merge 7 commits into
masterfrom
test-bcftools-ops
Open

test: Add bcftools pipeline equivalence tests#334
mwiewior wants to merge 7 commits into
masterfrom
test-bcftools-ops

Conversation

@mwiewior

@mwiewior mwiewior commented Mar 4, 2026

Copy link
Copy Markdown
Collaborator

Summary

  • Add test_bcftools.py with 12 tests verifying pb.sql() replicates a three-step bcftools pipeline:
    1. bcftools view --samples-file (sample subsetting via register_vcf(samples=...))
    2. bcftools filter --exclude 'QUAL<20 || AVG(FORMAT/GQ)<15 || AVG(FORMAT/DP)<15 || AVG(FORMAT/DP)>150' (SQL WHERE clause)
    3. bcftools filter --exclude 'FORMAT/GQ<10 | FORMAT/DP<10 | FORMAT/DP>200' --set-GTs '.' (vcf_set_gts + list_and/list_gte/list_lte)
  • Include precomputed bcftools output (100 variants × 2000 samples) for regression comparison
  • Verify row counts, positions, GT masking counts, per-sample GT/DP/GQ values, and sink_vcf roundtrip

Test plan

  • All 12 tests pass locally (pytest tests/test_bcftools.py -v)
  • CI passes

🤖 Generated with Claude Code

Comment thread tests/test_bcftools.py
assert result_df.filter(pl.col("qual") < 20).shape[0] == 0

def test_avg_gq_filter(self, result_df):
"""All variants have avg(GQ) >= 15."""

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

tutaj można by dodać assert walidujący że przed filtrowaniem liczba takich wierszy była większa od zera

Comment thread tests/test_bcftools.py
FILTER_SQL = """\
SELECT chrom, start, "end", id, ref, alt, qual, filter, \
named_struct(\
'GT', vcf_set_gts(genotypes."GT", \

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

w którym miejscu tutaj ustawiamy GT na missing ./. ?

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

spojrzałem do https://github.com/biodatageeks/datafusion-bio-formats/pull/87/changes#diff-445df8bbc29de5d41e0d0def0dab28ef2d17eeafe0d8de99359f22c41ff27081 i myślę, że niezbyt dobrze że ./. jest domyślnym replacement w vcf_set_gts - to nie jest czytelne ani oczywiste. Lepiej byłoby chyba usunąc opcję dwuargumentową dla vcf_set_gts i pozostawić vcf_set_gts(gt_list, mask, replacement) ze wszystkimi argumentami wymaganymi

Comment thread tests/test_bcftools.py
format_fields=["GT", "DP", "GQ"],
)
roundtrip_df = pb.sql("SELECT * FROM bcftools_roundtrip").collect()
assert roundtrip_df.shape[0] == 100

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

dobrze byłoby porównać wyprodukowane pliki wynikowe po dekompresji. Obecnie wydaje mi się że polars-bio wyprodukuje troche inny plik, m.in. header będzie obcięty. Ale już pomijając header przede wszystkim dobrze byłoby sprawdzić czy wiersze są takie same

@mwiewior
mwiewior force-pushed the test-bcftools-ops branch 2 times, most recently from c34a301 to d4240a7 Compare March 6, 2026 14:18
mwiewior and others added 6 commits March 6, 2026 15:20
Verify that pb.sql() with VCF UDFs (vcf_set_gts, list_avg, list_gte,
list_lte, list_and) replicates a three-step bcftools pipeline:
sample subsetting, quality filtering, and per-sample GT masking.

Includes precomputed bcftools output for regression comparison.

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
Verify that GT/DP/GQ values survive the full pb.sql → sink_vcf → re-read
roundtrip and match the precomputed bcftools output.

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
- register_vcf gains format_fields and samples parameters
- VCF UDFs: list_avg, list_gte, list_lte, list_and, vcf_set_gts,
  vcf_an, vcf_ac, vcf_af
- register_vcf_long_view for per-sample SQL queries
- Updated Rust bindings, write path, and metadata extractors
- Tests for format columns, UDFs, and write roundtrips

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
…fix)

Verify parsing of VCF with Number=. AD field, mixed FORMAT fields
across rows (GT:AB:AD:DP:GQ:PL vs GT:AD:DP:GQ:PL), and haploid
GTs on chrX. Uses existing head_106667_tail_6.vcf test data.

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
- Add vcf_process_genotypes Rust UDF for single-pass genotype processing
  (GT normalization, quality masking, PL correction, dosage calculation)
- Add list_normalize_gt Rust UDF for haploid GT normalization
- Register UDFs in DataFusion session context
- Add benchmark scripts comparing 3 approaches:
  - Explode/reaggregate (baseline)
  - Polars list ops + map_batches
  - Pure SQL with Rust UDFs (5.1x faster, 3.2x less memory)

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
@mwiewior
mwiewior force-pushed the test-bcftools-ops branch from d4240a7 to 26759a2 Compare March 6, 2026 14:23
Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
output_files = {}

for partitions in [1, 2, 4, 8, 12]:
out_path = tempfile.mktemp(suffix=f"_{partitions}p.vcf.gz")

Check failure

Code scanning / CodeQL

Insecure temporary file High

Call to deprecated function tempfile.mktemp may be insecure.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants