Skip to content

Commit debf062

Browse files
committed
test_read_phasing: add synthetic somatic+germline phasing test
Builds a synthetic mini-BAM (via MockAlignmentFile + make_pysam_read) with reads spanning an imagined somatic SNV at chr1:4 (T>G) and a nearby germline SNP at chr1:8 (T>A): three reads carry both alts, two carry the somatic alt only, two carry the germline alt only. Walks the full ReadCollector -> annotate_phased_variants -> IsovarReadPhasing stack and verifies the somatic and germline calls are reported as in-cis partners. This is the canonical motivating scenario for openvax/varcode#268's germline-aware effect prediction: a somatic mutation whose codon overlaps a germline polymorphism, where the protein outcome depends on whether the two are in cis -- and only RNA-read co-occurrence disambiguates.
1 parent c252810 commit debf062

1 file changed

Lines changed: 106 additions & 0 deletions

File tree

tests/test_read_phasing.py

Lines changed: 106 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -20,9 +20,12 @@
2020
from isovar.allele_read import AlleleRead
2121
from isovar.isovar_result import IsovarResult
2222
from isovar.phasing import annotate_phased_variants
23+
from isovar.read_collector import ReadCollector
2324
from isovar.read_evidence import ReadEvidence
2425

2526
from .common import eq_
27+
from .genomes_for_testing import grcm38
28+
from .mock_objects import MockAlignmentFile, make_pysam_read
2629
from .testing_helpers import data_path
2730

2831

@@ -245,3 +248,106 @@ def test_adapter_respects_min_shared_fragments_threshold():
245248
# Threshold filters out the single shared read, so no partners.
246249
eq_(phasing.partners_in_cis(v1), ())
247250
eq_(phasing.partners_in_cis(v2), ())
251+
252+
253+
# ---- synthetic somatic + germline phasing scenario -----------------------
254+
#
255+
# Builds a synthetic mini-BAM in-memory and walks the full
256+
# ReadCollector -> annotate_phased_variants -> IsovarReadPhasing
257+
# stack. The scenario mirrors the canonical motivating use case for
258+
# openvax/isovar#182 / openvax/varcode#268: a somatic SNV whose codon
259+
# overlaps a germline SNP, where reading direction alone can't tell
260+
# whether the two variants are in cis or trans -- only RNA-read
261+
# co-occurrence can.
262+
#
263+
# Imagined locus (mouse, b16-adjacent): chr1:1-10 reference is
264+
# `ACCTTGATCG`. A somatic SNV sits at chr1:4 (T>G), and an imagined
265+
# germline SNP sits at chr1:8 (T>A). The reads below are the kind of
266+
# RNA fragments Isovar's read collector would harvest from a tumor
267+
# BAM at this locus.
268+
#
269+
# ref : A C C T T G A T C G
270+
# pos (1-based) 1 2 3 4 5 6 7 8 9 10
271+
# somatic@4 : ^ T>G
272+
# germline@8: ^ T>A
273+
#
274+
# Three reads carry both alts (cis evidence), two carry only the
275+
# somatic alt, two carry only the germline alt. After phasing with
276+
# the default min_shared_fragments threshold of 2, the somatic and
277+
# germline calls should be reported as partners.
278+
279+
280+
def test_synthetic_somatic_and_germline_phased_on_rna_reads():
281+
chromosome = "1"
282+
somatic = Variant(
283+
chromosome, 4, "T", "G", grcm38, normalize_contig_names=False)
284+
germline = Variant(
285+
chromosome, 8, "T", "A", grcm38, normalize_contig_names=False)
286+
287+
# Reads cover positions 1-10 (reference_start=0, 10M cigar). The
288+
# MD tag encodes mismatches relative to the assumed reference
289+
# `ACCTTGATCG`:
290+
# * both alts: ACCGTGAACG -> mismatches at 0-based 3 and 7 -> MD=3T3T2
291+
# * somatic only: ACCGTGATCG -> mismatch at 0-based 3 -> MD=3T6
292+
# * germline only: ACCTTGAACG -> mismatch at 0-based 7 -> MD=7T2
293+
def both_alts(name):
294+
return make_pysam_read(
295+
seq="ACCGTGAACG", cigar="10M", mdtag="3T3T2",
296+
name=name, reference_start=0)
297+
298+
def somatic_only(name):
299+
return make_pysam_read(
300+
seq="ACCGTGATCG", cigar="10M", mdtag="3T6",
301+
name=name, reference_start=0)
302+
303+
def germline_only(name):
304+
return make_pysam_read(
305+
seq="ACCTTGAACG", cigar="10M", mdtag="7T2",
306+
name=name, reference_start=0)
307+
308+
reads = [
309+
both_alts("frag-both-1"),
310+
both_alts("frag-both-2"),
311+
both_alts("frag-both-3"),
312+
somatic_only("frag-somatic-1"),
313+
somatic_only("frag-somatic-2"),
314+
germline_only("frag-germline-1"),
315+
germline_only("frag-germline-2"),
316+
]
317+
samfile = MockAlignmentFile(references=(chromosome,), reads=reads)
318+
collector = ReadCollector()
319+
320+
somatic_evidence = collector.read_evidence_for_variant(
321+
variant=somatic, alignment_file=samfile)
322+
germline_evidence = collector.read_evidence_for_variant(
323+
variant=germline, alignment_file=samfile)
324+
325+
# Sanity-check the synthetic-data construction itself before
326+
# exercising the adapter -- if the read collector doesn't pick
327+
# the alts up as expected, downstream assertions would be
328+
# misleading.
329+
eq_(
330+
somatic_evidence.alt_read_names,
331+
{"frag-both-1", "frag-both-2", "frag-both-3",
332+
"frag-somatic-1", "frag-somatic-2"})
333+
eq_(
334+
germline_evidence.alt_read_names,
335+
{"frag-both-1", "frag-both-2", "frag-both-3",
336+
"frag-germline-1", "frag-germline-2"})
337+
338+
annotated = annotate_phased_variants([
339+
IsovarResult(
340+
variant=somatic,
341+
read_evidence=somatic_evidence,
342+
predicted_effect=None),
343+
IsovarResult(
344+
variant=germline,
345+
read_evidence=germline_evidence,
346+
predicted_effect=None),
347+
])
348+
phasing = IsovarReadPhasing(annotated)
349+
350+
assert phasing.has_evidence(somatic)
351+
assert phasing.has_evidence(germline)
352+
eq_(phasing.partners_in_cis(somatic), (germline,))
353+
eq_(phasing.partners_in_cis(germline), (somatic,))

0 commit comments

Comments
 (0)