Thirty Times Nothing

Zetobit · Bioinformatics Insight Series · No. 12

Pipeline Design

Thirty Times Nothing

Mean coverage is an average, and an average is a summary of a distribution nobody looks at. The 30× on your sequencing quote says almost nothing about whether the base you care about was ever read.

Kanna Nandakumar, PhD Founder & CEO, Zetobit LLC ~1,900 words

0.0003

Bases left unread at 30×
under the Poisson model.
Reality is far worse.

0.20%

Autosomal genome below
10% of mean coverage —
at 120× Illumina

395

Medically relevant genes
GIAB originally excluded
from benchmark sets

Every sequencing quote you have ever received specifies a number followed by a multiplication sign. Thirty times. Fifty times. A hundred times for tumor. The number is treated as a guarantee — as though purchasing 30× buys you a genome observed thirty times over, uniformly, everywhere. It does not. The number is a ratio of total bases sequenced to genome size, and it is silent on the only question that matters clinically: was the specific base underlying this patient's phenotype covered deeply enough to call?

This distinction is not academic. It is the difference between a negative result that means "no variant here" and a negative result that means "no data here." Those two are reported identically by most pipelines and are not the same finding. A laboratory that cannot separate them is issuing clinical reports with an unstated confidence interval.

The Arithmetic

Coverage is a distribution, not a value

The idealized model is old and still useful. Lander and Waterman treated read start positions as a Poisson process along the genome¹. Under that model, with mean coverage c, the probability that any particular base receives exactly k reads is P(k) = e−c · ck / k!. The probability it receives zero reads is e−c. At 30×, that is roughly 9 × 10⁻¹⁴ — about 0.0003 bases in a 3.1 Gb genome. Comfortable.

Except the model assumes reads land independently and uniformly. Neither assumption survives contact with a real library. Fragmentation is not random with respect to sequence. PCR during library construction and cluster generation amplifies GC-extreme templates inefficiently. Alignment is not possible where the reference is ambiguous. The result is not a Poisson distribution around 30 — it is a heavy-tailed distribution with a long left shoulder, and the left shoulder is where your clinically actionable genes live.

The Poisson model tells you what coverage would look like if the genome were a featureless string. The genome's features are precisely what you are sequencing it to find.

Ross and colleagues quantified this directly. In Illumina data at 120× mean coverage — four times the clinical standard — 0.20% of the autosomal genome sat below 10% of the genome-wide average². After excluding regions attributable to known bias motifs or sample-reference differences, 0.045% remained unexplained. Their analysis also established that GC-rich regions prone to low coverage disproportionately include human promoters, and they catalogued a thousand that proved exceptionally resistant to sequencing².

Read that again with a clinical lens. Promoters. Regulatory regions where variants alter expression of the gene downstream. At 120×. The distribution does not care about your mean.

Figure 1. The same mean, two different genomes. The idealized Poisson distribution (dashed) places almost no mass below the callable threshold. The observed distribution (solid) is left-skewed by GC bias, mappability, and library artifacts. The shaded region is the part of the genome you paid to sequence and did not receive. Schematic; not to scale.

The Substitution

Replace the mean with a callable fraction

The remedy is unglamorous and entirely within a laboratory's control: stop reporting mean coverage as the quality metric and start reporting the callable fraction — the proportion of your target region that met a defined per-base depth and mapping-quality floor, computed per sample, per gene, per exon.

This is a policy decision, not a technical one. The tooling has existed for years: mosdepth will produce per-base and per-region depth summaries quickly enough to run on every sample³; samtools depth and GATK's callable-loci logic cover the same ground. What is usually missing is the institutional commitment to define the threshold in advance, apply it uniformly, and let it fail samples.

Defining the floor

A defensible threshold has three parts, and each must be justified from your own validation data rather than inherited from a vendor slide:

• Minimum depth. The read count below which your caller's sensitivity for the variant classes you report falls beneath your validated limit. For germline heterozygous SNVs on short reads this commonly lands around 20×; for indels it is higher; for somatic variants at low allele fraction it is dramatically higher and depth alone is insufficient.
• Minimum mapping quality. Reads with MAPQ 0 are reads the aligner could not place uniquely. Counting them toward depth inflates coverage precisely in the segmental duplications and homologous gene families where you most need honest numbers.
• Minimum base quality. Depth composed of low-confidence base calls is depth in name only.

A base that clears all three is callable. A base that fails any is a gap, and gaps get reported. Not buried in a QC appendix — reported alongside the variant list, gene by gene.

Bioinformatics note

A common failure in exome and panel pipelines is computing callable fraction against the capture BED rather than the coding regions of the transcript actually reported. Capture kits routinely include probe-adjacent padding that inflates the denominator. The correct denominator is the clinically reported region — the CDS plus your defined splice-site window on the reference transcript — intersected with the target, not the target alone.

