Why Read Counts Inform Genotype Calls
A sequencing experiment does not necessarily report a perfectly reliable genotype. Instead, it produces reads that may support either the reference allele or an alternate allele. The number of alternate reads can therefore provide evidence for a genotype without determining it with certainty. A genotype likelihood quantifies how compatible the observed read data are with a specified genotype under a stated error model.
Consider a diploid site with two possible alleles: reference and alternate. The three genotypes contain zero, one, or two copies of the alternate allele. If reads are sampled from the two chromosome copies and allele calls can be wrong, each genotype implies a different probability of observing an alternate read. The read-count likelihood follows from that probability and a binomial sampling model.
From the Error Model to the Likelihood
The model assumes that a true reference allele is miscalled as alternate with probability \(e\), and a true alternate allele is miscalled as reference with the same probability \(e\). A genotype with \(g\) alternate copies has alternate-allele fraction \(g/2\). For a randomly sampled chromosome copy, the probability of an alternate call is therefore the probability of sampling a reference copy and miscalling it, plus the probability of sampling an alternate copy and calling it correctly:
Thus \(q_0=e\), \(q_1=1/2\), and \(q_2=1-e\). The heterozygote has alternate-read probability \(1/2\) under this model because the two allele types are equally represented and have symmetric error rates. The homozygote probabilities are near zero or one when the error rate is small, but they are not exactly zero or one unless \(e=0\).
To obtain the likelihood, first imagine the \(R\) reads in a fixed order. Under conditional independence, any particular ordered pattern with \(K\) alternate calls and \(R-K\) reference calls has probability \(q_g^K(1-q_g)^{R-K}\). The observed data are counts rather than a particular ordering. There are \(\binom{R}{K}\) different orderings with those counts, so their probabilities sum to the binomial likelihood in the definition. The binomial coefficient is the same for all candidate genotypes for the same observed data; it can be omitted when comparing their relative likelihoods, but it belongs in the full probability of the count.
Worked Example: Many Alternate Reads
Worked Example: Many Alternate Reads
Suppose a hypothetical sample has \(R=20\) reads, of which \(K=18\) are called alternate. Assume a known symmetric error rate \(e=0.01\), with independent reads and no allele-specific bias. Compare the likelihoods for \(g=0,1,2\).
The genotype-specific alternate-read probabilities are \(q_0=0.01\), \(q_1=0.50\), and \(q_2=0.99\). The count likelihood is \(\binom{20}{18}q_g^{18}(1-q_g)^2\), and \(\binom{20}{18}=190\).
For \(g=0\), $$ L(0;18,20)=190(0.01)^{18}(0.99)^2 \approx 1.862\times10^{-34}. $$ For \(g=1\), $$ L(1;18,20)=190(0.50)^{18}(0.50)^2 =190(0.50)^{20} \approx 0.0001812. $$ For \(g=2\), $$ L(2;18,20)=190(0.99)^{18}(0.01)^2 \approx 0.01586. $$
Among these three likelihoods, \(g=2\) gives the largest value: the observed count of 18 alternate reads is much more compatible with the alternate homozygote than with either of the other genotypes under the stated model. This comparison alone is not a probability that \(g=2\) is the true genotype; it evaluates the read data under each genotype separately.
Comparing Genotype Likelihoods
For two candidate genotypes \(g_a\) and \(g_b\), the likelihood ratio summarizes how much more compatible the observed count is with one candidate than the other. When \(0<e<0.5\), all three values \(q_0,q_1,q_2\) lie strictly between zero and one. The likelihoods are then positive for any possible count \(0\le K\le R\), so their ratio is well-defined:
The binomial coefficient cancels because it is common to both likelihoods. A ratio greater than one means the observed count has a larger likelihood under \(g_a\) than under \(g_b\); a ratio less than one favors \(g_b\) in this likelihood comparison. The ratio is not a posterior odds ratio: it contains no genotype-frequency information or other prior information.
Worked Example: Balanced Reads and a Heterozygote
Consider hypothetical data with \(R=10\) and \(K=5\), and take \(e=0.02\). Assume the same independent-read, symmetric-error model. Compare the heterozygote likelihood with the likelihood for either homozygote.
Here \(q_1=0.50\), \(q_0=0.02\), and \(q_2=0.98\). The heterozygote likelihood is $$ L(1;5,10)=\binom{10}{5}(0.50)^5(0.50)^5 =252(0.50)^{10} \approx 0.2461. $$ The reference homozygote likelihood is $$ L(0;5,10)=252(0.02)^5(0.98)^5 =252(0.0196)^5 \approx 7.29\times10^{-7}. $$
For the alternate homozygote, interchanging the alternate and reference counts leaves the probability unchanged: $$ L(2;5,10)=252(0.98)^5(0.02)^5 \approx 7.29\times10^{-7}. $$
The ratio of either homozygote likelihood to the heterozygote likelihood is $$ \frac{L(0;5,10)}{L(1;5,10)} =\left(\frac{0.02}{0.50}\right)^5 \left(\frac{0.98}{0.50}\right)^5 =(0.0784)^5 \approx 2.962\times10^{-6}. $$ The observed balanced count is therefore far more compatible with \(g=1\) than with either homozygote under these assumptions. This is a comparison of read-data likelihoods, not a claim that the heterozygote has a particular posterior probability.
Boundary Cases and Log-Likelihoods
The error-rate boundary \(e=0\) requires special care. In that case the model is deterministic for the two homozygotes: \(q_0=0\) and \(q_2=1\), while \(q_1=1/2\). Some observed counts then have zero likelihood under a candidate genotype. Expressions involving a ratio of likelihoods or logarithms must not be used as if every likelihood were positive.
Worked Example: Zero Error Rate
Suppose a hypothetical sample has \(R=4\), \(K=0\), and \(e=0\). The model gives \(q_0=0\), \(q_1=0.5\), and \(q_2=1\). For \(g=0\), all four reference reads are certain under the model: $$ L(0;0,4)=\binom{4}{0}(0)^0(1)^4=1, $$ where the zero-count factor is understood as contributing 1. For \(g=1\), $$ L(1;0,4)=\binom{4}{0}(0.5)^0(0.5)^4=0.0625. $$ For \(g=2\), no reference reads are possible: $$ L(2;0,4)=\binom{4}{0}(1)^0(0)^4=0. $$
The data have positive likelihood under \(g=0\) and \(g=1\), but zero likelihood under \(g=2\). A ratio with \(L(2;0,4)\) as its denominator is undefined; it is not a finite likelihood ratio. If two candidate genotypes both assign zero likelihood to an observation, their ratio is also undefined, and the data do not distinguish them by that ratio.
The likelihood itself is always the safest starting point at the boundary. For log-likelihood calculations, use \(\log L\), with \(\log(0)=-\infty\). Equivalently, in the expanded expression
a term whose count is zero contributes zero, including the limiting convention \(0\log(0)=0\). A positive count multiplying \(\log(0)\) contributes \(-\infty\), correctly indicating a zero likelihood. This convention avoids treating \(0\log(0)\) as an ordinary numerical product. At \(e=0.5\), all genotypes have \(q_g=0.5\), so the read counts have identical likelihoods under all three genotypes and provide no genotype discrimination in this model.
Practical Interpretation and a Common Pitfall
A genotype likelihood is a function of the candidate genotype with the read data held fixed. It answers: under the specified error model, how probable are these read counts if this genotype were the genotype? It is not the probability of the genotype given the reads. Converting likelihoods into posterior genotype probabilities requires additional information, such as genotype prior probabilities, and is a separate calculation.
A common mistake is to say, “The genotype with the largest likelihood has probability one,” or to treat the likelihood values as probabilities that the candidate genotypes are true. The correction is that the largest likelihood identifies the maximum-likelihood genotype among the candidates under this read model, but it does not by itself express uncertainty about the genotype. The likelihood comparison also depends on the model assumptions and the specified error rate.
Check Your Understanding
Use the symmetric-error binomial model described above.
- For \(e=0.03\), calculate \(q_0\), \(q_1\), and \(q_2\). Why does the heterozygote probability equal one-half?
- For \(R=8\) and \(K=6\), write the full likelihood for \(g=2\) when \(e=0.02\), including the binomial coefficient.
- Explain why the binomial coefficient cancels when taking a likelihood ratio for two candidate genotypes with the same observed \(K\) and \(R\).
- At \(e=0\), what is \(L(2;0,4)\), and why is a likelihood ratio with that value as denominator not defined?
- In one sentence, distinguish a genotype likelihood from the probability of a genotype given the read data.