Coverage and callability: what your file let us read, and what it did not

Depth is a count, callability is a judgement

Stack a genome's reads on the reference and count the reads over one base: that is its depth. It is cheap to measure and easy to misread, because a deep stack is not automatically a useful one. A read carries its own estimate of how likely it is to sit in the wrong place: the mapping quality, defined in the SAM specification as minus ten times the logarithm of the chance that the position is wrong^[3]. A mapping quality of 20 means one chance in a hundred. A read that fits equally well in more than one place gets a very low score, and most variant callers ignore such reads^[3].

Callability folds those things together. The Medical Genome Initiative, a consortium of North American healthcare and research organisations, defines a callable region as one where a single-nucleotide genotype can be reliably derived, and says an assessment of callable regions should use depth, base quality and mapping quality^[1]. Depth is one input. A base can have plenty of it and still be uncallable.

Where a 30x genome thins out

If reads landed at random, a 30x genome would almost never fall below ten reads at a base: under that assumption the chance is about seven in a million (our arithmetic, treating depth as a Poisson count with a mean of 30). Real genomes are not. Sequencing under-represents DNA fragments that are very rich in G and C, and also those very rich in A and T, most likely because of the PCR step, and the size of the effect is not the same from one sample to the next^[4]. A stretch that is hard to amplify comes back thin in one person's genome and fine in the next.

The second failure is structural and the same in everyone. Ebbert and colleagues went looking for it in whole-genome data from ten men, defining a region as dark when it had five reads or fewer, or when at least nine in ten of its reads had a mapping quality under 10^[5]. They found dark regions in thousands of genes.

Many of those regions are what the authors call camouflaged: a gene with a near-identical copy elsewhere, so that a short read fits both places equally well. SMN1 and SMN2, the genes behind spinal muscular atrophy, were completely or nearly completely camouflaged in their data^[5]. Sequencing deeper does not fix that, because the extra reads are just as ambiguous. Longer reads help: in the same study, linked-read and long-read technologies cut the dark protein-coding stretches to between about half and a tenth of the short-read figure^[5].

Where to draw the line is a choice

Different groups draw the line in different places, each for a reason. Ebbert and colleagues note that 20 to 30 reads is often considered a reasonable minimum to identify a heterozygous change with confidence, and set their own dark-region line far lower, at five, to stay conservative^[5]. gnomAD, the public catalogue of human variation, counts a genotype as high quality from a depth of 10 with a genotype quality of at least 20^[6].

Our floor of 10 sits at the permissive end of that range, level with gnomAD's depth rule^[7]. It is a floor for saying a base was read, which is a weaker statement than saying a heterozygous change there would have been caught every time. That is why the report never turns a callable gene into a statement that you do not carry something.

A plain VCF has no record of where it looked

A plain VCF lists the positions where a genome differs from the reference and says nothing about the rest. When it has no line for a position, the reads may have matched the reference there, or there may have been too few reads to call anything. A gVCF closes that gap with reference blocks, records that cover a run of positions where the reads matched and give the depth behind them^[8]. Our guide to the genome file formats covers how the two are written.

#CHROM  POS   REF  ALT     INFO      FORMAT              Sample
1       4370  G    <*>     END=4383  GT:DP:GQ:MIN_DP:PL  0/0:25:60:23:0,60,900
1       4384  C    <*>     END=4388  GT:DP:GQ:MIN_DP:PL  0/0:25:45:25:0,42,630
1       4389  T    TC,<*>  .         GT:DP:GQ:PL         0/1:23:99:51,0,36,93,92,86
1       4390  C    <*>     END=4390  GT:DP:GQ:MIN_DP:PL  0/0:26:0:26:0,0,315
1       4391  C    <*>     END=4395  GT:DP:GQ:MIN_DP:PL  0/0:27:63:27:0,63,945

A gVCF can also say no outright. The specification writes a position where no call could be made as ./.^[8], and a block can carry that genotype with a depth of thirty. We have seen a real whole-genome gVCF do this at every position our report reads. Counting depth alone would have marked almost all of them covered, so the code skips any record whose genotype is a no-call, whatever its depth.