Where It Breaks

The genome is not uniformly sequenceable

Low coverage is not randomly scattered. It concentrates in exactly the regions that are hardest to sequence and, by an unkind coincidence of biology, often the most clinically consequential. The Genome in a Bottle Consortium made this explicit when constructing its benchmarks: nearly 400 medically relevant genes were excluded from earlier GIAB variant benchmark sets because of repetitiveness or polymorphic complexity⁴.

Wagner and colleagues subsequently characterized 273 of those 395 challenging autosomal genes using a haplotype-resolved whole-genome assembly, producing curated benchmarks against both GRCh37 and GRCh38 for HG002⁴. The most instructive finding was not the benchmark itself but a diagnosis of the reference: false duplications in GRCh37 and GRCh38 caused reference-specific missed variants for both short- and long-read technologies in genes including CBS, CRYAA, and KCNE1. When those false duplications were masked, variant recall improved from 8% to 100% in affected regions⁴.

Eight percent to one hundred. That is not a coverage problem being fixed by sequencing deeper. That is an alignment problem masquerading as a coverage problem, and no amount of additional reads would have solved it. The reads were present. They were being placed in a duplicated copy that does not exist in the sample, dropping their mapping quality to zero and, under any sane callable definition, removing the region from the callable set.

If your only lever is depth, every gap looks like a sequencing failure. Most are not.

Caution

Increasing mean coverage from 30× to 60× will not rescue a region that is uncallable because of mappability. Both the numerator and the denominator scale. Doubling the depth of MAPQ 0 reads produces twice as many reads you have agreed not to trust. Before authorizing deeper sequencing to close a gap, determine whether the gap is a depth gap or a mappability gap. The remedies are entirely different — the first is more reads, the second is a better reference, a different aligner, a long-read orthogonal assay, or a documented limitation.

Making It Operational

What the report should say

The practical output of this argument is a change to two documents: the validation summary and the clinical report.

Metric	What it claims	What it actually establishes
Mean coverage	The genome was sequenced 30 times over.	Total bases sequenced ÷ target size. Silent on distribution. Inflated by high-coverage outlier regions.
% target ≥ 20×	Nearly all bases are callable.	A depth floor only. Includes MAPQ 0 pileups unless explicitly filtered.
Callable fraction	These bases met depth, MAPQ, and base-quality floors.	A defensible denominator for sensitivity. Falsifiable, per-sample, per-gene.
Per-gene gap list	These specific coordinates were not interrogated.	Converts an unstated limitation into a disclosed one. The only honest form of a negative result.

The gap list is the deliverable that changes clinical behavior. When an ordering clinician sees that exon 7 of the gene matching the phenotype was not covered, they order a confirmatory assay. When they see "mean coverage 34.2×, 98.1% of target ≥ 20×," they conclude the gene was ruled out. One of those two reports caused a diagnostic odyssey to end. The other extended it.

Three changes worth making this quarter

None of these require new instrumentation, and all of them are auditable:

• Compute callable fraction with mapping-quality filtering enabled, and against the reported transcript CDS rather than the capture target. Compare the two numbers on twenty historical samples. The delta is the size of the problem you have been carrying.
• Emit a per-gene gap list for every clinical report, including when it is empty. An empty gap list is a positive assertion of complete interrogation and is far more valuable than a coverage percentage.
• Classify each recurrent gap as depth-limited or mappability-limited before proposing a remedy. Recurrent gaps that appear in every sample at the same coordinates are almost never depth problems.

Thirty times is a purchasing specification. It belongs on an invoice. It does not belong in the sentence that tells a family whether the gene was examined.

Work with Zetobit

Coverage policy is validation policy.

Zetobit designs and validates CAP/CLIA-compliant NGS pipelines for clinical laboratories, biotech, and pharma — including callable-region definition, per-gene gap reporting, and the validation evidence that supports both.

Start a conversation

References

1. Lander, E. S. & Waterman, M. S. Genomic mapping by fingerprinting random clones: a mathematical analysis. Genomics 2, 231–239 (1988). doi:10.1016/0888-7543(88)90007-9
2. Ross, M. G. et al. Characterizing and measuring bias in sequence data. Genome Biology 14, R51 (2013). doi:10.1186/gb-2013-14-5-r51
3. Pedersen, B. S. & Quinlan, A. R. Mosdepth: quick coverage calculation for genomes and exomes. Bioinformatics 34, 867–868 (2018). doi:10.1093/bioinformatics/btx699
4. Wagner, J. et al. Curated variation benchmarks for challenging medically relevant autosomal genes. Nature Biotechnology 40, 672–680 (2022). doi:10.1038/s41587-021-01158-1

Zetobit, LLC
710 E Main St, Lexington, KY 40502
Bioinformatics Insight Series — original technical writing on genomics, clinical pipelines, and multi-omics.
zetobit.com