Estimation of Divergence Times for Major Lineages of Primate Species (original) (raw)

Abstract

Although the phylogenetic relationships of major lineages of primate species are relatively well established, the times of divergence of these lineages as estimated by molecular data are still controversial. This controversy has been generated in part because different authors have used different types of molecular data, different statistical methods, and different calibration points. We have therefore examined the effects of these factors on the estimates of divergence times and reached the following conclusions: (1) It is advisable to concatenate many gene sequences and use a multigene gamma distance for estimating divergence times rather than using the individual gene approach. (2) When sequence data from many nuclear genes are available, protein sequences appear to give more robust estimates than DNA sequences. (3) Nuclear proteins are generally more suitable than mitochondrial proteins for time estimation. (4) It is important first to construct a phylogenetic tree for a group of species using some outgroups and then estimate the branch lengths. (5) It appears to be better to use a few reliable calibration points rather than many unreliable ones. Considering all these factors and using two calibration points, we estimated that the human lineage diverged from the chimpanzee, gorilla, orangutan, Old World monkey, and New World monkey lineages approximately 6 MYA (with a range of 5–7), 7 MYA (range, 6–8), 13 MYA (range, 12–15), 23 MYA (range, 21–25), and 33 MYA (range 32–36).

Introduction

In recent years a large number of authors have investigated the evolutionary relationships of primate species, using both molecular and paleontological data, and we now have a rough picture of the phylogenetic relationships of the major lineages of primate species (e.g., Goodman et al. 1998). However, the times of divergence of these lineages are still controversial (e.g., Horai et al. 1995; Takahata and Satta 1997; Arnason, Gullberg, and Janke 1998; Arnason et al. 2000; Cao et al. 2000; Chen and Li 2001). For example, the estimate of the time of divergence between humans and chimpanzees varies from 3.6 MYA (Easteal and Herbert 1997) to 13 MYA (Arnason, Gullberg, and Janke 1998). This controversy has occurred because different authors have used different types of molecular data (e.g., nuclear genes, mitochondrial genes, and noncoding DNA regions), different statistical methods, and different calibration points.

Estimation of divergence time is generally more difficult than reconstruction of a phylogenetic tree, because, strictly speaking, no gene would evolve at a constant rate. For this reason, recent authors have used many independently evolving genes to estimate divergence times in the hope of reducing the effect of rate variation (e.g., Doolittle et al. 1996; Wray, Levinton, and Shapiro 1996; Kumar and Hedges 1998). The traditional method of using information from many different genes is to compute an estimate of divergence time between two species or two groups of species for each gene and then take the average of all the estimates (individual gene [IG] or individual protein [IP] approach). Nei, Xu, and Glazko (2001) showed that this method tends to give biased estimates of divergence times, particularly overestimates when the calibration date is smaller than the time to be estimated. They then proposed the “concatenated distance” method, in which some sorts of concatenated distances for all genes are first computed and the divergence time is then estimated from the distances for all pairs of species. In particular, they suggested the use of a gamma distance for concatenated sequences (CS) for all genes (multigene or multiprotein gamma distance).

It has been customary to use protein sequences rather than DNA sequences for time estimation, because the former are generally more conserved than the latter and can be handled by simpler mathematical models (e.g., Doolittle et al. 1996; Kumar and Hedges 1998; Nei and Kumar 2000, chapter 10). Some authors have argued that because noncoding regions of DNA sequences are not direct targets of natural selection, they should give more reliable estimates (Goodman et al. 1998; Chen and Li 2001). However, noncoding regions are subject to insertion and deletion more often than coding regions, and therefore they may not necessarily give reliable estimates. Nevertheless, for estimating relatively short evolutionary times, as in the present case, DNA sequences both for coding and noncoding regions may be more informative than protein sequences. Because of abundant availability, mitochondrial (mt) DNA have also been used extensively for time estimation in the past (e.g., Horai et al. 1995; Arnason et al. 1996, 1998, 2000). However, the estimates obtained from mt genes are controversial because the evolutionary rate of mt genes apparently varies rather extensively among different groups of mammals (Gissi et al. 2000). We have therefore decided to compare estimates of divergence times obtainable from nuclear protein-coding genes, noncoding DNA sequences, and mt genes.

One of the important factors that determine the accuracy of estimates of divergence times is reliability of the calibration point used for producing the time scale of the phylogenetic tree constructed. In this study we use the times of divergence between humans and orangutans (about 13 MYA) and between primates and artiodactyls (90 MYA) as calibration points. We are interested in finding whether these two calibration points give similar time estimates for other branch points in the tree.

Materials and Methods

Species and Genes Used

Our intention in this study was to use as many genes as possible for a group of species under consideration. However, the number of gene sequences available for primate species is quite limited except for a few primate species. Furthermore, to estimate divergence times, we needed species that would provide reasonably good calibration points. One of the calibration points we used, as already mentioned, is the time of divergence between the orangutan and the human lineages (13 MYA). This fossil dating based on Sivapithecus was once questioned (Pilbeam et al. 1990), but recent statistical analyses of cranial and postcranial characters (Begun, Ward, and Rose 1997; Ward 1997) suggest that the “_Sivapithecus_-Pongo (orangutan) clade remains the strongest phylogenetic hypothesis” (Ward 1997). We therefore decided to use this calibration point, and for this reason inclusion of orangutan genes was essential. Another fossil-based calibration point we used was the time of divergence between archaic fossil ungulates (clade Ungulatomorpha) and primates. The fossils indicate that Ungulatomorpha appeared about 85–90 MYA (Archibald 1996; Archibald, Averianev, and Ekdale 2001). Because fossils generally give a minimum estimate of splitting time, we assumed that primates and artiodactyls diverged 90 MYA and used primarily cattle (Bos taurus) genes for artiodactyls.

For estimating the divergence times for a group of species, we have to have outgroups to determine the root of the tree for the group (fig. 1). For this purpose, we used mice and rats (see Results for the justification of using rodents as outgroups). However, the number of genes available varied considerably with species. We therefore conducted separate analyses for the following four species groups.

Species Group 1

For this group, we could use 29 shared nuclear genes. The primate species used were humans (Homo sapiens), chimpanzees (Pan troglodytes), gorillas (Gorilla gorilla), and orangutans (Pongo pygmaeus), and their topological relationships are given in figure 1_A_. Artiodactyls (primarily cattle genes) were used for calibrating evolutionary times, whereas mice (Mus musculus) and rats (Rattus norvegicus) were used for determining the root of the tree for the remaining species. The names of the genes used are listed in table S1 of the online Supplementary Material and the Web site http://www.bio.psu.edu/people/faculty/nei/lab/databases.htm. The smallest number of genes available from GenBank was for orangutans, and we compiled all genes that are orthologous to the orangutan genes. To avoid paralogous genes, we constructed a Neighbor-Joining (NJ) tree (Saitou and Nei 1987) for each gene using uncorrected p distance (Nei and Kumar 2000, p. 18) and eliminated all the genes that produced the incorrect topology except for humans, chimpanzees, and gorillas, or the genes that showed zero distances for all pairs of humans and African ape species. Apparently because humans, chimpanzees, and gorillas diverged in a short period of evolutionary time and the ancestral populations were polymorphic, different genes are known to show different topologies for the three species (Saitou and Nei 1986; Satta, Klein, and Takahata 2000; Chen and Li 2001; Klein and Takahata 2002; O'hUigin et al. 2002). Therefore, we cannot eliminate paralogous genes for these species by comparing the gene tree with the species tree. In this case we used the most plausible orthologous genes, although there were only a few such cases. The final number of genes chosen in this way was 29 protein-coding genes with a total of 6,966 codons.

Species Group 2

For estimating the divergence times for Old World (OW) monkeys (mostly macaque [_Macaca mulatta_] genes) and New World (NW) monkeys (mostly marmoset [_Callithrix jacchus_] genes), we could obtain only 13 nuclear genes with 2,425 codons. Therefore, we conducted a separate analysis for the nine species listed in figure 2_A_. The genes used in this study were a subset of the genes used for species group 1 (see online Supplementary Material).

Species Group 3

Species group 3 was chosen primarily for studying time estimates obtainable from noncoding regions of DNA sequences. In this group we estimated the times of divergence of the human lineage from the chimpanzee, gorilla, orangutan, and OW monkey (macaque) lineages and used NW monkeys (marmoset) as the outgroup (fig. 3). The noncoding DNA regions used were the flanking and intergenic regions of the ϵ and the γ1-γ2 globin genes (9,818 bp) (http://cmmg.biosci.wayne.edu/lgross/). The noncoding sequences of the globin gene regions are available even for loris, lemurs, and others (Goodman et al. 1998). However, because the comparison of sequences from distantly related species showed many deletions and insertions, we did not include these species in this study.

Species Group 4

For some unknown reasons, the evolutionary rate of mt genes varies considerably from species to species (Gissi et al. 2000), and therefore these genes may not give reliable estimates of divergence times. However, a large number of authors have used these genes for studying primate evolution. We have therefore examined time estimates obtainable from mt genes. In this study we used humans, chimpanzees, gorillas, orangutans, gibbons, baboons (Papio hamadryas; OW monkeys), capuchins (Cebus albifrons; NW monkeys), slow loris (Nycticebus coucang; strepsirhines), artiodactyls (Bos taurus), and two rodent species (fig. 4_A_). Each of these species has 13 coding genes, and we used all the genes, regardless of the transcription direction. Some authors (e.g., Arnason et al. 1996, 2000) excluded genes NAD6 and COX2 because of the difference in transcription direction or a higher rate of evolution in some species groups. However, our preliminary study showed that inclusion or exclusion of these genes has little effect on time estimates.

Statistical Methods

As already mentioned, we are primarily interested in the CS approach in this study. When protein sequences were used, we first constructed a NJ tree using MEGA2 (Kumar et al. 2001) with multiprotein gamma distance (_d_MG), which is a Poisson-correction (PC) gamma distance obtained for the concatenated amino acid sequences for all the proteins used. The gamma parameter a was estimated by Gu and Zhang's (1997) method. In practice, however, the a value obtained by this method is sometimes too small for the purpose of time estimation (Nei, Xu, and Glazko 2001), and the classical Dayhoff distance, which can be obtained by a PC gamma distance with a = 2.25 (Nei and Kumar 2000, p. 23), often gives reasonable time estimates when the divergence time considered is relatively short. Actually, it is known that even PC distance with a = ∞ gives reasonable time estimates when conserved proteins such as cytochrome c and hemoglobins are used (Dickerson 1971). Because the proteins we used were quite conserved, we used Dayhoff and PC distances as well. When we analyzed DNA sequences, the Jukes-Cantor, Kimura, and Kimura gamma distances (Nei and Kumar 2000, chapter 3) were used. Although the phylogenetic tree of the primates species used here is reasonably well established (Goodman et al. 1998), we also constructed maximum parsimony (MP) trees using the branch-and-bound algorithm of MEGA2 and maximum likelihood (ML) trees using the Poisson model of PROTML (Adachi and Hasegawa 1996) for protein data and PAUP* (Swofford 1998) for DNA data. The topologies of these trees were always the same as those of NJ trees.

Once the topology of the species was determined, the branch lengths of the tree were estimated by the least squares method. We then used Takezaki, Rzhetsky, and Nei's (1995) two-cluster and branch-length tests to examine the molecular clock hypothesis (computer program LINTREE; see http://mep.bio.psu.edu). Only when these tests were significant at the 1% level did we consider the deviation biologically meaningful, because reasonably good time estimates are known to be obtained even if the deviation is considerably large (Nei and Kumar 2000, chapter 10). When the molecular clock hypothesis was acceptable, we constructed a linearized tree to estimate the times of species divergence. When it was rejected, we used the stem-lineage method proposed by Nei, Xu, and Glazko (2001) and Nei and Glazko (2002).

As mentioned earlier, the traditional method of time estimation is first to construct a linearized tree for each gene with a gamma distance and estimate the divergence times for this tree. Let us consider the linearized tree in figure 1_B_ and assume that the divergence time (_T_1) between primates and artiodactyls is known (90 MYA). The divergence time (t) between humans and chimpanzees can then be estimated by

where a and b are branch length estimates for the human and artiodactyl lineages, respectively, in figure 1_B_. If t̂ is computed for all genes examined, t is estimated by their average (t̄). That is,

where i stands for the _i_th gene and k is the number of genes used. If a and b are normally distributed and k is large, the distribution of t̄ is known to be given by the following Cauchy distribution:

(Johnson and Kotz 1970, p. 154). This has a wider distribution than the normal distribution, and the theoretical variance is infinite (fig. 5). One way to avoid the problem of large variances is to eliminate so-called outliers (say 5% of uppermost and 5% lowermost t̂ i values) (Kumar and Hedges 1998). In this approach, however, it is unclear how t̂ is affected by the elimination, and some subjective judgment may enter into the computation. Nei, Xu, and Glazko (2001) and Nei and Glazko (2002) also showed that t̄ tends to be an overestimate when the divergence time to be estimated is older than the calibration point.

By contrast, if we use a multigene gamma distance for concatenated sequences of many genes and estimate _t_1 by equation (1) using a and b obtained from multigene distances, the bias inherent in equation (2) virtually disappears. The standard errors of _t̂_1 is computed by the bootstrap test using genes (rather than amino acids or nucleotides) as the units of resampling. This computation can be made by our computer program TIMER (http://mep.bio.psu.edu).

Results

Species Group 1

We first analyzed protein sequences following Nei, Xu, and Glazko (2001). The phylogenetic tree presented in figure 1_A_ is an NJ tree obtained by using multiprotein gamma distance (_d_MG) with a = 0.47. The topology of this tree is well supported by the bootstrap test. It was also supported by MP and ML analyses. The topology of the primate portion of the tree is identical with the generally accepted molecular topology (Goodman et al. 1998). However, the remaining portion of the tree is somewhat controversial. There is no fossil record that resolves the evolutionary relationships of primates, artiodactyls, and rodents (Bromham, Phillips, and Penny 1999). Most molecular studies in the past have supported a close relationship between primates and artiodactyls rather than between primates and rodents or between artiodactyls and rodents (e.g., Li et al. 1990; Cao et al. 1998; Arnason et al. 2000; Reyes, Pesole, and Saccone 2000). Recently, however, Murphy et al. (2001) constructed a phylogenetic tree of about 42 placental mammalian species using 19 nuclear and 3 mitochondrial genes and suggested that primates and rodents are phylogenetically closer to each other than to artiodactyls. However, our phylogenetic analysis of 71 nuclear proteins with 24,952 amino acids for humans, artiodactyls, rodents, chicken, and Xenopus has generated the traditional molecular topology with a high level of statistical support (Nei and Glazko 2002). We have therefore decided to use the phylogenetic tree given in figure 1 for our study of divergence times of hominoid species. For estimating divergence times for primate species, however, the topological relationships among primates, artiodactyls, and rodents do not matter very much, because the linearized tree is constructed without rodent species.

To estimate the divergence times, we first tested the molecular clock hypothesis by using Takezaki, Rzhetsky, and Nei's (1995) method. The results of this test showed that the deviation from the clock hypothesis is not statistically significant at the 1% level. We therefore constructed the linearized tree (fig. 1_B_). The evolutionary time scale given for this tree was obtained by using the divergence time between primates and artiodactyls (_T_1 = 90 MYA) as the calibration point and a rate of amino acid substitution of 1.2 × 10−9 per year per lineage. Estimates of the divergence times for humans and ape species are presented in table 1. The estimates obtained by using the human/orangutan divergence (_T_2 = 13 MYA) as the calibration point are also presented in table 1. It is interesting to note that the two sets of estimates are close to each other and that the estimate of divergence time between humans and chimpanzees (about 5.5–5.7 MYA) is close to the ages of the recently discovered oldest hominid fossils (5.4–7.0 MYA; Aiello and Collard 2001; Haile-Selassie 2001; Brunet et al. 2002). When the second calibration point (_T_2 = 13 MYA) is used, we obtained 93 MYA for the divergence time between primates and artiodactyls. This estimate is close to _T_1 = 90 MYA. Table 1 also includes estimates obtained by using Dayhoff and PC distances. These distances give reasonably good time estimates for humans and apes if we assume that the lower limit of the human/chimpanzee divergence time is about 5.4 MYA and the upper limit of the primate/artiodactyl divergence time is about 90 MYA. These results show that when a ranges from 0.47 to ∞ the time estimates remain nearly the same for these relatively closely related species. (_d_MG with a = 0.47 is better for estimating the primate/artiodactyls divergence time.)

Table 1 includes the time estimates obtained by the IG approach with gamma distances. In this case the time estimate for the human/chimpanzee divergence is apparently too small when _T_1 = 90 MYA was used as the calibration point, and the estimate for the primate/artiodactyl divergence is too large when _T_2 = 13 MYA was used (75% overestimate of the calibration date). This type of underestimation and overestimation occurs even when Dayhoff distance or PC distance is used.

In the present data set, the extent of protein divergence among the human and ape species was rather small, so we suspected that DNA sequences might give more reliable results. We therefore estimated the divergence times using the concatenated DNA sequences for the species in figure 1_A_, for which the DNA sequences of 24 genes were available. (No DNA sequences were available for five genes for some of the species used.) In this study we first used the Kimura gamma distance for the entire sequences (18,272 bp), and then for the sequences of first and second codon positions only. The results obtained are presented in table 2. When we used all three codon positions of DNA sequences and _T_1 = 90 MYA as the calibration point, the estimates of divergence times among hominoid species appeared to be too low, but the primate/artiodactyl divergence time was apparently overestimated when _T_2 = 13 MYA was used as the calibration point. This tendency did not change, even when we used the sequence data for first and second codon positions and Kimura gamma distance. These results suggest that protein data generally give more reliable estimates than DNA data even for closely related species.

Interestingly, however, DNA data gave reasonably good estimates of the human/orangutan and the primate/artiodactyl divergence times when Kimura distance rather than Kimura gamma distance was used for first and second codon position data. This finding might suggest that the parameter for the gamma distance a is underestimated for closely related sequences and that Kimura distance for first and second codon positions is a good option in this case. However, this is somewhat illogical, because there must be rate heterogeneity for concatenated sequences. Therefore, we are not sure whether we should give much weight to this result.

Species Group 2

The time estimates obtained by multiprotein gamma distance (a = 0.61) with _T_1 = 90 MYA and _T_2 = 13 MYA for this species group are presented in table 3. They are similar to each other and are rather close to the estimates in table 1, whenever comparable estimates are available. The estimate for the divergence between humans and orangutans (12.6 MYA) is also close to the paleontological estimate (13 MYA). Similarly, the estimate for the divergence time between primates and artiodactyls (92.8 MYA) is close to the paleontological data (90 MYA). Actually, these statements hold true even with Dayhoff and PC distances, although the primate/artiodactyls divergence time tends to be underestimated when _T_2 = 13 MYA is used. The average estimates of the time of divergence between humans and OW monkeys and NW monkeys are approximately 23 MYA and 33 MYA, respectively. These estimates are close to the rough estimates obtained by Goodman et al. (1998).

By contrast, the IG approach again gives unduly low estimates for the human/chimp divergence but gives an unduly high estimate for the primates/artiodactyl divergence, as expected theoretically (Nei, Xu, and Glazko 2001; Nei and Glazko 2002). When _T_2 = 13 MYA is used, the estimates of the time of divergence of humans from OW monkeys and NW monkeys are also considerably higher than those obtained by the CS approach.

We also used all three codon position data and first and second codon positions of DNA sequences to estimate divergence times. In this case we could use only nine genes. The results were quite similar to those presented in table 2. That is, when Kimura distance or Kimura gamma distance for all three codon positions was used, the calibration point of _T_1 = 90 MYA gave too low estimates for the human/chimp and the human/gorilla divergence, whereas _T_2 = 13 MYA gave a too high estimate of the primate/artiodactyl divergence (table S1 of the online Supplementary Material). By contrast, when Kimura distance was used for first and second codon position data, both _T_1 = 90 MYA and _T_2 = 13 MYA gave reasonable estimates. However, the estimates obtained were again similar to those obtained from multiprotein gamma distance with a = 0.61.

Species Group 3

The phylogenetic tree for the noncoding regions of the ϵ and the γ1-γ2 gene clusters is presented in figure 3. When the marmoset (NW monkey) was used as the outgroup, the evolutionary change of hominoid and macaque sequences did not deviate significantly from the molecular clock. We therefore constructed a linearized tree for the hominoids and OW monkeys and estimated the divergence times for these species using the human/orangutan divergence time as the calibration point (table 4). When Kimura gamma distance with a = 0.26 was used, the estimated times of divergence of humans from chimpanzees and gorillas were slightly smaller than those obtained by multiprotein gamma distance (tables 1 and 3), but the estimate time for OW monkeys was slightly higher. When Kimura or Jukes-Cantor distance was used, the results hardly changed. However, this study is not very informative, because we could not use the calibration point of _T_1.

Species Group 4

Figure 4_A_ shows the NJ tree for eight primate and three nonprimate species obtained by 13 mtDNA-encoded proteins with 3752 amino acids. Essentially the same tree topology was obtained by MP and ML analyses. In this case the molecular clock obviously does not work, and the orangutan, baboon, and capuchin monkey sequences evolved significantly faster than the average sequence (1% level), whereas the slow loris and cattle sequences evolved significantly slower.

Therefore, the linearization of the tree cannot be justified unless we eliminate all deviant species. However, if we eliminate the deviant species, we have no more calibration points. So, despite this clear violation of the molecular clock, we attempted to construct a linearized tree using all species. We again used the primate/artiodactyl divergence time (_T_1 = 90 MYA) and the human/orangutan divergence times (_T_2 = 13 MYA) as the calibration points (table 5). The results were unexpectedly interesting, because when _T_1 = 90 MYA was used, the estimates of times of divergence of the human lineage from the chimpanzee, gorilla, orangutan, gibbon, baboon, capuchin, and slow loris lineages were approximately 11, 15, 32, 33, 63, 90, and 90 MYA, and these estimates are very similar to those obtained by Arnason, Gullberg, and Janke (1998) (13, 16, 30, 35, 52, 70, and 90 MYA, respectively). Arnason, Gullberg, and Janke (1998) actually used some kind of rate-adjustment method, because they were aware of the slow rate of evolution of the artiodactyl sequence, but their rate adjustment was probably insufficient, because we obtained essentially the same results without any rate adjustment. The inadequacy of the linearized tree method in this case is also clear from the fact that the use of the calibration point _T_2 = 13 MYA gives estimates very different from those obtained by using _T_1 = 90 MYA.

Nei, Xu, and Glazko (2001) and Nei and Glazko (2002) suggested that the stem-lineage method for estimating divergence times may be used when the evolutionary rate varies with exterior branch and the number of amino acids or nucleotides used is very large. In this method the stem lineage that keeps generating exterior branches is assumed to evolve at a constant rate, and the divergence times are estimated by using only the stem lineage. As in figure 4_A_, denoting the lengths of the exterior and interior branches of the primate stem lineage by a, b, c, d, e, f, g, h, and i, we can assume that the sum of stem branch lengths, S = (a + b)/2 + c + d + e + f + g + h + i, corresponds to the primate/cattle divergence time, i.e., _T_1 = 90 MYA, and that the divergence time for a pair of species is proportional to the appropriate sum of branch lengths in the stem lineage. For example, we can estimate the divergence time between humans and baboons by [(a + b)/2 + c + d + f] _T_1/S. In the case of _d_MG with a = 0.51 we have a = 0.022, b = 0.024, c = 0.007, d = 0.026, e = 0.010, f = 0.031, g = 0.059, h = 0.075, i = 0.006, and S = 0.237. Therefore, we obtain 37.2 MYA. By contrast, if we use _T_2 = 13 MYA as the calibration point, the human/baboon divergence time is estimated by [(a + b)/2 + c + d + _f_] _T_2/S, where _T_2 = 13 MYA and S = (a + b)/2 + c + d = 0.056. Therefore, we have 22.5 MYA.

The estimates of divergence times obtained for other species are presented in table 5. The time estimates obtained by using _T_1 = 90 MYA are again much higher than those in tables 1 and 3, but those obtained by using _T_2 = 13 MYA are rather close to those in the latter tables except for artiodactyls. Nevertheless, because the stem-lineage method depends on an unproven assumption, we had better not to give much weight to these estimates.

We conducted a similar statistical analysis for the DNA sequences of 13 coding genes using first, second, and third codon positions. The time estimates obtained were even more divergent from those given in tables 1 and 2, and the estimates obtained under the assumptions of _T_1 = 90 MYA and _T_2 = 13 MYA were even more inconsistent than those for protein data (data not shown). These results again suggest that mt DNA sequence data are less suitable for time estimation than nuclear proteins.

Discussion

In this article we have presented several problems related to the estimation of divergence times. First, we have shown that the IG method of time estimation often gives biased estimates of divergence times, as was shown theoretically by Nei, Xu, and Glazko (2001) and Nei and Glazko (2002). This method tends to give overestimates of divergence times when the calibration point is younger than the estimated time and to give underestimates when the calibration point is older than the estimated time. A similar conclusion was obtained by Rodriguez-Trelles, Tarrio, and Ayala (2002) in their computer simulation. Our observations support this theoretical prediction. These biases are caused mainly by the variances and covariances of a and b in equation (1). Therefore, if we use the concatenated sequences of many genes, these variances and covariances are reduced, and consequently the estimation bias is reduced, as is clear from tables 1 and 3. For this reason, we recommend that the CS method rather than the IG method be used in time estimation.

For estimating divergence times for distantly related organisms, it has been customary to use protein sequences rather than DNA sequences, because DNA sequences are usually less conserved than protein sequences and the substitution pattern in DNA sequences varies extensively with codon position (Nei and Kumar 2000, chapters 2 and 3). In the present study, PC gamma distance for nuclear proteins gave reasonably good time estimates for simian primates, and the gamma parameter value (a) did not affect the results seriously as long as a was greater than the estimated one. In the case of DNA sequences the estimates depended on the codon positions and the a value used. In general, Kimura distance with a = ∞ for first and second codon position data gave reasonable estimates, but Kimura gamma distance with the estimated a value for all three codon positions or first and second codon positions did not. Therefore, it appears that protein sequences give more robust estimates than DNA sequences, even for relatively closely related species as long as many genes are used. However, this problem should be studied in more detail considering both distantly related and closely related species.

Many authors have used mt genes for time estimation. The present study shows that even for relatively closely related species such as simian primates these genes do not give good results because the evolutionary rate varies extensively among different evolutionary lineages. It is therefore preferable to use nuclear genes rather than mitochondrial genes for time estimation.

As the genome sequencing project proceeds in various organisms, we will have many genes that can be used for time estimation. However, because the number of gene sequences available usually varies from species to species, the number of genes shared by all species is often rather small. For this reason, some authors (e.g., Stauffer et al. 2001) have used only three species or three groups of species at a time for time estimation without using outgroups. This approach certainly increases the number of genes available, but the molecular clock cannot be tested properly unless the root of the three-species tree is determined by using outgroups (Takezaki, Rzhetsky, and Nei 1995). Therefore, this method may not give accurate estimates even if a large number of genes is used. The different sets of genes used for different triplet of species or species groups may also give different time estimates. It is generally advisable first to construct a phylogenetic tree for all the species involved and then estimate divergence times.

In the present article we used two calibration dates for estimating divergence times to see whether the two different dates give similar estimates or not. This investigation was useful in identifying the right kind of molecular data and the right kind of statistical methods. However, once we know the best data set or the best statistical method, we can now use the two or more calibration dates to obtain a more reliable evolutionary rate by using the regression method, as was done by Hughes and Nei (1990) and Takahashi, Rooney, and Nei (2000). If the calibration dates are reliable, this method would give more reliable time estimates. In practice, however, there are cases in which one calibration point is more reliable than others (Kumar and Hedges 1998). If this is the case, use of one or two reliable calibration points may be preferable.

A number of authors (e.g., Sanderson 1997; Rambaut and Bromham 1998; Kishino, Thorne, and Bruno 2001; Soltis et al. 2002) have used sophisticated statistical methods to take care of rate variation in the presence of multiple calibration dates. However, even these methods do not seem to work well when the extent of rate variation is large (Soltis et al. 2002). It appears that the most important thing for obtaining reliable time estimates is to use molecular data that follow the molecular clock more closely than others.

Some authors estimated divergence times by using several local evolutionary rates in different parts of the tree (e.g., Yoder and Yang 2000). If we know local evolutionary rate very accurately, this approach is expected to give reliable time estimates. In practice, however, the local (relative) rates are determined intuitively by looking at the branch lengths of the original phylogenetic tree. Therefore, the results obtained are often different depending on the number of local rates and the calibration points used. For example, Yoder and Yang (2000) estimated the human/chimpanzee divergence time using mitochondrial genes from 31 mammalian species. They used different substitution models, different data types (protein sequences, and different codon positions of DNA sequences), different numbers of local clocks, and different calibration points and obtained various estimates ranging from 2.68 to 10.12 MYA. It was therefore difficult to choose the most likely estimate from the statistical analysis alone, and they chose estimates that were most likely to agree with the available fossil record. We also analyzed our mt gene data using the Yoder and Yang method with two or four local clocks, but the results were no better than those obtained by the stem-lineage method shown in table 5 (see online Supplementary Material). This again suggests that mt gene data are not appropriate for time estimation in primates, whatever method is used. What is important in time estimation is to use genes that follow the global molecular clock as much as possible.

Note that the estimation of divergence times from molecular data is not to fit molecular data to the fossil record available. Fossil records are usually very poor in providing divergence time estimates as mentioned below, and the utility of molecular clocks is to provide time estimates that are difficult to obtain from the fossil record. Therefore, a global clock that applies at least to a group of species is necessary.

In this article, we used several different data sets and distance measures each with a single global clock. We obtained relatively close but slightly different time estimates using different distance measures. For example, our estimate of the time of divergence between humans and chimpanzees obtained by the CS method varied from 5.5 MYA to 7.4 MYA in table 1 and table 3 (excluding the standard errors), the average of the 12 observations being 6.3 MYA. (Although the estimates in table 3 are based on a subset of the genes used in table 1, we treated them as independent estimates because the estimates depend on the genes and distance measures used and we wanted to know only a crude magnitude of variation of the estimates without consideration of standard errors.) Therefore, we can probably say that the divergence between humans and chimpanzees occurred about 6 MYA with a rough range of 5–7 MYA. If we use this crude approach, the times of divergence of the human lineage from the gorilla, orangutan, OW monkey, and NW monkey lineages become 7 MYA (range, 6–8), 13 MYA (range, 12–15), 23 MYA (range, 21–25), and 33 MYA (range, 32–36). Here we included the fossil estimate (13 MYA) for the computation of the human/orangutan divergence and used only the estimates in table 3 for the computation of the human/OW monkey and the human/NW monkey divergence.

Note that the above estimates were obtained without consideration of uncertainty of fossil dating. We used _T_1 = 90 MYA for the primate/artiodactyl divergence, but the actual dating of Ungulatomorpha varies from 85 MYA to 90 MYA. The dating of Sivapithecus also varies from 6.8 MYA to 12.7 MYA (Ward 1997). Furthermore, these dates do not necessarily indicate the actual time of species divergence (Easteal 1999). Therefore, the actual time of divergence may deviate even more from our estimates. In the presence of this uncertainty, what kind of estimates should we trust? In our opinion, the best way would be to construct linearized trees for a group of species (many different species of primates in the present case) using several different groups of genes and examine the consistency among time estimates obtained from different sets of genes. If different genes give similar estimates, we can accept them until they are rejected by other new sets of genes.

If we know the uncertainty of calibration points, it is clear that the standard error computed here is only a small portion of the uncertainty of time estimates. The magnitude of a standard error also depends on the method of computation. In this study, we computed the standard error of a time estimate from concatenated sequences using genes as the units of bootstrap resampling. Theoretically, it is possible to compute this standard error using amino acids (or nucleotides) as the units of resampling. However, the latter method gives an unduly small standard error, because the unit of evolution is a gene rather than an amino acid. In the IG method the magnitude of standard error is determined in part by the extent of elimination of outliers as mentioned earlier. Because the extent of elimination of outliers is subjective, the reliability of standard errors is difficult to evaluate. For these reasons, the standard error attached to a time estimate does not give a real extent of uncertainty, and we should not place much emphasis on it.

Jeffery Long, Associate Editor

Fig. 1.

(A) Neighbor-Joining tree for hominoid 4 species and artiodactyls constructed by using multiprotein gamma distance (_d_MG) with a = 0.47 for 29 protein sequences. Two rodent species were used as outgroups. (B) Linearized tree of the above tree. The timescale does not apply to rodents, because these species are outgroups

Fig. 2.

(A) Neighbor-Joining tree for simian primates and artiodactyls constructed by using multiprotein gamma distance with a = 0.61 for 13 protein sequences. (B) Linearized tree of the above tree. The time scale does not apply to rodents, because these species are outgroups

Fig. 3.

Neighbor-Joining tree for simian primates constructed by using Kimura distance for the noncoding regions of the ϵ and the γ1-γ2 genes

Fig. 4.

(A) Neighbor-Joining tree for simian and prosimian primates and non-primate mammalian species obtained by multiprotein gamma distances (_d_MG) for 13 mitochondrial proteins (3,752 amino acids) with a = 0.51. *Species that evolved significantly faster the average. †Species that evolved significantly more slowly. Two rodent species were used as outgroups. (B) Enforced linearized tree of tree A.

Fig. 5.

Cauchy distribution compared with the normal distribution. p(t̂) refers to the probability density of t̂

Table 1

Estimates (± Standard Errors) of Divergence Times of the Human Lineage from Other Primate Species and Artiodactyls Obtained by Using Different Approaches (29 Nuclear Proteins).

Table 1

Estimates (± Standard Errors) of Divergence Times of the Human Lineage from Other Primate Species and Artiodactyls Obtained by Using Different Approaches (29 Nuclear Proteins).

Table 2

Estimates of Divergence Times (± Standard Errors) of the Human Lineage from Other Primate Species and Artiodactyls Obtained by Using Nucleotide Sequences of 24 Genes.

Table 2

Estimates of Divergence Times (± Standard Errors) of the Human Lineage from Other Primate Species and Artiodactyls Obtained by Using Nucleotide Sequences of 24 Genes.

Table 3

Estimates of Divergence Times (± Standard Errors) of the Human Lineage from Other Primate Species and Artiodactyls (13 Nuclear Proteins).

Table 3

Estimates of Divergence Times (± Standard Errors) of the Human Lineage from Other Primate Species and Artiodactyls (13 Nuclear Proteins).

Table 4

Estimates (± Standard Errors) of Divergence Times (MYA) of the Human Lineage from Other Primate Species Obtained from Noncoding DNA Sequences of the ϵ and γ1-γ2 Globin Gene Regions.

Note.—The total number of nucleotides used is 9,818. The human/orangutan divergence time was used as the calibration point.

Table 4

Estimates (± Standard Errors) of Divergence Times (MYA) of the Human Lineage from Other Primate Species Obtained from Noncoding DNA Sequences of the ϵ and γ1-γ2 Globin Gene Regions.

Note.—The total number of nucleotides used is 9,818. The human/orangutan divergence time was used as the calibration point.

Table 5

Estimates (± Standard Errors) of Divergence Times (MYA) of the Human Lineage from Other Primate Species and Artiodactyls Obtained by Using Multiprotein Gamma Distances (13 Mitochondrial Proteins with a = 0.51).

Table 5

Estimates (± Standard Errors) of Divergence Times (MYA) of the Human Lineage from Other Primate Species and Artiodactyls Obtained by Using Multiprotein Gamma Distances (13 Mitochondrial Proteins with a = 0.51).

Literature Cited

Adachi, J., and M. Hasegawa.

1996

. MOLPHY: programs for molecular phylogenetics. Institute of Statistical Mathematics, Tokyo.

Aiello, L. C., and M. Collard.

2001

. Palaeoanthropology: our newest oldest ancestor?

Nature

410

:

526

-527.

Archibald, J. D.

1996

. Fossil evidence for a late Cretaceous origin of “hoofed” mammals.

Science

272

:

1150

-1152.

Archibald, J. D., A. O. Averianov, and E. G. Ekdale.

2001

. Late cretaceous relatives of rabbits, rodents, and other extant eutherian mammals.

Nature

414

:

62

-65.

Arnason, U., A. Gullberg, A. S. Burguete, and A. Janke.

2000

. Molecular estimates of primate divergences and new hypotheses for primate dispersal and the origin of modern humans.

Hereditas

133

:

217

-228.

Arnason, U., A. Gullberg, S. Gretarsdottir, B. Ursing, and A. Janke.

2000

. The mitochondrial genome of the sperm whale and a new molecular reference for estimating eutherian divergence dates.

J. Mol. Evol.

50

:

569

-578.

Arnason, U., A. Gullberg, and A. Janke.

1998

. Molecular timing of primate divergences as estimated by two nonprimate calibration points.

J. Mol. Evol.

47

:

718

-727.

Arnason, U., A. Gullberg, A. Janke, and X. Xu.

1996

. Pattern and timing of evolutionary divergences among hominoids based on analyses of complete mtDNAs.

J. Mol. Evol.

43

:

650

-661.

Begun, D. R., C. V. Ward, and M. D. Rose.

1997

. Events in hominoid evolution. Pp. 389–415 in D. R. Begun, C. V. Ward, and M. D. Rose, eds. Function, phylogeny and fossils: Miocene hominoid evolution and adaptation. New York, Plenum.

Bromham, L., M. J. Phillips, and D. Penny.

1999

. Growing up with dinosaurs: molecular dates and the mammalian radiation.

Trends Ecol. Evol.

14

:

113

-118.

Brunet, M., F. Guy, D. Pilbeam, H. T. Mackaye, A. Likius, D. Ahounta, A. Beauvilain, C. Blondel, H. Bocherens, and J. R. Boisserie, et al. (39 co-authors).

2002

. A new hominid from the Upper Miocene of Chad, Central Africa.

Nature

418

:

145

-151.

Cao, Y., M. Fujiwara, M. Nikaido, N. Okada, and M. Hasegawa.

2000

. Interordinal relationships and timescale of eutherian evolution as inferred from mitochondrial genome data.

Gene

259

:

149

-158.

Cao, Y., A. Janke, P. T. Waddell, M. Westerman, O. Takenaka, S. Murata, N. Okada, S. Pääbo, and M. Hasegawa.

1998

. Conflict among individual mitochondrial proteins in resolving the phylogeny of eutherian orders.

J. Mol. Evol.

47

:

307

-322.

Chen, F. C., and W.-H. Li.

2001

. Genomic divergences between humans and other hominoids and the effective population size of the common ancestor of humans and chimpanzees.

Am. J. Hum. Genet.

68

:

444

-456.

Dickerson, R. E.

1971

. The structures of cytochrome c and the rates of molecular evolution.

J. Mol. Evol.

1

:

26

-45.

Doolittle, R. F., D.-F. Feng, S. Tsang, G. Cho, and E. Little.

1996

. Determining divergence times of the major kingdoms of living organisms with a protein clock.

Science

271

:

470

-477.

Easteal, S.

1999

. Molecular evidence for the early divergence of placental mammals.

Bioessays

21

:

1052

-1058.

Easteal, S., and G. Herbert.

1997

. Molecular evidence from the nuclear genome for the time frame of human evolution.

J. Mol. Evol.

44

:(Suppl.):

121

-132.

Gissi, C., A. Reyes, G. Pesole, and C. Saccone.

2000

. Lineage-specific evolutionary rate in mammalian mtDNA.

Mol. Biol. Evol.

17

:

1022

-1031.

Goodman, M., C. A. Porter, J. Czelusniak, S. L. Page, H. Schneider, J. Shoshani, G. Gunnell, and C. P. Groves.

1998

. Toward a phylogenetic classification of primates based on DNA evidence complemented by fossil evidence.

Mol. Phylogenet. Evol.

9

:

585

-598.

Gu, X., and J. Zhang.

1997

. A simple method for estimating the parameter of substitution rate variation among sites.

Mol. Biol. Evol.

14

:

1106

-1113.

Haile-Selassie, Y.

2001

. Late Miocene hominids from the Middle Awash, Ethiopia.

Nature

412

:

178

-181.

Horai, S., K. Hayasaka, R. Kondo, K. Tsugane, and N. Takahata.

1995

. Recent African origin of modern humans revealed by complete sequences of hominoid mitochondrial DNAs.

Proc. Natl. Acad. Sci. USA

92

:

532

-536.

Hughes, A. L., and M. Nei.

1990

. Evolutionary relationships of class II major-histocompatibility-complex genes in mammals.

Mol Biol Evol.

7

:

491

-514.

Johnson, N., and S. Kotz.

1970

. Continuous univariate distributions—1. Houghton Mifflin, Boston.

Kishino, H., J. L. Thorne, and W. J. Bruno.

2001

. Performance of a divergence time estimation method under a probabilistic model of rate evolution.

Mol. Biol. Evol.

18

:

352

-361.

Klein, J., and N. Takahata.

2002

. Where do we come from? The molecular evidence for human descent. Springer-Verlag, New York.

Kumar, S., and S. B. Hedges.

1998

. A molecular timescale for vertebrate evolution.

Nature

392

:

917

-920.

Kumar, S., K. Tamura, I. B. Jakobsen, and M. Nei.

2001

. MEGA2: molecular evolutionary genetics analysis software.

Bioinformatics

17

:

1244

-1245.

Li, W.-H., M. Gouy, P. M. Sharp, C. O'hUigin, and Y. W. Yang.

1990

. Molecular phylogeny of rodentia, lagomorpha, primates, artiodactyla, and carnivora and molecular clocks.

Proc. Natl. Acad. Sci. USA

87

:

6703

-6707.

Murphy, W. J., E. Eizirik, S. J. O'Brien, O. Madsen, M. Scally, C. J. Douady, E. Teeling, O. A. Ryder, M. J. Stanhope, and W. W. de Jong, et al. (11 co-authors).

2001

. Resolution of the early placental mammal radiation using Bayesian phylogenetics.

Science

294

:

2348

-2351.

Nei, M., and S. Kumar.

2000

. Molecular evolution and phylogenetics. Oxford University Press, Oxford.

Nei, M., and G. V. Glazko.

2002

. Estimation of divergence times for a few mammalian and several primate species.

J. Hered.

93

:

157

-164.

Nei, M., P. Xu, and G. Glazko.

2001

. Estimation of divergence times from multiprotein sequences for a few mammalian species and several distantly related organisms.

Proc. Natl. Acad. Sci. USA

98

:

2497

-2502.

O'hUigin, C., Y. Satta, N. Takahata, and J. Klein.

2002

. Contribution of homoplasy and of ancestral polymorphism to the evolution of genes in anthropoid primates.

Mol. Biol. Evol.

19

:

1501

-1513.

Pilbeam, D., M. D. Rose, J. C. Barry, and S. M. Shah.

1990

. New Sivapithecus humeri from Pakistan and the relationship of Sivapithecus and Pongo.

Nature

348

:

237

-239.

Rambaut, A., and L. Bromham.

1998

. Estimating divergence dates from molecular sequences.

Mol. Biol. Evol.

15

:

442

-448.

Reyes, A., G. Pesole, and C. Saccone.

2000

. Long-branch attraction pheonomenon and the impact of among-site rate variation on rodent phylogeny.

Gene

259

:

177

-187.

Rodriguez-Trelles, F., R. Tarrio, and F. J. Ayala.

2002

. A methodological bias toward overestimation of molecular evolutionary time scales.

Proc. Natl. Acad. Sci. USA

99

:

8112

-8115.

Satta, Y., J. Klein, and N. Takahata.

2000

. DNA archives and our nearest relative: the trichotomy problem revisited.

Mol. Phylogenet. Evol.

14

:

259

-275.

Saitou, N., and M. Nei.

1986

. The number of nucleotides required to determine the branching order of three species, with special reference to the human-chimpanzee-gorilla divergence.

J. Mol. Evol.

24

:

189

-204.

Saitou, N., and M. Nei.

1987

. The neighbor-joining method: a new method for reconstructing phylogenetic trees.

Mol. Biol. Evol.

4

:

406

-425.

Sanderson, M. J.

1997

. A nonparametric approach to estimating divergence times in the absence of rate constancy.

Mol. Biol. Evol.

14

:

1218

-1231.

Soltis, P. S., D. E. Soltis, V. Savolainen, P. R. Crane, and T. G. Barraclough.

2002

. Rate heterogeneity among lineages of tracheophytes: integration of molecular and fossil data and evidence for molecular living fossils.

Proc. Natl. Acad. Sci. USA

99

:

4430

-4435.

Stauffer, R. L., A. Walker, O. A. Ryder, M. Lyons-Weiler, and S. B. Hedges.

2001

. Human and ape molecular clocks and constraints on paleontological hypotheses.

J. Hered.

92

:

469

-474.

Swofford, D. L.

1998

. PAUP*: phylogenetic analysis using parsimony (*and other methods). Sinauer Associates, Sunderland, Mass.

Takahashi, K., A. P. Rooney, and M. Nei.

2000

. Origins and divergence times of mammalian class II MHC gene clusters.

J. Hered.

19

:

189

-204.

Takahata, N., and Y. Satta.

1997

. Evolution of the primate lineage leading to modern humans: phylogenetic and demographic inferences from DNA sequences.

Proc. Natl. Acad. Sci. USA

94

:

4811

-4815.

Takezaki, N., A. Rzhetsky, and M. Nei.

1995

. Phylogenetic test of the molecular clock and linearized trees.

Mol. Biol. Evol.

12

:

823

-833.

Ward, S.

1997

. The taxonomy and phylogenetic relationships of Sivapithecus revisited. Pp. 269–290 in Begun, D. R., C. V. Ward, and M. D. Rose, eds. Function, phylogeny and fossils: Miocene hominoid evolution and adaptation. New York, Plenum.

Wray, G. A., J. S. Levinton, and L. H. Shapiro.

1996

. Molecular evidence for deep precambrian divergences among metazoan phyla.

Science

274

:

568

-573.

Yoder, A. D., and Z. Yang.

2000

. Estimation of primate speciation dates using local molecular clocks.

Mol. Biol. Evol.

17

:

1081

-1090.

Society for Molecular Biology and Evolution