Figure 1. Orientations and phases of gene overlap. Genes can overlap on the same strand or on the opposite strand. The reference gene in a pair of overlapping genes is called phase 0. Same-strand overlaps can be in two phases (1 and 2); opposite-strand overlaps can be in three phases (0, 1, and 2). First and second codon positions, in which ,5% and 0% of the changes are synonymous, are marked in red. Third codon positions, in which ,70% of the changes are synonymous, are marked in blue. doi:10.1371/journal.pone.0003996.g001 overlap phases (Figure 1). To understand the consequences of estimating selection pressures on overlapping genes as if they are independent genes, let us consider a simplified view of the genetic code, in which all changes in first and second codon positions are nonsynonymous and all changes in third codon position are synonymous. (In reality, the proportions of changes that are synonymous are ,5%, 0%, and ,70% for the first, second, and third codon positions, respectively). From Figure 1 we see that in all overlap types, but one (opposite-strand phase 2), all synonymous changes in one gene are nonsynonymous in the overlapping gene, while half of the nonsynonymous changes are synonymous in the overlapping gene. Since the rate of synonymous substitutions is usually higher than that of nonsynon- ymous substitutions, ignoring overlap constraints would result in the underestimation of the rate of synonymous substitutions. (In the case of opposite-strand phase-2 overlaps, ignoring the overlap would result in the underestimation of nonsynonymous substitu- tions rate.) The bias in the estimation would be correlated with the strength of purifying selection on the overlapping gene. Thus, a false inference of positive selection is likely for genes under relaxed purifying selection when the overlapping gene is under strong purifying selection. Goldman and Yang’s [37,38] method for the estimation of selection intensity in non-overlapping coding sequences The most commonly used method for estimating selection intensity on protein coding genes fits a Markov model of codon substitution to data of two homologous sequences [37,38]. The codon-based model of nucleotide substitution is specified by the substitution-rate matrix, Qcodon = {qij}, where qij is the instantaneous rate of change from codon i to codon j. 8 0, if i and j differ at two or three codon positions, > > >>> pj, if i and j differ by a synonymous transversion, > > >< kpj, if i and j differ by a synonymous transition, qij~ 1 ð Þ > vpj, if i and j differ by a nonsynonymous > >>>>> transversion, >: vkpj, if i and j differ by a nonsynonymous transition: PLoS ONE | www.plosone.org 2 December 2008 | Volume 3 | Issue 12 | e3996 Selection in Overlapping Genes Here, k is the transition/transversion rate, v is the nonsynon- ymous/synonymous rate ratio (dN/dS), and pj is the equilibrium frequency of codon j, which can be estimated from the sequence data by several models [Fequal, F164, F364, and F61, reviewed in 38]. Parameters pj and k characterize the pattern of mutations, whereas v characterizes selection on nonsynonymous mutations. Qcodon is used to calculate the transition-probability matrix 8 0, > > > > > > > >>> pj, > > > > > > >>> kpj, > > > > > > > > >>> v1pj, > > > > > > > > > > > > >>> v2pj, > > > < qij~ if i and j differ at two or three codon positions, if i and j differ by a synonymous transversion in both genes, if i and j differ by a synonymous transition in both genes, if i and j differ by a nonsynonymous transversion in gene A and synonymous in gene B, if i and j differ by a nonsynonymous transversion in gene B and synonymous in gene A, if i and j differ by a nonsynonymous > > > >>> v1kpj, > >>>> transition in gene A and synonymous > >>>> in gene B, > >>>> v2kpj, if i and j differ by a nonsynonymous > >>>> transition in gene B and synonymous > >>>>> in gene A, >>>> v1v2pj, if i and j differ by a nonsynonymous > >>>>> transversion in both genes, >>>> v1v2kpj, if i and j differ by a nonsynonymous > > : transition in both genes: PðtÞ~ pijðtÞ ~eQcodont , where pij(t) is a probability that a given codon i will become j after time t. Parameters k, t, and v are estimated by maximization of the log-likelihood function ‘ðtÞ~ X X nij log pipijðtÞ , i j where nij is the number of sites in the alignment consist of codons i and j. The estimated parameters are then used to calculate dN and dS [38]. A new method for the simultaneous estimation of selection intensities in overlapping genes We follow the maximum likelihood approach of Goldman and Yang [37,38] to construct a model that accounts for different selection pressures on the genes in the overlap. We start with the simplest case, that of opposite-strand phase-0 overlaps. The reason this is the simplest case is that each codon overlaps only one codon in the overlapping gene. The substitution of nucleotides in opposite-strand phase-0 overlaps is specified by the substitution- rate matrix, Qcodon = {qij}, where qij is the instantaneous rate of change from codon i to codon j. 2 ð Þ 3 ð Þ 5 ð Þ