Skip to content

Synthetic expansion ladder: at what size does sizing break, and does it fail loudly? #34

Description

@wdecoster

STRdust's worst failure mode is reporting a real expansion short rather than not at all: a no-call is visible and prompts a second look, a confidently wrong small number does not. Nothing has ever tested for it directly, and the HG002 benchmark structurally cannot.

The mechanism

find_insertions keeps only read.is_primary. minimap2's map_ont defaults (zdrop=400, max_gap=5000) can split a sufficiently large insertion across alignments — the non-primary piece is then discarded and the allele comes back truncated, not missing. flanking = 5000 bounds it from the other side, since the repeat-compressed reference only carries that much context.

Why the truth-set benchmark cannot answer it

On 5,072 loci carrying a real expansion over 200 bp, detection was 92.2% found with 1.0% reported "short". But that measures the expansions GIAB HG002 happens to contain. It cannot distinguish "the truth set has few enormous alleles" from "we size enormous alleles wrongly", because the loci that would discriminate are the ones the truth set lacks.

The experiment

A synthetic ladder: reads carrying insertions of 100 / 500 / 1000 / 2000 / 4000 / 8000 bp at a single locus, genotyped with defaults, in both modes.

Read out, per rung: is the allele sized correctly, sized short, or no-called? The answer wanted is where sizing breaks, and whether it fails loudly or quietly. A clean transition from correct to no-call is acceptable behaviour; a region where alleles come back plausibly-but-wrongly short is not.

Worth also varying read length (10 kb vs 20 kb) to separate "the read does not span it" from "the aligner split it".

Background: KNOWN_ISSUES §27.4 and §30.1 as of 3fba652.

Activity

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions