*Department of Biology, University College London, England;
and
Department of Organismic and Evolutionary Biology, Harvard University
![]() |
Abstract |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() |
Introduction |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
The simplest problem in this regard is estimation of the numbers of synonymous (dS) and nonsynonymous (dN) substitutions per site between two sequences. In the past two decades, a number of intuitive methods have been suggested for this estimation. They involve ad hoc treatments that cannot be justified rigorously, and they will be referred to here as approximate methods. In common, they involve three steps. First, the numbers of synonymous (S) and nonsynonymous (N) sites in the sequences are counted. Second, the numbers of synonymous and nonsynonymous differences between the two sequences are counted. Third, a correction for multiple substitutions at the same site is applied to calculate the numbers of synonymous (dS) and nonsynonymous (dN) substitutions per site between the two sequences. Here, we use the notation of Nei and Gojobori (1986
, referred to later as "NG"), which appears to be the most commonly used approximate method; definitions of symbols are given in Table 1
. The reader is referred to Ina (1995, 1996)
for a recent discussion of important concepts. While the above strategy appears simple, well-known features of DNA sequence evolution, such as unequal transition and transversion rates and unequal nucleotide or codon frequencies, make it a real challenge to count sites and differences correctly.
|
A maximum-likelihood (ML) method for estimating dS and dN between two sequences was developed by Goldman and Yang (1994)
based on an explicit model of codon substitution. The ML method does not involve ad hoc approximations. Furthermore, the ML method is flexible in that knowledge of the substitution process such as transition/transversion bias, codon usage biases, and even chemical differences between amino acids can easily be incorporated into the model.
Relying on insights gained through previous methods, particularly ML estimation, we propose in this paper an approximate method for estimating dS and dN that accounts for two major features of DNA sequence evolution: the transition/transversion bias and the base (codon) frequency bias. We examine the similarities and differences among the ML method, the approximate method of this paper, and two other approximate methods by a consistency analysis of infinitely long sequences and computer simulation of finite data. A real data set is also analyzed using different estimation methods.
![]() |
Methods for Estimating Synonymous and Nonsynonymous RatesMaximum-Likelihood Estimation |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Parameter is the (mutational) transition/transversion rate ratio, and
= dN/dS is the nonsynonymous/synonymous rate ratio. The equilibrium codon frequencies (
j) are calculated using the nucleotide frequencies at the three codon positions; that is, codon frequencies are proportional to the products of nucleotide frequencies at the three codon positions. This approach was found to produce results similar to those obtained using all codon frequencies as free parameters, although codon frequencies are often quite different from those expected from nucleotide frequencies at codon positions (e.g., Goldman and Yang 1994
; Pedersen, Wiuf, and Christiansen 1998
; Yang and Nielsen 1998
). We note that different assumptions about codon frequencies can be easily incorporated into the ML method and the approximate method of this paper. Equation (1) is similar to the simulation model of Gojobori (1983)
and the likelihood model of Muse and Gaut (1994)
, although those authors did not consider the transition/transversion rate bias or different base frequencies at the three codon positions.
The diagonal elements of the rate matrix, Q = {qij}, are determined by the mathematical requirement (e.g., Grimmett and Stirzaker 1992
, p. 241) that the row sums are zero:
Because time and rate are confounded, we multiply the rate matrix by a scaling factor so that the expected number of nucleotide substitutions per codon is one:
This scaling means that time t (or, equivalently, branch length or sequence distance) is measured by the expected number of (nucleotide) substitutions per codon.
The transition probability matrix is calculated as
![]() |
![]() |
The log-likelihood function is given by the multinomial probability
where nij is the number of sites occupied by codons i and j in the two sequences. The codon frequencies i are estimated using the observed nucleotide frequencies at the three codon positions in the two sequences. Parameters t,
, and
are estimated by maximizing the likelihood function numerically, and the estimates are used to calculate dS and dN. Specifically, the proportions of synonymous and nonsynonymous substitutions are given as
and *N = 1 -
*S, respectively. The summation is taken over all codon pairs i and j (i
j) that code for the same amino acid, and aai is the amino acid encoded by codon i. The numbers of synonymous and nonsynonymous substitutions per codon are then t
*S and t
*N, respectively. The proportions of synonymous and nonsynonymous sites are defined as the proportions of synonymous and nonsynonymous "mutations" before the operation of natural selection at the amino acid level (Goldman and Yang 1994
; Ina 1995
). These can be calculated in a manner similar to equation (7) as
1S and
1N (equivalent to
S and
N in Goldman and Yang 1994
), using the ML estimate of
but with
= 1 fixed. The numbers of synonymous and nonsynonymous sites per codon are 3
1S and 3
1N, respectively. The numbers of synonymous and nonsynonymous substitutions per site are then dS = t
*S/(3
1S) and dN = t
*N/(3
1N), respectively (see Table 1
).
A New Approximate Method
We suggest an approximate method for estimating dS and dN using the strategy adopted by previous authors: counting sites, counting differences, and correcting for multiple hits. In all three steps, we take into account the transition/transversion rate bias and the base (codon) frequency bias. Approximately, our method is based on the HKY85 nucleotide mutation (substitution) model (Hasegawa, Kishino, and Yano 1985
). Although this is not the most general model available, it accounts for two most important features of the mutation process, that is, the transition/transversion bias and unequal base frequencies. Previous results (e.g., Yang 1994a
) suggest that adding further complication is often unnecessary.
Estimating the Transition/Transversion Rate Ratio ()
We use the fourfold-degenerate sites at the third codon positions and the nondegenerate sites to estimate . Mutations at the fourfold-degenerate sites do not change the amino acid, and thus the transition/transversion rate bias at those sites should reflect the mutational bias. Mutations at nondegenerate sites all lead to amino acid changes and can also be used to estimate
(see eq. 1). Here, we assume that different nonsynonymous substitutions have the same rate irrespective of the pair of amino acids involved, although the assumption is unrealistic (Yang, Nielsen, and Hasegawa 1998
). We calculate an average of
, weighted by the numbers of nucleotide sites in the two site classes. Since no simple formula is available for estimating
under the HKY85 model, we use the formula for the F84 model (Yang 1994 b
) instead, relying on the similarity of the two models. We calculate
where T and V are proportions of transitional and transversional differences, respectively, and Y =
T +
C and
R =
A +
G. Then,
![]() |
(Yang 1994b
). The estimated
F84 is then transformed to
(i.e.,
HKY85) using the following formula (see Goldman 1993
)
For data of multiple sequences, we suggest estimating a common by averaging estimates from all pairwise comparisons and using the combined estimate of
in the calculation of pairwise dN and dS rates.
Counting Synonymous and Nonsynonymous Sites
Inas (1995)
Table 1
for counting synonymous and nonsynonymous sites in each codon is correct for mutation models more general than that of Kimura (1980)
, although Inas Table involves minor errors for codons that can change to stop codons in one step. In general, synonymous and nonsynonymous sites can be counted as in the ML method discussed above for any codon-substitution model (Goldman and Yang 1994
). We count sites using codons in the two compared sequences, rather than the equilibrium codon frequencies expected from the model (see discussion in Yang and Nielsen 1998
). As there should be about 4% loss of sites due to mutations to stop codons, this scaling means that we are slightly underestimating dS and dN, although the
ratio is not affected (Yang and Nielsen 1998
). The numbers of sites (S and N) are scaled so that S + N = 3Lc, where Lc is the number of codons. Nucleotide frequencies at synonymous and nonsynonymous sites are recorded and used later for multiple-hit corrections.
Counting Synonymous and Nonsynonymous Differences
Observed nucleotide differences between the two sequences are classified into four categories: synonymous transitions, synonymous transversions, nonsynonymous transitions, and nonsynonymous transversions. When the two compared codons differ at one position, the classification is obvious. When they differ at two or three positions, there will be two or six parsimonious pathways along which one codon could change into the other, and all of them should be considered. Since different pathways may involve different numbers of synonymous and nonsynonymous changes, they should be weighted differently. Miyata and Yasunaga (1980)
and Li, Wu, and Luo (1985)
made valuable attempts to weight pathways. They also pointed out that equal weighting of pathways, which is used in later methods such as those of Nei and Gojobori (1986)
and Ina (1995)
, biases estimates of
toward 1; that is, equal weighting tends to overestimate
when
< 1 and to underestimate
when
> 1.
The most appropriate weights are the relative probabilities of pathways, which we will use in the new approximate method. The probabilities depend on the parameters being estimated: the sequence divergence level (t), the transition/transversion rate ratio (), and the dN/dS ratio (
). For given values of t,
, and
, it is easy to construct the rate matrix Q and calculate the transition probability matrix P(t) (see eqs. 1 and 4). We use the Taylor expansion in this case for its fast speed.
The number of terms used is determined by a preset accuracy level. The weight, that is, the probability for each pathway, is calculated as the product of the probabilities of all changes involved in the pathway. Pathways involving stop codons are given weight 0. If all pathways involve stop codons (for example, between AAG and TGG in the mammalian mitochondrial code), ad hoc decisions have to be made.
An example is given in Figure 1
using a pair of codons in the mitochondrial genes of the human and the orangutan. The concatenated sequence of 12 genes on the H-strand has 3,331 codons (see below). Between the two sequences, 1,198 codons are different at one position, 151 are different at two positions, and 21 are different at all three positions. The 329th codon is TTA in the human and CTC in the orangutan. Transition probabilities for changes involved in each of the two pathways are given in Figure 1
, calculated using the estimates obtained by the new method (t = 0.873, = 10.61, and
= 0.057). The probability for the first path is then pTTA,TTC(t) x pTTC,CTC(t) = 0.00541 x 0.04974 = 0.00027, while that for the second path is 0.29219 x 0.08679 = 0.02536. Weights for the two pathways are thus 0.011 and 0.989, and there are 0.022 nonsynonymous and 1.978 synonymous differences between the two codons. Since the nonsynonymous rate is much lower than the synonymous rate (
< 1), the first path is much less likely. Equal weighting of pathways would give one synonymous and one nonsynonymous difference between the two codons. Other codon pairs may not be so extreme, as the different pathways may involve the same numbers of different types of changes.
|
The Algorithm
Our method for estimating dS and dN can be summarized in the following iterative algorithm.
1. Estimate from the fourfold-degenerate sites and the nondegenerate sites under the HKY85-F84 model using base (codon) frequencies from the real data. The estimated
is used in later steps.
2. Count the numbers of synonymous and nonsynonymous sites (S and N, respectively) using the estimated and the observed base (codon) frequencies.
3. Choose starting values for t and (e.g., using estimates from the NG method).
4. Count the numbers of synonymous and nonsynonymous differences (both transitions and transversions) using , the codon frequencies, and the current values of t and
. The transition probability matrix P(t) is calculated by equation (12) and used to weight pathways when the two codons differ at more than one position. This step generates the proportions of transitional (T) and transversional (V) differences for each of the synonymous and nonsynonymous site classes.
5. Correct for multiple hits to calculate dS and dN using counts of sites and differences and base frequencies at synonymous and nonsynonymous sites. This step updates t and : t = dS x 3S/(S + N) + dN x 3N/(S + N), and
= dN/dS.
6. Repeat steps 45 until the algorithm converges.
In general, two or three rounds of iteration are sufficient. Some variations of the above algorithm are possible. For example, one may use an estimate of obtained externally. Furthermore, no iteration is needed if pathways are weighted equally when counting differences.
![]() |
Comparison of Methods for Estimating dS and dN |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Two approaches are taken to evaluate the methods. The first examines infinitely long sequences and may be termed a "consistency analysis." Instead of the observed codons in the two sequences, the data consist of the expected frequencies (fij) of all 61 x 61 codon "site patterns," calculated using equation (5) for given parameters t, ,
, and codon frequencies
j. ML estimates are known to be statistically consistent when the model is correct (Stuart, Ord, and Arnold 1999
, chapter 18). Since the dN and dS rates (and their ratio) are defined as functions of parameters t,
,
, and
j, ML estimates of dN and dS will also be consistent. Approximate methods, including the method of this paper, involve ad hoc approximations and in general do not give the true values as estimates. They are statistically inconsistent. However, a good approximate method should not deviate too far from the truth with an infinite amount of data. The second approach we take is computer simulation. Finite data sets are generated by simulation and then analyzed by different methods to examine their biases and sampling variances.
We examine effects of the transition/transversion rate ratio (), base (codon) frequencies, and the selective pressure on the gene reflected in parameter
. We initially fix t = 1 nucleotide substitution per codon, although the effect of sequence divergence is examined later. For a neutral gene (
= 1), this translates to 1/3 synonymous and 1/3 nonsynonymous substitutions per site. Three values of
are considered:
= 1 (no selection),
= 0.3 (purifying selection) and
= 3 (positive selection). Estimates of
(dN/dS) from real data vary widely from gene to gene, and
= 0.3 appears to represent moderate purifying selection (see, e.g., Ohta 1995
; Li 1997
; Yang and Nielsen 1998
; Eyre-Walker and Keightley 1999
). There are not many genes under positive selection, but estimates at about
= 3 are found in real data (e.g., Lee, Ota, and Vacquier 1995
; Messier and Stewart 1997
). Three sets of base frequencies at codon positions are used. The first set has equal base (codon) frequencies. The second set is from primate mitochondrial protein-coding genes and has very biased base frequencies. The third set is from HIV env genes (Table 2
). The universal genetic code is used in both the consistency analysis and the computer simulation.
|
Estimates of by different methods are plotted against the transition/transversion rate ratio
for different values of the dN/dS ratio (
). Results for the three sets of codon frequencies are shown in Figure 2AI.
|
The case of = 10 is explored in Table 3
, which shows that NG substantially underestimates
and gives 0.669, 0.216, and 1.812 for
= 1, 0.3, and 3, respectively. In this case, the proportion of synonymous sites (S%) should be 33%, but NG (assuming
= 1) gives 26% (Yang and Nielsen 1998
, Fig. 3
). Use of equal weighting (assuming
= 1) in counting differences by NG tends to bias the estimate of
toward 1. However, compared with the bias in counting sites, the bias in counting differences is much less important because there may not be many pairs of codons that are different at two or three positions and because different pathways may involve the same numbers of synonymous and nonsynonymous changes. As mentioned above, NG underestimates S considerably (25.5/32.5 = 0.785) when
= 10, and almost all of this underestimation is translated into the overestimation of dS (0.3333/0.4134 = 0.806 for
= 1). For similar reasons, the underestimation of
by the NG method is more serious for
= 3 than for
= 0.3 (Fig. 2B and C
). In the latter case, equal weighting assuming
= 1 in the NG method counterbalances the effect of ignoring
in counting sites, while in the former case, the two biases are in the same direction (Table 3
).
|
|
Primate Mitochondrial Codon Frequencies
Results obtained using base frequencies at the three codon positions from the primate mitochondrial genes (see Table 2
) are shown in Figure 2DF.
Estimates of S, dS, and dN for = 10 are shown in Table 3
. The results are very different from those of Figure 2AC
under equal codon frequencies. Except for small values of
(
< 2), Inas method performs more poorly than NG. The two methods give very different counts of sites (S). While transition bias always leads to more synonymous sites, the effect of base frequency bias is more complicated. Extreme codon usage bias can cause the proportion of synonymous sites to range from 0% (e.g., when only codons TTC and TTA are present in the sequences) to 100% (e.g., when only codons CTT, CTC, CTA, and CTG are present). There tend to be more synonymous sites if the two most frequent nucleotides at third positions are both purines or both pyrimidines.
For the mitochondrial genes, S is much smaller than expected under equal base frequencies, causing NG to overestimate rather than underestimate S. For example, NG gives S% = 26.7%, which is higher than the correct value at either = 1 (23.1%) or
= 10 (23.7%) (Table 3
). The overestimation of S caused by ignoring the base frequency bias more than compensates for the underestimation caused by ignoring the transition bias. As a result, NG overestimates
when
= 1 (with estimates from 1.1 to 1.3; Fig. 2D
) or
= 0.3 (with estimates from 0.42 to 0.46; Fig. 2E
). When
= 3, equal weighting of pathways (assuming
= 1) in counting differences combined with the assumption of no transition bias (
= 1) in counting sites cancels the effect of the base frequency bias in counting sites, such that NG produces a quite reliable estimate of
(Fig. 2F
). Inas (1995)
method, by considering the transition bias alone and ignoring the base/codon frequency bias, substantially overestimates the proportion of synonymous sites and overestimates
(Table 3
and Fig. 2DF
). Nevertheless, it should be noted that the observed pattern depends on the particular set of codon frequencies. For frequencies from other data sets, NG may be considerably worse than Inas method.
The method of this paper is slightly better than NG for = 1 although the two methods have opposite biases. When
= 0.3, the new method has little bias. When
= 3, the new method underestimates
, with estimates from 2.5 to 2.6. Since the new method counts sites correctly (see Table 3
), the bias must be due to counting of differences and correction for multiple hits. Table 3
suggests that the new method underestimates both dS and dN, but the underestimation of dN is more serious, leading to underestimation of the
ratio. Apart from the case in which
= 3, the new method is better than both NG and Inas method.
HIV env Codon Frequencies
Figure 2GI
shows estimates of when base/codon frequencies from the HIV envelope genes (see Table 2
) are used. Base frequencies in this gene are less biased than are those in the mitochondrial genes, and the effect of ignoring the base frequency bias is minor. For example, the correct proportion of synonymous sites at
= 1 is 21.9%, while NG gives 24.0%, with very slight overestimation. Patterns in Figure 2GI
are quite similar to those for equal codon frequencies (Fig. 2AC
). Exactly at
= 1, NG gives the estimates 1.105, 0.371, and 2.554 when the true
= 1, 0.3, and 3, respectively. The estimates are biased toward 1, mainly due to the use of equal weighting in counting differences. When
= 10, NG underestimates the proportion of synonymous sites (24.0% vs. the correct value, 28.5%). The bias is not as extreme as that in the case of equal codon frequencies, as unequal base/codon frequencies appear to counterbalance the effect of transition bias to some extent (Tables 3 and 4
).
Inas (1995)
method overestimates the
ratio because it ignores the base frequency bias and thus overestimates the number of synonymous sites. The bias is not as extreme as it is for mitochondrial genes. The new method gives estimates very close to the true values for
= 1 and
= 0.3. When
= 3, the new method slightly underestimates the ratio, as in the case of mitochondrial codon frequencies.
Computer Simulations
The data consist of a pair of codon sequences and are simulated by sampling codon site patterns from the multinomial distribution specified by the site pattern probabilities fij (eq. 5). The sequence has Lc = 100 or 500 codons. Three values of are used: 1 (no bias), 2 (small bias), and 20 (large bias). Most estimates of
from nuclear genes are in the range (1.5, 5), so a value of 2 is typical. Estimates from mitochondrial genes vary considerably among data sets, from 2 or 3 to over 100.
The averages of the estimates among simulated replicates are listed in Table 4
for the three sets of codon frequencies. Standard errors for the ML estimates are also presented, while those for other methods (not shown) are very small due to the use of many more replicates. Averages of dS and dN are calculated for all methods but not shown. We note that the simulation results are highly consistent with those found for infinite data, discussed above. For example, if a method gives estimates smaller than the true value in infinite data, it tends to have negative biases in finite samples as well.
|
The NG method has little bias if codon frequencies are equal and if there is no transition/transversion bias ( = 1). When
1, the method tends to bias toward 1 due to its use of equal weighting in counting sites. The bias is nevertheless small. These results agree well with previous simulations by Ota and Nei (1994)
and Muse (1996)
, who used similar simple models to examine the performance of NG. However, NG is biased in most other parameter combinations. The biases in general agree with findings of the consistency analysis (Fig. 2
). In particular, ignoring the transition/transversion bias leads to underestimates of
and ignoring unequal base frequencies leads to overestimates of
. For equal codon frequencies and no selection (
= 1), NG gives severe underestimates of
when
is large. The effects of the transition/transversion rate bias and base frequency biases tend to cancel each other, such that NG has smaller biases than ML when
= 3 for the mitochondrial codon frequencies. In almost all other cases, ML has smaller biases than NG.
Inas (1995)
method has small biases when base frequencies are equal. The method tends to overestimate
when
< 1 and to underestimate
when
> 1, probably due to its use of equal weighting in counting sites. This is the same pattern as that found in infinite data (Fig. 2AC
). Inas method considerably overestimates
for all values of
under the mitochondrial codon frequencies, probably because it overestimates the number of synonymous sites. For the HIV env codon frequencies, the method overestimates
when
1 and underestimates
when
> 1, as noted for infinite data (Fig. 2GI
).
The new method appears to have little bias over most of the parameter space examined (Table 4
). When
1, it is less biased than NG or Inas (1995)
method. When
= 3, it tends to overestimate
in small samples, like the ML method, but the bias seems smaller than that of ML. The new method appears to provide a close approximation of ML over the range of parameter values examined.
Since all methods are biased for at least some parameter combinations, the mean squared error (MSE) may be an appropriate criterion by which to compare methods. The MSE of a parameter estimator is defined as MSE() = E( - )2, where
is the true value. Since MSE() = Var() + [E() -
]2, this measures both bias and variance. The square root of the MSE is plotted in Figure 3
against the sequence length (number of codons). Two parameter combinations are considered. The first is for mitochondrial genes with a strong transition bias (
= 20) and purifying selection (
= 0.3), and the second is for the HIV env gene with moderate transition bias (
= 2) and positive selection (
= 3). In the first case (Fig. 3A
), ML and the method of this paper performed very similarly, and both have the smallest MSEs, while Inas (1995)
method has very large MSEs due to the positive bias of the method (Fig. 2E
). In the second case (Fig. 3B
), ML performed much worse than other methods for short genes (Lc < 300 codons) due to its large positive bias at
= 3, while for large genes (Lc > 500), ML and the new method are better than NG and Inas method. In both cases, the new method lies between NG and ML.
We also performed a small-scale simulation to examine the effect of sequence divergence level (t). The results are shown in Figure 4 . We examine two sets of parameter values, with the sequence length fixed at Lc = 500. In the first case (Fig. 4A
), equal codon frequencies are used with = 2 and
= 0.3. The new method of this paper overestimates
at small divergences but underestimates
at large divergences. Other methods are insensitive to sequence divergence level. ML and the new method are less biased than NG and Inas method. In the second case (Fig. 4B
), mitochondrial codon frequencies are used with
= 20 and
= 0.3. In this case, ML and the new method have little bias over the whole range of the sequence divergence level. Note that the synonymous rate is quite high, with dS = 0.71 at t = 1, and dS = 1.1 at t = 1.5. NG and Inas method involve positive biases, and the biases become more serious when the sequences are more divergent. Although average estimates of both dN and dS by NG increase with t, dN increases at a faster rate, such that the
ratio increases with the increase of t. Muse (1996)
discussed the fact that at high sequence divergences, NG does not produce distance estimates linear with time.
|
![]() |
Comparison of Human and Orangutan Mitochondrial Genes |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
|
![]() |
Discussion |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
The ML method for pairwise comparison is less biased and has a lower MSE than the approximate methods for almost all parameter combinations. Only for short sequences and high ratios does it involve a positive bias and perform more poorly than some of the approximate methods. We suggest that, in general, the ML method, which accounts for both the transition bias and the codon usage bias, should be the preferred method for estimating dS and dN between two sequences. Only in the case of very short sequences may it be advantageous to use simpler models. In the course of this study, we realized that correcting for biases involved in the NG method is extremely complicated, despite the fact that the method is well known for its simplicity. In contrast, ML is conceptually much simpler, mainly because the probability theory employed by the method takes care of the difficult tasks of weighting evolutionary pathways and correcting for multiple hits, with no need for ad hoc approximations. Specifically, the Chapman-Kolmogorov theorem (e.g., Grimmett and Stirzaker 1992
, p. 239) states that pij(t) =
k pik(s)pkj(t - s) for any 0
s
t; that is, the probability that codon i changes to codon j over time t is a sum over all possible codons (k) at any intermediate time point s. This obvious result ensures that the likelihood calculation (eqs. 46) accounts for all possible pathways of changes between the two codons, weighting them appropriately according to their relative probabilities of occurrence.
The major advantage of ML appears to lie in its flexibility in simultaneous comparison of multiple sequences, taking into account their phylogenetic relationship. Hypotheses concerning variable dN/dS ratios among lineages (Yang 1998
; Yang and Nielsen 1998
) or among sites (Nielsen and Yang 1998
) can be tested using the likelihood ratio test. The ML model can easily be extended to include important features of DNA sequence evolution such as the dependence of nonsynonymous rates on the chemical properties of the amino acids (Yang, Nielsen, and Hasegawa 1998
).
![]() |
Program Availability and Performance |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() |
Acknowledgements |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() |
Footnotes |
---|
1 Keywords: synonymous rate,
nonsynonymous rate,
approximate methods,
maximum likelihood,
molecular evolution,
adaptive evolution,
positive selection.
2 Address for correspondence and reprints: Ziheng Yang, Department of Biology, 4 Stephenson Way, London NW1 2HE, England. E-mail: z.yang{at}ucl.ac.uk
![]() |
References |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Akashi, H. 1995. Inferring weak selection from patterns of polymorphism and divergence at "silent" sites in Drosophila DNA. Genetics 139:10671076.
Comeron, J. M. 1995. A method for estimating the numbers of synonymous and nonsynonymous substitutions per site. J. Mol. Evol. 41:11521159.[ISI][Medline]
Crandall, K. A., and D. M. Hillis. 1997. Rhodopsin evolution in the dark. Nature 387:667668.
Eyre-Walker, A., and P. D. Keightley. 1999. High genomic deleterious mutation rates in hominoids. Nature 397:344347.
Gillespie, J. H. 1991. The causes of molecular evolution. Oxford University Press, Oxford, England.
Gojobori, T. 1983. Codon substitution in evolution and the "saturation" of synonymous changes. Genetics 105:10111027.
Goldman, N. 1993. Statistical tests of models of DNA substitution. J. Mol. Evol. 36:182198.[ISI][Medline]
Goldman, N., and Z. Yang. 1994. A codon-based model of nucleotide substitution for protein-coding DNA sequences. Mol. Biol. Evol. 11:725736.
Grimmett, G. R., and D. R. Stirzaker. 1992. Probability and random processes. 2nd edition. Clarendon Press, Oxford, England.
Hasegawa, M., H. Kishino, and T. Yano. 1985. Dating the human-ape splitting by a molecular clock of mitochondrial DNA. J. Mol. Evol. 22:160174.[ISI][Medline]
Ina, Y. 1995. New methods for estimating the numbers of synonymous and nonsynonymous substitutions. J. Mol. Evol. 40:190226.[ISI][Medline]
. 1996. Pattern of synonymous and nonsynonymous substitutions: an indicator of mechanisms of molecular evolution. J. Genet. 75:91115.[ISI]
Jukes, T. H., and C. R. Cantor. 1969. Evolution of protein molecules. Pp. 21123 in H. N. Munro, ed. Mammalian protein metabolism. Academic Press, New York.
Kimura, M. 1980. A simple method for estimating evolutionary rates of base substitutions through comparative studies of nucleotide sequences. J. Mol. Evol. 16:111120.[ISI][Medline]
. 1983. The neutral theory of molecular evolution. Cambridge University Press, Cambridge, England.
Lee, Y. H., T. Ota, and V. D. Vacquier. 1995. Positive selection is a general phenomenon in the evolution of abalone sperm lysin. Mol. Biol. Evol. 12:231238.[Abstract]
Lewontin, R. 1989. Inferring the number of evolutionary events from DNA coding sequence differences. Mol. Biol. Evol. 6:1532.[Abstract]
Li, W.-H. 1993. Unbiased estimation of the rates of synonymous and nonsynonymous substitution. J. Mol. Evol. 36:9699.[ISI][Medline]
. 1997. Molecular evolution. Sinauer, Sunderland, Mass.
Li, W.-H., C.-I. Wu, and C.-C. Luo. 1985. A new method for estimating synonymous and non-synonymous rates of nucleotide substitutions considering the relative likelihood of nucleotide and codon changes. Mol. Biol. Evol. 2:150174.[Abstract]
Messier, W., and C.-B. Stewart. 1997. Episodic adaptive evolution of primate lysozymes. Nature 385:151154.
Miyata, T., and T. Yasunaga. 1980. Molecular evolution of mRNA: a method for estimating evolutionary rates of synonymous and amino acid substitutions from homologous nucleotide sequences and its applications. J. Mol. Evol. 16:2336.[ISI][Medline]
Moriyama, E. N., and J. R. Powell. 1997. Synonymous substitution rates in Drosophila: mitochondrial versus nuclear genes. J. Mol. Evol. 45:378391.[ISI][Medline]
Muse, S. V. 1996. Estimating synonymous and nonsynonymous substitution rates. Mol. Biol. Evol. 13:105114.[Abstract]
Muse, S. V., and B. S. Gaut. 1994. A likelihood approach for comparing synonymous and nonsynonymous nucleotide substitution rates, with application to chloroplast genome. Mol. Biol. Evol. 11:715724.
Nei, M., and T. Gojobori. 1986. Simple methods for estimating the number of synonymous and nonsynonymous nucleotide substitutions. Mol. Biol. Evol. 3:418426.[Abstract]
Nielsen, R., and Z. Yang. 1998. Likelihood models for detecting positively selected amino acid sites and applications to the HIV-1 envelope gene. Genetics 148:929936.
Ohta, T. 1995. Synonymous and nonsynonymous substitutions in mammalian genes and the nearly neutral theory. J. Mol. Evol. 405663.
Ota, T., and M. Nei. 1994. Variance and covariances of the numbers of synonymous and nonsynonymous substitutions per site. Mol. Biol. Evol. 11:613619.[Abstract]
Pamilo, P., and N. O. Bianchi. 1993. Evolution of the Zfx and Zfy genesrates and interdependence between the genes. Mol. Biol. Evol. 10:271281.[Abstract]
Pedersen, A.-M. K., C. Wiuf, and F. B. Christiansen. 1998. A codon-based model designed to describe lentiviral evolution. Mol. Biol. Evol. 15:10691081.[Abstract]
Perler, F., A. Efstratiadis, P. Lomedica, W. Gilbert, R. Kolodner, and J. Dodgson. 1980. The evolution of genes: the chicken preproinsulin gene. Cell 20:555566.
Stuart, A., K. Ord, and S. Arnold. 1999. Kendalls advanced theory of statistics. Vol. 2a, 6th edition. Arnold, London.
Tateno, Y., N. Takezaki, and M. Nei. 1994. Relative efficiencies of the maximum-likelihood, neighbor-joining, and maximum-parsimony methods when substitution rate varies with site. Mol. Biol. Evol. 11:261277.[Abstract]
Yang, Z. 1994a. Estimating the pattern of nucleotide substitution. J. Mol. Evol. 39:105111.
. 1994b. Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: approximate methods. J. Mol. Evol. 39:306314.
. 1998. Likelihood ratio tests for detecting positive selection and application to primate lysozyme evolution. Mol. Biol. Evol. 15:568573.[Abstract]
. 1999. Phylogenetic analysis by maximum likelihood (PAML). Version 2. University College London, England.
Yang, Z., and R. Nielsen. 1998. Synonymous and nonsynonymous rate variation in nuclear genes of mammals. J. Mol. Evol. 46:409418.[ISI][Medline]
Yang, Z., R. Nielsen, and M. Hasegawa. 1998. Models of amino acid substitution and applications to mitochondrial protein evolution. Mol. Biol. Evol. 15:16001611.