The main difference between this model and the single-gene model is that here we distinguish between two dN/dS ratios (v1 and v2 for gene 1 and gene 2, respectively). Another difference is the estimation of codon-equilibrium frequencies. Since the parameters of codon frequencies characterize processes that are independent of the selection on overlapping regions, we estimate these frequencies using the non-overlapping regions of each gene. The calculation of the transition-probability matrix and the log-likelihood function is done in the same way as in the single-gene model (equations 2 and 3). The above model is a simple expansion of the single-gene model to account for opposite-strand overlaps in phase 0. However, this model cannot be used in the other four overlap cases, same-strand phase-1 and phase-2 overlaps and opposite-strand phase-1 and phase-2 overlaps, because in all these cases a codon overlaps two codons of the second gene. Therefore, we set the unit of evolution to be a codon (the reference codon) and its two overlapping codons, which together constitute a sextet (Figure 2). The sextet is, therefore, the smallest unit of evolution in overlapping genes. In our model, each gene constitutes a set of sextets and within each sextet, only the reference codon is allowed to evolve. Changes in this codon affect the two overlapping codons. For example, consider the red and blue overlapping genes in Figure 2a. A change from G to A in position five (Figure 2a, bold) is illustrated in Figure 2b for the red gene as a reference and in Figure 2c for the blue gene as a reference. Restricting changes to the reference codon only is essential for the model, since changes outside the reference codon will require the consideration of other overlap- ping codons outside of the sextet, and so ad infinitum. In addition, this restriction allows the model to maintain the assumption that each reference codon evolves independently. For gene A as the A reference gene, we specify the substitution-rate matrix, Q sextet = {qAuv} where q A uv is the instantaneous rate from sextet u to sextet v with the codons of gene A as the reference codons: 8 0, if u and v dif f er at two or three codon > > > > positions or at a position outside the > > > > > ref erence codon, > > > >>>> pv, if u and v dif f er by a synonymous > > > transversion in both genes, > > > >>> kpv, if u and v dif f er by a synonymous > > > > transition in both genes, > > > >>> v1pv, if u and v dif f er by a nonsynonymous > > > > transversion in gene A and synonymous > > > > > in gene B, > > > >>>> v2pv, if u and v dif f er by a nonsynonymous < quAv~ transversion in gene B and synonymousð6Þ >>>> in gene A, > >>>> v1kpv, if u and v dif f er by a nonsynonymous > >>>> transition in gene A and synonymous > >>>> in gene B, > >>>> v2kpv, if u and v dif f er by a nonsynonymous > >>>> transition in gene B and synonymous > >>>> in gene A, > >>>> v1v2pv, if u and v dif f er by a nonsynonymous > >>>> transversion in both genes, > >>>> v1v2kpv, if u and v dif f er by a nonsynonymous > :> transition in both genes: Selection in Overlapping Genes Figure 2. a. An overlapping gene pair (red and blue). b. The codon that is allowed to evolve is marked in red. The substitution in the second- codon position affects the overlapping codon in blue. c. The opposite situation in which only the codon marked in blue is allowed to change. doi:10.1371/journal.pone.0003996.g002 Similarly, we specify the substitution-rate matrix, Q B sextet = {qBuv} B for gene B as the reference gene, where q uv is the instantaneous rate from sextet u to sextet v with gene B codons as the reference A B codons. These substitution-rate matrixes, Q sextet and Q sextet, can be used to calculate transition-probability matrixes (equation 2). However, these transition-probability matrixes cannot be used directly in the maximization of a log-likelihood function (equation 3) because they do not allow changes between any two sextets (as required in a Markov process). For example, the transition probability between sextets AAAAAA and CAAAAA (where the reference codons at positions 3-5 are underlined) would be zero for any given time t, because changes at a position outside of the reference codon are not allowed. A similar difficulty led Pedersen and Jensen [36] to use a complicated, computationally-expensive, simulation procedure A B to estimate model parameters. Hence, we use Q sextet and Q sextet to construct codon-based substitution-rate matrixes QcAodon~ qAij and QcBodon~ qBij by summing the rates over all sextets that share the same reference codon. Similar approach was used by Yang et al. [39] to construct an amino acid substitution-rate matrix from a codon substitution-rate matrix. Let I and J represent the sets of sextets whose reference codons are i and j, respectively, than, the substitution rate from codon i to codon j is qij~ X u[I,v[J quv: 7 ð Þ A B Q codon and Q codon are used to calculate a transition-probability matrix for each of the genes as in equation 2. PAðtÞ~ p A ijðtÞ ~eQcAodon t and PBðtÞ~ p B ij ðtÞ ~eQcBodon t: 8 ð Þ The new transition-probability matrixes are suitable for a maximization of a log-likelihood function since they allow transition between each two codons. PA(t) and PB(t) can be used separately to estimate model parameters in a log-likelihood function for each gene (equation 3). However, in order to use all the information in the data, we combine the two transition- probability matrixes to create the following log-likelihood function: PLoS ONE | www.plosone.org 3 December 2008 | Volume 3 | Issue 12 | e3996