Nor does the format guarantee a complete file. Of the two real 30x gVCFs we measured, one confirmed 99.0 percent of the 25.5 million bases our report examines and the other only 14.2 percent, because the pipeline that wrote it left most of the panel without records^[9]. The two plain VCFs confirmed 0.17 and 0.16 percent, only the bases where they listed a variant^[9].

Chips: on the array or not examined

A chip export from 23andMe, AncestryDNA or a similar service has no depth at all. The chip measures a fixed list of positions chosen before your DNA arrived and reports your two letters at each. A position off the list was never looked at, and nothing in the file can turn that into a result. We counted how many of the 1,572 positions our report scores ten real exports could read, counting a position only when our reader would use it: no-calls, insertions, deletions and the ambiguous A/T and C/G pairs do not count^[9].

The row count is a poor guide. The 23andMe v3 export had the most rows of the ten, 955,958, and read fewer scored positions than the 2024 AncestryDNA file with 674,131; the 2018 MyHeritage file had more rows than that AncestryDNA file and read 308^[9]. What matters is which positions a chip was built to measure, and that varies by product and year.

Being on the array is also not the same as being read correctly. Chips call common variants well and very rare ones badly: in UK Biobank, only 16 percent of 4,757 heterozygous chip calls for variants with a frequency under 0.001 percent were confirmed by sequencing^[10]. That is a separate problem from coverage, and our piece on chip false positives deals with it.

What the report says, file by file

For a position the report needs and your file does not list as a variant, what can be said depends on the file, and the report uses different words for each case.

On a plain VCF the report opens with a banner headed 'Read this as a variant-only file', which says an absence means no variant was reported and could be the reference letter or a position your file did not cover. The trait and similar modules report how many of their positions were read from your file, and calls the rest read as the common (reference) genotype, 'an inference, not a measurement'. As on every sequenced file, the clinical-findings card says "Not found" means no matching variant was present in your file, not "examined and absent". There is no per-gene table, because there is nothing to fill it with.

On a gVCF the banner gives the share of panel bases callable at read depth 10 or more, and a table gives each gene's callable percentage with a band. A second table goes through the actionable conditions one by one: screened, nothing found; screened except a named gene your file could not read; or could not be read at all. Trait and other modules say how many of their positions were examined, 'proven callable by your file's own coverage'. Here the 90 percent rule applies: a reference position counts as examined only when every panel region of its gene is at least 90 percent callable. Positions outside the measured regions stay inferences, and the report says how many.

On a chip export there is no coverage table. Positions are on the array or they are not. A single-variant finding the chip could not read is named and labelled 'not examined', and a trait whose marker is missing from the array is left unresolved rather than shown as the common version. Medication cards list each defining position as on your array, on it but not resolvable from chip data, or not on your array.

Deep Read measures it from the reads

With a BAM or CRAM, Deep Read runs mosdepth, a coverage tool^[11], over each gene on its pharmacogene panel, counts the bases covered by at least ten reads with mapping quality 20 or more, and divides by the length of the gene's region. Each gene's card carries the result as a badge, such as '94% callable', flagged when it is under 90 percent. The badge's own note says what it is: the share of the gene's region read at 10x or more, 'measured across the whole gene, not at the positions that define your alleles'.

Below 30 percent the report withholds the phenotype. The card says the coverage over that gene was low, gives the percentage, and says we did not assess a phenotype from your reads rather than risk a false-confident call. That threshold is deliberately low: it catches a gene that essentially could not be read, and leaves everything above it to the diplotype call.

Two limits apply. The per-gene map across the clinical panel and the deletion screen run only on reads aligned to GRCh38; on GRCh37 the pharmacogene panel still runs, and the report says the other two were not measured rather than estimating them. And a camouflaged gene stays hard however good the file: callability says how much of it could be read, and does not make the rest readable. The sample report shows the gVCF tables and the Deep Read badges on an illustrative file.