Introduction to Bioinformatics and Computational Genomics

Week 6Variant calling and structural variants

Once you have reads aligned to a reference, the question becomes which differences are real. This is the unit where the tools you run are the actual tools the field uses.

Questions this week answers

  • Why is a variant call a probability rather than a lookup?
  • Your mean coverage is 20x. Is that enough to call a heterozygote?
  • Every filter in a VCF exists because of a specific artefact. Which one catches what?

By the end of this week you can

  • Explain why variant calling is a probability, not a lookup
  • Read a VCF and know what each field is claiming
  • Apply the standard filters and say what artefact each one catches
  • Distinguish the four ways structural variants show up in read data
0 of 4 done

  1. Mean coverage is an average, and variants are called at the tail

    Set the mean to 20x. Then look at how much of the genome is still below the depth you would need to trust a call.

    A mean of 10.0x leaves 2.9% of the genome below 5x, or about 87.8 million bases. Depth at a base is Poisson around the mean, so the tail never disappears, it only thins.

    Mean coverage10x30.0 Gb of reads
    Genome below 5x2.9%87.8 million bases
    Nothing at all4.5e-3%no read covers it
    Heterozygote seen at 5x93.8%both alleles sampled
    A bacterium is a few Mb; the human genome is about 3,000.
    Below this you would not report a genotype.
    Things to try0 of 3

    The 150 bp x 400 M over 3 Gb example is from the week 4 sequencing lecture.

What should survive this week

  • You never observe a genotype, only reads, so what you want is Pr(G | R). That needs a prior, and transitions outnumber transversions roughly two to one.
  • Coverage is Poisson around its mean, so the low tail thins but never disappears. A healthy average and millions of uncallable bases coexist comfortably.
  • Missing a true heterozygote requires every read to draw the same allele, which is 2 x (1/2)^d. Eight reads gets you past 99%.
  • Strand bias, allele balance and Hardy-Weinberg each catch a different artefact. VCF is tab-separated text, so grep and awk work on it like anything else.